在此之前我們再回顧下上一篇文章中的系數求解公式。二次拟合曲線的方程為:
系數行列式為:
另外:
所求系數為:
使用C實現最小二乘法為:
#include<stdio.h>
#include <math.h>
#include<stdbool.h>
#define DATA_NUM (6)
#define DOUBLE_PRECISION (1e-15)
double x[DATA_NUM]={0,2,4,6,8,10};
double y[DATA_NUM]={0,6,25,42,70,110};
//y = a*x^2 b*x c
bool LeastSquares(double *x, double *y, unsigned int data_num, double *a, double *b, double *c) {
double sumx=0,sumx2=0,sumx3=0,sumx4=0,sumy=0,sumxy=0,sumx2y=0;
double D=0;
if(!data_num) return false;
for(int i=0;i<data_num;i ) {
sumx =x[i];sumy =y[i];
sumx2 =pow (x[i],2); sumxy =x[i]*y[i];
sumx3 =pow(x[i],3); sumx2y =pow(x[i],2)*y[i];
sumx4 =pow(x[i],4);
}
D = sumx2*sumx2*sumx2 sumx*sumx*sumx4 data_num*sumx3*sumx3 - data_num*sumx2*sumx4 - 2*sumx*sumx2*sumx3;
if (fabs(D) < DOUBLE_PRECISION) {
return false;
}
*a = (sumy*(sumx2*sumx2-sumx*sumx3) sumxy*(data_num*sumx3-sumx*sumx2) sumx2y*(sumx*sumx-data_num*sumx2))/D;
*b = (sumy*(sumx*sumx4-sumx2*sumx3) sumxy*(sumx2*sumx2-data_num*sumx4) sumx2y*(data_num*sumx3-sumx*sumx2))/D;
*c = (sumy*(sumx3*sumx3-sumx2*sumx4) sumxy*(sumx*sumx4-sumx2*sumx3) sumx2y*(sumx2*sumx2-sumx*sumx3))/D;
return true;
}
int main() {
double a,b,c;
LeastSquares(x, y, DATA_NUM, &a, &b, &c);
printf("a=%9.8f,\nb=%9.8f,\nc=%9.8f\n",a,b,c);
printf ("y=%9.6fx*x %9.6fx %9.6f",a,b,c);
return 0;
}
上述在求解行列式值的時候使用的是代數餘子式方法求解。即用某一個列與代數餘子式相乘。
往期推薦
最小二乘法系數求解
最小二乘法曲線拟合推導過程
,更多精彩资讯请关注tft每日頭條,我们将持续为您更新最新资讯!