C++ Linear Regression: A Complete Guide

Linear regression is a fundamental statistical method used to model the relationship between a dependent variable and one or more independent variables. While Python and R are popular for data science, C++ offers unmatched performance for high-frequency trading, embedded systems, and real-time analytics. In this guide, you'll learn how to calculate linear regression using C++ from scratch, covering both the closed-form solution and gradient descent.

Why Implement Linear Regression in C++?

C++ is not typically the first language for data science, but it shines when speed and memory efficiency are critical. A c++ regression implementation can be embedded into production systems, game engines, or IoT devices. Moreover, understanding the math behind linear regression deepens your grasp of machine learning algorithms.

The Mathematics of Linear Regression

Simple linear regression models the relationship as:

y = β₀ + β₁x + ε

where β₀ is the intercept, β₁ is the slope, and ε is the error term. For multiple regression, we extend to:

y = Xβ + ε

We estimate β using the least squares method, which minimizes the sum of squared residuals. The closed-form solution is:

β = (XᵀX)⁻¹Xᵀy

This is known as the normal equation. Alternatively, we can use gradient descent for large datasets.

C++ Least Squares Implementation

Let's implement simple linear regression using the least squares approach. We'll compute the slope and intercept directly from the data.

#include <iostream>
#include <vector>
#include <numeric>
#include <cmath>

class LinearRegression {
private:
    double slope;
    double intercept;
public:
    LinearRegression() : slope(0), intercept(0) {}
    
    void fit(const std::vector<double>& x, const std::vector<double>& y) {
        size_t n = x.size();
        double sum_x = std::accumulate(x.begin(), x.end(), 0.0);
        double sum_y = std::accumulate(y.begin(), y.end(), 0.0);
        double sum_xy = 0.0, sum_xx = 0.0;
        for (size_t i = 0; i < n; ++i) {
            sum_xy += x[i] * y[i];
            sum_xx += x[i] * x[i];
        }
        double denom = n * sum_xx - sum_x * sum_x;
        slope = (n * sum_xy - sum_x * sum_y) / denom;
        intercept = (sum_y - slope * sum_x) / n;
    }
    
    double predict(double x) const {
        return intercept + slope * x;
    }
};

This implementation uses the standard formulas for slope and intercept. It's efficient and works well for small to medium datasets.

Multiple Linear Regression with Gradient Descent

For multiple features or large datasets, gradient descent is more practical. Here's a C++ class for multiple linear regression using gradient descent:

#include <vector>
#include <cmath>

class MultipleLinearRegression {
private:
    std::vector<double> weights;
    double bias;
    double learning_rate;
    int epochs;
public:
    MultipleLinearRegression(double lr = 0.01, int ep = 1000) 
        : learning_rate(lr), epochs(ep), bias(0) {}
    
    void fit(const std::vector<std::vector<double>>& X, const std::vector<double>& y) {
        size_t n_samples = X.size();
        size_t n_features = X[0].size();
        weights.assign(n_features, 0.0);
        
        for (int epoch = 0; epoch < epochs; ++epoch) {
            std::vector<double> dw(n_features, 0.0);
            double db = 0.0;
            for (size_t i = 0; i < n_samples; ++i) {
                double y_pred = predict(X[i]);
                double error = y_pred - y[i];
                for (size_t j = 0; j < n_features; ++j) {
                    dw[j] += error * X[i][j];
                }
                db += error;
            }
            for (size_t j = 0; j < n_features; ++j) {
                weights[j] -= learning_rate * dw[j] / n_samples;
            }
            bias -= learning_rate * db / n_samples;
        }
    }
    
    double predict(const std::vector<double>& x) const {
        double result = bias;
        for (size_t i = 0; i < weights.size(); ++i) {
            result += weights[i] * x[i];
        }
        return result;
    }
};

This code implements batch gradient descent. You can adjust the learning rate and number of epochs for convergence.

Evaluating Your Model

After training, evaluate performance using metrics like Mean Squared Error (MSE) or R-squared. Here's a quick MSE function:

double mean_squared_error(const std::vector<double>& y_true, const std::vector<double>& y_pred) {
    double sum = 0.0;
    for (size_t i = 0; i < y_true.size(); ++i) {
        sum += std::pow(y_true[i] - y_pred[i], 2);
    }
    return sum / y_true.size();
}

Practical Considerations

• Feature Scaling: Gradient descent converges faster when features are normalized.
• Regularization: Add L1 or L2 penalties to prevent overfitting.
• Matrix Libraries: For large-scale problems, use Eigen or Armadillo for efficient matrix operations.
• Numerical Stability: The normal equation can be unstable if XᵀX is ill-conditioned; use QR decomposition instead.

Conclusion

Implementing c++ linear regression from scratch is a rewarding exercise that combines programming with statistical learning. Whether you choose the closed-form c++ least squares or gradient descent, you now have a solid foundation. For production, consider using optimized libraries, but understanding the core algorithm is invaluable. Happy coding!