2、=0对数据点S（1,1）到S（1，M）区间按0.01步长作图计算δrms=M1δi2;δmax=max1<=i<=
3、δi
4、.end程序代码：%多项式最小二乘法拟合数据N=input('pleaseputinhowmanytimesthepowerwillyouoverfit:');M=input('howmanycouplesofstatisticsarethereinthetable:');%读入数据文件f=fopen('mytxt.txt','r');S=fscanf(f,'%g',[M2]);fcl
5、ose(f);S=S';%显示数据文件，确保正确输入disp('S(x,y)=');disp(S);%建立多项式系数法方程组中间矩阵C=zeros(N+1,M);fori=1:N+1forj=1:MifS(1,j)==0C(i,j)=0;elseC(i,j)=S(1,j).^(i-1);endendend%建立法方程组A=C*C';Y=zeros(M,1);fori=1:MY(i,1)=S(2,i);endb=C*Y;%用列主元高斯消元法接法方程组A=[A,b];fori=1:N+1max=abs(A(i
6、,i));forj=i+1:N+1ifabs(A(j,i))>maxflag=j;max=A(j,i);endendfork=i:N+2B=A(flag,k);A(flag,k)=A(i,k);A(i,k)=B;endforkh=i+1:N+1m=-A(kh,i)/A(i,i);A(kh,i)=0;forkl=i+1:N+2A(kh,kl)=A(kh,kl)+m*A(i,kl);endendendX=zeros(N+1,1);fori=N+1:-1:1forj=i-1:-1:1m=-A(j,i)/A(i,
7、i);A(j,N+2)=A(j,N+2)+m*A(i,N+2);endX(i,1)=A(i,N+2)/A(i,i);enddisp(X);%根据系数求得待定曲线symsx;expr=0;fori=1:N+1expr=expr+X(i,1)*x.^(i-1);end%输出得到的曲线表达式disp(expr);%计算偏差bias=zeros(M,1);forj=1:Mfori=1:N+1bias(j,1)=bias(j,1)+X(i,1)*S(1,j)^(i-1);endbias(j,1)=bias(j,1)
8、-S(2,j);end%寻找最大偏差max=abs(bias(1,1));flag=1;fori=2:Mifabs(bias(i,1))>maxflag=i;max=abs(bias(i,1));endenddisp('themaximunabsoulutevalueis:');disp(max);%计算均方误差rms=0;fori=1:Mrms=rms+bias(i,1)^2;endrms=sqrt(rms);disp('thesquarebiasis:');disp(rms);%制图a=S(1,1):
9、0.01:S(1,M);y=subs(expr,x,a);plot(a,y);holdon;gridon;fori=1:Mx=S(1,i);y=S(2,i);plot(x,y,'*');holdon;end运行结果：表达式中分式难以化简，但在表达式前给出了次幂前的四位有效数字的系数
