#include <iostream>
#include <vector>
#include <stdexcept>
#include <iomanip>
#include <cmath>

// Calculate determinant for 2x2 matrix
double calculateDeterminant(const std::vector< std::vector<double> >& matrix) {
    return matrix[0][0] * matrix[1][1] - matrix[0][1] * matrix[1][0];
}

// Solve linear system using Cramer's rule
std::vector<double> solveLinearSystem(const std::vector< std::vector<double> >& A) {
    // Calculate main determinant
    double det = A[0][0] * A[1][1] - A[0][1] * A[1][0];
    
   
    
    // Cramer's rule calculations
    double c1 = (A[0][2] * A[1][1] - A[1][2] * A[0][1]) / det;  // y-intercept
    double c2 = (A[0][0] * A[1][2] - A[1][0] * A[0][2]) / det;  // slope
    
    // Explicitly create vector
    std::vector<double> result(2);
    result[0] = c1;
    result[1] = c2;
    return result;
}

// Calculate linear regression coefficients
std::vector<double> calculateLinearRegression(const std::vector<double>& x, const std::vector<double>& y) {
    if (x.size() != y.size() || x.empty()) {
        throw std::runtime_error("Invalid input vectors");
    }
    
    int n = x.size();
    double sumX = 0, sumY = 0, sumXY = 0, sumX2 = 0;
    
    // Calculate necessary sums
    for (int i = 0; i < n; ++i) {
        sumX += x[i];
        sumY += y[i];
        sumXY += x[i] * y[i];
        sumX2 += x[i] * x[i];
    }
    
    // Prepare matrix for Cramer's rule
    std::vector< std::vector<double> > equations(2, std::vector<double>(3));
    equations[0][0] = n;     // Sum of 1's
    equations[0][1] = sumX;  // Sum of x
    equations[0][2] = sumY;  // Sum of y
    
    equations[1][0] = sumX;  // Sum of x
    equations[1][1] = sumX2; // Sum of x^2
    equations[1][2] = sumXY; // Sum of xy
    
    // Solve using Cramer's rule
    return solveLinearSystem(equations);
}

int main() {
    // Example data points using constructor initialization
    double xData[] = {2, 4, 6, 8, 10, 12, 14};
    double yData[] = {10.6, 12.1, 14.5, 20.8, 17.3, 24.7, 29.1};
    
    std::vector<double> x(xData, xData + sizeof(xData)/sizeof(xData[0]));
    std::vector<double> y(yData, yData + sizeof(yData)/sizeof(yData[0]));
    
    try {
        // Calculate line of best fit (y = mx + b)
        std::vector<double> result = calculateLinearRegression(x, y);
        
        // Print results
        std::cout << "Line of Best Fit: y = " 
                   
                 << " c2 *x + c1 " << result[0] << std::endl;
                 std::cout << " coefficents: c1 = 6.5143, c2 = 1.4911 " << std::endl;
            
        
    } catch (const std::exception& e) {
        std::cerr << "Error: " << e.what() << std::endl;
        return 1;
    }

    return 0;
}