【C言語】最小二乗法による単回帰と重回帰

【C言語】最小二乗法による単回帰と重回帰

最小二乗法とは

最小二乗法(Least Squares Method)は、与えられたデータ点に最も適した直線や超平面を求めるための方法です。

目標は、以下の誤差関数(残差平方和, Residual Sum of Squares: RSS)を最小化することです。

$$ RSS = \sum_{i=1}^{n} (y_i – \hat{y}_i)^2 $$

ここで、$\hat{y}_i$ は予測値、$y_i$ は実際のデータ点です。

単回帰の理論とC言語実装

単回帰モデル

単回帰(Simple Linear Regression)は、次のような形で表されます:

$$ y = \beta_0 + \beta_1 x + \epsilon $$

最小二乗法を用いて、$\beta_0$(切片)と$\beta_1$(傾き)を求めます。

$$ \beta_1 = \frac{\sum (x_i – \bar{x})(y_i – \bar{y})}{\sum (x_i – \bar{x})^2} $$

$$ \beta_0 = \bar{y} – \beta_1 \bar{x} $$

単回帰のC言語実装

以下は、単回帰をC言語で実装するコードです。


#include <stdio.h>

void linear_regression(double x[], double y[], int n, double *b0, double *b1) {
    double sum_x = 0, sum_y = 0, sum_xy = 0, sum_x2 = 0;
    for (int i = 0; i < n; i++) {
        sum_x += x[i];
        sum_y += y[i];
        sum_xy += x[i] * y[i];
        sum_x2 += x[i] * x[i];
    }
    double x_mean = sum_x / n;
    double y_mean = sum_y / n;
    
    *b1 = (sum_xy - n * x_mean * y_mean) / (sum_x2 - n * x_mean * x_mean);
    *b0 = y_mean - (*b1) * x_mean;
}

int main() {
    double x[] = {1, 2, 3, 4, 5};
    double y[] = {2, 4, 5, 4, 5};
    int n = 5;
    double b0, b1;

    linear_regression(x, y, n, &b0, &b1);
    printf("y = %.2f + %.2f * x\n", b0, b1);
    return 0;
}

重回帰の理論とC言語実装

重回帰モデル

重回帰(Multiple Linear Regression)は、以下のように表されます:

$$ y = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \dots + \beta_m x_m + \epsilon $$

通常、重回帰の係数は、行列を用いて求められます。

$$ \mathbf{b} = (\mathbf{X}^T \mathbf{X})^{-1} \mathbf{X}^T \mathbf{y} $$

重回帰のC言語実装

以下は、ガウスの消去法を用いた重回帰のC言語実装です。


#include <stdio.h>

#define N 3 // 説明変数の数(例として2つ)

void solve_linear_system(double A[N+1][N+2], double b[N+1]) {
    for (int i = 0; i <= N; i++) {
        double pivot = A[i][i];
        for (int j = i; j <= N+1; j++) {
            A[i][j] /= pivot;
        }
        for (int k = 0; k <= N; k++) {
            if (k != i) {
                double factor = A[k][i];
                for (int j = i; j <= N+1; j++) {
                    A[k][j] -= factor * A[i][j];
                }
            }
        }
    }
    for (int i = 0; i <= N; i++) {
        b[i] = A[i][N+1];
    }
}

int main() {
    double X[N+1][N+2] = {
        {1, 1, 2, 5},
        {1, 2, 3, 8},
        {1, 3, 5, 14}
    };
    double b[N+1];

    solve_linear_system(X, b);

    printf("回帰係数: ");
    for (int i = 0; i <= N; i++) {
        printf("%.2f ", b[i]);
    }
    printf("\n");
    return 0;
}

応用例

  • 株価予測
  • 住宅価格予測
  • 広告費と売上の関係分析

コメントは受け付けていません。