用C语言实现简单的多元线性回归算法（二）

上一篇我们贴上了简单粗暴的线性回归的代码，里面各种参数都设置的是固定参数，不具有可扩展性，今天我们在上一篇的基础上做了部分改进，当然对于熟悉C语言的大侠来说可能这篇博客会太low了，您完全可以跳过。我们在这里只是讲如何用C语言自己实现一个有实用性的线性回归。有人会说用python不是很简单吗，干嘛费劲巴拉的用C语言实现。首先这里是为了满足某些只支持C语言的环境，例如IOT的部分设备；另外也是为了加深对算法的理解和思考。至于线性回归原理网上到处是，我们在这里只贴一下参数更新的推导过程公式，

线性回归模型的代价函数是下面这样子：

其中，f(x)为：

我们的目标是求最小J(w)

即在当J(w)取最小值时，对应的W的值，这时候梯度下降就用上了，用公式可表示为如下样式：

我们把J(w)对w的偏导求出来：

把(2)带入（1）可得：

以上就是线性回归中的梯度下降算法。

下面我们贴上代码，代码中会有比较详细的注释，在这里就不多废话了。

MultipleLinearRegressionTest.h：

//
// Created by lenovo on 2020/7/26.
//

#ifndef MULTIPLELINEARREGRESSION_MULTIPLELINEARREGRESSIONTEST_H
#define MULTIPLELINEARREGRESSION_MULTIPLELINEARREGRESSIONTEST_H

//初始化参数
void init_test(double learning_rate, long int X_Len,long int X_arg_count,int channel,double *y_pred_pt,double* theta_pt,double* temp_pt);

//开始训练
void fit_test(double *X_pt,double *Y_pt);

#endif //MULTIPLELINEARREGRESSION_MULTIPLELINEARREGRESSIONTEST_H

MultipleLinearRegressionTest.c：

//
// Created by lenovo on 2020/7/26.
//


#include <stdio.h>
#include "MultipleLinearRegressionTest.h"


// 样本个数，这里其实是bachsize的意思，只不过我们只有一个bach
long int g_data_len = 10;
//学习率
double g_learning_rate = 0.01;
//参数个数
long int g_param_count = 6;
//通道数
int g_channel = 1;


double g_out_Y_pt = 0;
//预测输出值指针
double *g_pred_pt = 0;
//参数指针
double *g_theta_pt = 0;

double *temp  = 0;


//损失值
double loss_val = 1.0;


//初始化各个参数
void init_test(double learning_rate, long int data_len,long int param_count,int channel,double *y_pred_pt,double* theta_pt,double*temp_pt){
    g_learning_rate = learning_rate;
    g_data_len = data_len;
    g_param_count = param_count;
    g_channel = channel;
    g_pred_pt = y_pred_pt;

    g_theta_pt = theta_pt;
    temp = temp_pt;
}


//更新参数值
void train_step_test(const double *temp_step, double *theta) {
    for (int i = 0; i < g_param_count; i++) {
        //这里的temp_step[i] 其实就是算出来的梯度值
        theta[i] = theta[i] - temp_step[i];
        printf(" theta[%d] = %f\n",i, theta[i]);
    }
}


//计算一个样本的输出值  y = WX + b
double f_test(const double *train_x_item,const double *theta){

    g_out_Y_pt = -1;
    for (int i = 0; i < g_param_count; ++i) {
        g_out_Y_pt += theta[i]*(*(train_x_item+i));

    }
    return g_out_Y_pt;
}


//预测
double* predict_test(double *train_x,double *theta){
    for (int i = 0; i < g_data_len; ++i) {
        g_pred_pt[i] = f_test(train_x+i*g_param_count,theta);
    }
    return g_pred_pt;
}

//求损失 
double loss_test(double *train_x,const double *train_y,
                 double *theta,double loss_val){
    predict_test(train_x,theta);
    loss_val = -1;
    for (int i = 0; i < g_data_len; i++) {
        loss_val += (train_y[i] - g_pred_pt[i]) * (train_y[i] - g_pred_pt[i]);
    }
    //次数计算的是一个bachsize 的平均损失(避免个别离群值引起过拟合)
    loss_val = loss_val / g_data_len;

    printf(" loss_val = %f\n", loss_val);
    return loss_val;
}

//求梯度
void calc_gradient_test(const double *train_x, const double *train_y,
                        double *temp_step, double *theta) {

    //对每个参数分别子梯度
    for (int i = 0; i < g_param_count; i++) {
        double sum = 0;
        //用一个bachsize的平均梯度表示当前参数的梯度(temp_step[i])
        for (int j = 0; j < g_data_len; j++) {
            double hx = 0;
            for (int k = 0; k < g_param_count; k++) {
                hx += theta[k] * (*(train_x+j*g_param_count+k));
            }
            sum += (hx - train_y[j]) * (*(train_x+j*g_param_count+i));
        }
        temp_step[i] = sum / g_data_len * 0.01;
    }
    printf("--------------------\n");
    train_step_test(temp_step, theta);
}

// 开始训练 (次数迭代60000次)
void fit_test(double *X_pt,double *Y_pt) {
    for (int i = 0; i < 60000; ++i) {
        printf("\nstep %d: \n", i);
        calc_gradient_test(X_pt,Y_pt,temp,g_theta_pt);
        loss_val = loss_test(X_pt,Y_pt,g_theta_pt,loss_val);
    }
}

main.c :

#include "src/MultipleLinearRegressionTest.h"

int main() {



    //该数据由 Y=              4*X1 + 9*X2 + 10*X3 + 2*X4 + 1*X5 + 6 产生(数据中尚为加入噪声)
    //数据总的第一列都设置为1，是为了表示截距(Y=WX+B)的参数
    double X_pt[10][6] = {{1, 7.41,  3.98,  8.34,   8.0,    0.95},
                          {1, 6.26,  5.12,  9.57,   0.3,    7.79},
                          {1, 1.52,  1.95,  4.01,   7.96,   2.19},
                          {1, 1.91,  8.58,  6.64,   2.99,   2.18},
                          {1, 2.2,   6.88,  0.88,   0.5,    9.74},
                          {1, 5.17,  0.14,  4.09,   9.0,    2.63},
                          {1, 9.13,  5.54,  6.36,   9.98,   5.27},
                          {1, 1.17,  4.67,  9.02,   5.14,   3.46},
                          {1, 3.97,  6.72,  6.12,   9.42,   1.43},
                          {1, 0.27,  3.16,  7.07,   0.28,   1.77}};

    double Y_pt[10] = {171.81, 181.21, 87.84, 165.42, 96.26, 89.47, 181.21, 156.65, 163.83, 108.55};
    double pred_pt[10] ={1,2,3,4,5,6,7,8,9,0};

    // 参数
    double theta_pt[6] = {1.0,1.0,1.0,1.0,1.0,2.0};
    //临时变量
    double temp_pt[6] =  {1.0,1.0,1.0,1.0,1.0,2.0};
    // 学习率(次数我们还是指定为一个固定值)
    double learning_rate = 0.02;

    // 样本个数
    int data_len = 10;
    // 参数个数
    int param_count = 6;
    //通道数
    int channel = 1;


    init_test(learning_rate, data_len,param_count,channel,pred_pt,theta_pt,temp_pt);

    fit_test((double *) X_pt, Y_pt);

    return 0;
}

以上便是线性归回模型梯度下降求解的代码实现。各位大侠如有疑问或者问题欢迎交流。