五阶最小二乘+迭代方法曲线拟合

在最小二乘的基础上，对误差进行拟合，将误差拟合的系数回代，逼近真实值，这个方法是错的，只是按客户要求做一下。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#define M 1000
#define N 6
double A[N]= {1.0,2.0,3.0,4.0,5.0,6.0},A_fit[N],A_e_fit[N];
double X[M],Y[M],YY[M],Y_e[M];//原始数据
const char PATH[64]="C:/Users/ren/Desktop/j2453/error/";
/*** 说明：5阶拟合函数* x：x坐标值* y：y坐标值* num：数据个数* a,b,c,d,e,f：返回拟合系数*/
void Fit5(double *x, double *y, int num,double *a, double *b, double *c, double *d, double *e, double *f)
{double sum_x = 0;double sum_x2 = 0;double sum_x3 = 0;double sum_x4 = 0;double sum_x5 = 0;double sum_x6 = 0;double sum_x7 = 0;double sum_x8 = 0;double sum_x9 = 0;double sum_x10 = 0;double sum_y = 0;double sum_xy = 0;double sum_x2y = 0;double sum_x3y = 0;double sum_x4y = 0;double sum_x5y = 0;for (int i = 0; i < num; ++i){sum_x += x[i];sum_x2 += x[i] * x[i];sum_x3 += x[i] * x[i] * x[i];sum_x4 += x[i] * x[i] * x[i] * x[i];sum_x5 += x[i] * x[i] * x[i] * x[i] * x[i];sum_x6 += x[i] * x[i] * x[i] * x[i] * x[i] * x[i];sum_x7 += x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i];sum_x8 += x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i];sum_x9 += x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i];sum_x10 += x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i] * x[i];sum_y += y[i];sum_xy += x[i] * y[i];sum_x2y += x[i] * x[i] * y[i];sum_x3y += x[i] * x[i] * x[i] * y[i];sum_x4y += x[i] * x[i] * x[i] * x[i] * y[i];sum_x5y += x[i] * x[i] * x[i] * x[i] * x[i] * y[i];}sum_x /= num;sum_x2 /= num;sum_x3 /= num;sum_x4 /= num;sum_x5 /= num;sum_x6 /= num;sum_x7 /= num;sum_x8 /= num;sum_x9 /= num;sum_x10 /= num;sum_y /= num;sum_xy /= num;sum_x2y /= num;sum_x3y /= num;sum_x4y /= num;sum_x5y /= num;double m01 = (sum_x8 * sum_x10 - sum_x9 * sum_x9) / (sum_x9 * sum_x10);double n01 = (sum_x7 * sum_x10 - sum_x8 * sum_x9) / (sum_x9 * sum_x10);double i01 = (sum_x6 * sum_x10 - sum_x7 * sum_x9) / (sum_x9 * sum_x10);double j01 = (sum_x5 * sum_x10 - sum_x6 * sum_x9) / (sum_x9 * sum_x10);double k01 = (sum_x4 * sum_x10 - sum_x5 * sum_x9) / (sum_x9 * sum_x10);double l01 = (sum_x4y * sum_x10 - sum_x5y * sum_x9) / (sum_x9 * sum_x10);double m02 = (sum_x7 * sum_x10 - sum_x9 * sum_x8) / (sum_x8 * sum_x10);double n02 = (sum_x6 * sum_x10 - sum_x8 * sum_x8) / (sum_x8 * sum_x10);double i02 = (sum_x5 * sum_x10 - sum_x7 * sum_x8) / (sum_x8 * sum_x10);double j02 = (sum_x4 * sum_x10 - sum_x6 * sum_x8) / (sum_x8 * sum_x10);double k02 = (sum_x3 * sum_x10 - sum_x5 * sum_x8) / (sum_x8 * sum_x10);double l02 = (sum_x3y * sum_x10 - sum_x5y * sum_x8) / (sum_x8 * sum_x10);double m03 = (sum_x6 * sum_x10 - sum_x9 * sum_x7) / (sum_x7 * sum_x10);double n03 = (sum_x5 * sum_x10 - sum_x8 * sum_x7) / (sum_x7 * sum_x10);double i03 = (sum_x4 * sum_x10 - sum_x7 * sum_x7) / (sum_x7 * sum_x10);double j03 = (sum_x3 * sum_x10 - sum_x6 * sum_x7) / (sum_x7 * sum_x10);double k03 = (sum_x2 * sum_x10 - sum_x5 * sum_x7) / (sum_x7 * sum_x10);double l03 = (sum_x2y * sum_x10 - sum_x5y * sum_x7) / (sum_x7 * sum_x10);double m04 = (sum_x5 * sum_x10 - sum_x9 * sum_x6) / (sum_x6 * sum_x10);double n04 = (sum_x4 * sum_x10 - sum_x8 * sum_x6) / (sum_x6 * sum_x10);double i04 = (sum_x3 * sum_x10 - sum_x7 * sum_x6) / (sum_x6 * sum_x10);double j04 = (sum_x2 * sum_x10 - sum_x6 * sum_x6) / (sum_x6 * sum_x10);double k04 = (sum_x * sum_x10 - sum_x5 * sum_x6) / (sum_x6 * sum_x10);double l04 = (sum_xy * sum_x10 - sum_x5y * sum_x6) / (sum_x6 * sum_x10);double m05 = (sum_x4 * sum_x10 - sum_x9 * sum_x5) / (sum_x5 * sum_x10);double n05 = (sum_x3 * sum_x10 - sum_x8 * sum_x5) / (sum_x5 * sum_x10);double i05 = (sum_x2 * sum_x10 - sum_x7 * sum_x5) / (sum_x5 * sum_x10);double j05 = (sum_x * sum_x10 - sum_x6 * sum_x5) / (sum_x5 * sum_x10);double k05 = (sum_x10 - sum_x5 * sum_x5) / (sum_x5 * sum_x10);double l05 = (sum_y * sum_x10 - sum_x5y * sum_x5) / (sum_x5 * sum_x10);double n11 = (m01 * n02 - n01 * m02) / (m01 * m02);double i11 = (m01 * i02 - i01 * m02) / (m01 * m02);double j11 = (m01 * j02 - j01 * m02) / (m01 * m02);double k11 = (m01 * k02 - k01 * m02) / (m01 * m02);double l11 = (m01 * l02 - l01 * m02) / (m01 * m02);double n12 = (m01 * n03 - n01 * m03) / (m01 * m03);double i12 = (m01 * i03 - i01 * m03) / (m01 * m03);double j12 = (m01 * j03 - j01 * m03) / (m01 * m03);double k12 = (m01 * k03 - k01 * m03) / (m01 * m03);double l12 = (m01 * l03 - l01 * m03) / (m01 * m03);double n13 = (m01 * n04 - n01 * m04) / (m01 * m04);double i13 = (m01 * i04 - i01 * m04) / (m01 * m04);double j13 = (m01 * j04 - j01 * m04) / (m01 * m04);double k13 = (m01 * k04 - k01 * m04) / (m01 * m04);double l13 = (m01 * l04 - l01 * m04) / (m01 * m04);double n14 = (m01 * n05 - n01 * m05) / (m01 * m05);double i14 = (m01 * i05 - i01 * m05) / (m01 * m05);double j14 = (m01 * j05 - j01 * m05) / (m01 * m05);double k14 = (m01 * k05 - k01 * m05) / (m01 * m05);double l14 = (m01 * l05 - l01 * m05) / (m01 * m05);double i21 = (n11 * i12 - i11 * n12) / (n11 * n12);double j21 = (n11 * j12 - j11 * n12) / (n11 * n12);double k21 = (n11 * k12 - k11 * n12) / (n11 * n12);double l21 = (n11 * l12 - l11 * n12) / (n11 * n12);double i22 = (n11 * i13 - i11 * n13) / (n11 * n13);double j22 = (n11 * j13 - j11 * n13) / (n11 * n13);double k22 = (n11 * k13 - k11 * n13) / (n11 * n13);double l22 = (n11 * l13 - l11 * n13) / (n11 * n13);double i23 = (n11 * i14 - i11 * n14) / (n11 * n14);double j23 = (n11 * j14 - j11 * n14) / (n11 * n14);double k23 = (n11 * k14 - k11 * n14) / (n11 * n14);double l23 = (n11 * l14 - l11 * n14) / (n11 * n14);double j31 = (i21 * j22 - j21 * i22) / (i21 * i22);double k31 = (i21 * k22 - k21 * i22) / (i21 * i22);double l31 = (i21 * l22 - l21 * i22) / (i21 * i22);double j32 = (i21 * j23 - j21 * i23) / (i21 * i23);double k32 = (i21 * k23 - k21 * i23) / (i21 * i23);double l32 = (i21 * l23 - l21 * i23) / (i21 * i23);double k33 = (j31 * k32 - k31 * j32) / (j31 * j32);double l33 = (j31 * l32 - l31 * j32) / (j31 * j32);*f = l33 / k33;*e = (l32 - k32 * (*f)) / j32;*d = (l23 - k23 * (*f) - j23 * (*e)) / i23;*c = (l14 - k14 * (*f) - j14 * (*e) - i14 * (*d)) / n14;*b = (l05 - k05 * (*f) - j05 * (*e) - i05 * (*d) - n05 * (*c)) / m05;*a = (sum_y - (*f) - sum_x * (*e) - sum_x2 * (*d) - sum_x3 * (*c) - sum_x4 * (*b)) / sum_x5;
}double Fun_N(double x,const int n,double *a)
{double ret=0.0;for(int i=0; i<n; i++){ret+=(pow(x,i)*a[i]);}return ret;
}
void save_erro(int x,double err[],int n)
{FILE* f;char path[64],num[16];strcpy(path,PATH);itoa(x,num,10);strcat(path,num);strcat(path,".txt");f = fopen(path, "w");for(int i=0; i<n; i++){fprintf(f,"%f\r",err[i]);}fclose(f);
}
int main()
{//给数据for(int i=0; i<M; i++){X[i]=0.02*i;Y[i]=Fun_N(X[i],N,A);YY[i]=Y[i]+0.0003*(double)(rand()%100-50);}//拟合Fit5(X,YY,M,&A_fit[5],&A_fit[4],&A_fit[3],&A_fit[2],&A_fit[1],&A_fit[0]);for(int i=0; i<50; i++){//求误差相for(int j=0; j<M; j++){double fit_val=Fun_N(X[j],N,A_fit);Y_e[j]=YY[j]-fit_val;}save_erro(i,Y_e,M);//误差行你和Fit5(X,Y_e,M,&A_e_fit[5],&A_e_fit[4],&A_e_fit[3],&A_e_fit[2],&A_e_fit[1],&A_e_fit[0]);//拟合系数回代for(int j=0; j<N; j++){A_fit[j]=A_fit[j]+A_e_fit[j];}}return 0;
}

五阶最小二乘+迭代方法曲线拟合相关推荐

说说牛顿迭代 -- 方法篇
说说牛顿迭代 – 方法篇写这个笔记主要是最近老在考虑最优化问题.今天刚好发现一个不错的手写公式的工具,加上前几天又发现Win10的Windows Ink比我想象得好用,于是来描几笔.主要是想试试这样 ...
协方差公式性质证明过程_论文推荐 | 刘志平：等价条件平差模型的方差-协方差分量最小二乘估计方法...
<测绘学报> 构建与学术的桥梁拉近与权威的距离等价条件平差模型的方差-协方差分量最小二乘估计方法刘志平1, 朱丹彤1, 余航1, 张克非1,2 1. 中国矿业大学环境与测绘学院, 江 ...
【论文速递】ISPRS2018 :基于增强极线几何约束以及自适应窗最小二乘匹配方法的立体SAR山区DSM
[论文速递]ISPRS2018 :基于增强极线几何约束以及自适应窗最小二乘匹配方法的立体SAR山区DSM [论文原文]:Radargrammetric DSM generation in mounta ...
机器学习笔记——2 简单线性模型及局部加权线性模型的基本原理和python实现（参数估计的两个基本角度：几何直观和概率直观。函数最值问题的两大基本算法：梯度方法与迭代方法）
简单线性模型及局部加权线性模型的基本原理和python实现(参数估计的两个基本角度:几何直观和概率直观.函数最值问题的两大基本算法:梯度方法与迭代方法) 线性模型是什么? 线性模型是监督学习的各类学习 ...
javascript 数组对象中的迭代方法
/* javascript 数组对象中的迭代方法 * ECMAScript5为数组定义了5个迭代方法.每个方法都接受两个参数,第一个是进行迭代的函数,第二个是该函数的作用域对象[可选]. * 进行迭代 ...
JavaScript数据迭代方法差别
js有很多总迭代方法,ES6之后又新增了几个: 这里主要讨论数组迭代遍历的方法所以不会细讲for...in... ES5.ES6数组迭代方法有: forEach map filter some eve ...
C语言试题二十之利用以下的简单迭代方法求方程cos(x)-x=0的一个实根。
1. 题目编写函数function,它的功能是:利用以下的简单迭代方法求方程cos(x)-x=0的一个实根. 迭代步骤如下: (1)取x1初值为0.0: (2)x0=x1,把x1的值赋各x0; (3 ...
JavaScript数组(2)---遍历/迭代方法 8种
最近工作中经常涉及到数据的处理,数组尤其常见,经常需要对其进行遍历.转换操作,网上的文章零零散散,不得已自己又找出红宝书来翻出来看,顺便记一笔,便于以后查询. 数组常用的方法 ECMAScript5为 ...
ols线性回归_普通最小二乘[OLS]方法使用于机器学习的简单线性回归变得容易
ols线性回归 Hello Everyone! 大家好! I am super excited to be writing another article after a long time sinc ...

五阶最小二乘+迭代方法曲线拟合

五阶最小二乘+迭代方法曲线拟合相关推荐

最新文章

热门文章