MATLAB在数据统计中的应用

MATLAB在数据统计中的应用______________________________________________目录：1、一元线性回归的matlab实现（含检验）【更新】2、一维数据滑动平均的matlab实现3、多元线性回归的matlab实现4、K阶自回归拟合及二阶自回归预测的Matlab实现5、一次指数平滑预测的matlab实现6、n次指数平滑及其预测7、一维数据移动平滑的matlab实现8、K阶自相关系数的matlab实现（含置信度检验）说明：1.正文中命令部分可以直接在Matlab中运行，作者(Yangfd09)在MATLABR2009a(7.8.0.347)中运行通过。2.限于作者水平问题，文中难免疏漏和错误，如蒙赐教，不胜感激！3.原创作品，仅供学习交流之用，会有不定期更新。一元线性回归的matlab实现（含检验）【更新】%求一元线性回归方程%数据要求：两行。第一行存放x的观察值，第二行存放y的观察值%数据文件名：data_yyhg.mat；变量名：test%loaddata_yyhg.matN=length(test(1,:));%注：也可以用[M,N]=size(test)%但不能用N=size(test(1,:))sx=0;sx2=0;sy=0;sy2=0;sxy=0;Lxy=0;Lyy=0;fori=1:Nsx=sx+test(1,i);sx2=sx2+test(1,i)^2;sy=sy+test(2,i);sy2=sy2+test(2,i)^2;sxy=sxy+test(1,i)*test(2,i);Lxy=Lxy+(test(1,i)-sum(test(1,:))/N)*(test(2,i)-sum(test(2,:)/N));Lyy=Lyy+(test(2,i)-sum(test(2,:))/N)^2;endr=[N,sx;sx,sx2]\[sy;sxy];a=r(1);b=r(2);%F分布检验U=b*Lxy;Q=Lyy-U;F=(N-2)*U/Q;%拟合优度检验x=test(1,:);y=a+b*x;eq=sum(test(2,:))/N;ssd=0;ssr=0;fori=1:Nssd=ssd+(test(2,i)-y(i))^2;ssr=ssr+(y(i)-eq)^2;endsst=ssd+ssr;RR=ssr/sst;%命令窗口中显示回归方程str=[blanks(5),'y=','(',num2str(a),')','+','(',num2str(b),')','*x'];disp('')disp('回归方程为:')disp(str)disp('R^2拟合优度检验：')strin=['R^2=',num2str(RR)];disp(strin)disp('F-分布显著性检验：')stri=['F计算值：',num2str(F),blanks(4),'自由度：f1=1,f2=',num2str(N-2)];disp(stri)disp('注：请对照F-分布表找到所需置信水平下的F临界值Fa，若FFa，则通过检验。')%绘制x-y散点图和回归直线yy=a+b*test(1,:);plot(test(1,:),test(2,:),'r.'),holdonplot(test(1,:),yy,'b-'),holdofftitle(str)附（可以直接粘贴到.mat文件中）：3.845.8811.314.416.516.213.810.86.74.777.751.260.154.155.456.84555.367.573.376.679.6一维数据滑动平均的matlab实现%滑动平均%数据格式：单行(按时间序列排序)%数据文件名：data_ph.mat,变量名：testloaddata_ph.matM=length(test);disp('请输入单侧平滑点数（时距）')k=input('（输入1对应于三点平滑，2对应五点平滑）：');y=zeros(1,M);if2*k+1=Mfori=1:M-2*kforj=i:i+2*ky(i+k)=y(i+k)+test(j);endy(i+k)=y(i+k)/(2*k+1);endy([1:k,M-k+1:M])=NaN;str=[int2str(k),'点滑动平均结果如下：'];disp(str)formatcompactdata=test,result=yformatelsedisp('Error:数据个数不足！')end附（直接复制到.mat文件中即可）：（某城市1999-2004年用水量数据）211.3260.18209.1248.79241250多元线性回归的matlab实现%多元线性回归%数据要求：M行N列。第一行为y的观察值，其余行分别为x1,x2,...,x(M-1)%数据文件名：data_dyhg.mat;变量名：testloaddata_dyhg.mat[M,N]=size(test);Y=test(1,:)';X=cat(2,ones(N,1),test(2:end,:)');A=X'*X;B=X'*Y;b=A\B;formatlonggstr=['系数按升幂排列为：',num2str(b')];disp(str)%F-分布检验eqx=sum(test(2:M,:),2)/N;eqy=sum(Y)/N;Liy=zeros(1,M-1);Lyy=0;fori=1:NLyy=Lyy+(Y(i)-eqy)^2;end%求回归平方和（方法一）：用y的拟合值y=zeros(1,N);UU=0;fori=1:Nforj=1:My(i)=y(i)+b(j)*X(i,j);endUU=UU+(y(i)-eqy)^2;endQ0=Lyy-UU;F0=(UU/(M-1))/(Q0/(N-M));%求回归平方和（方法二）：用偏回归系数fori=1:M-1forj=1:NLiy(i)=Liy(i)+(test(i+1,j)-eqx(i))*(test(1,j)-eqy);endendU=0;fori=1:M-1U=U+b(i+1)*Liy(i);endQ1=Lyy-U;F1=(U/(M-1))/(Q1/(N-M));%命令窗口中显示回归方程disp('')disp('F-分布检验：')stri='F计算值：';strin=['方法一：',num2str(F0),'方法二：',num2str(F1)];string=['自由度：','f1=',num2str(M-1),',f2=',num2str(N-M)];st=cat(2,stri,strin);disp(string)disp(st)disp('注：1.请对照F-分布表找到所需置信水平下的F临界值Fa，若FFa，则通过检验。')disp('2.方法一求回归平方和时用到了因变量的拟合值，方法二用的是偏回归系数。')附（数据较多，显示不便，直接复制到.mat文件中即可）：48.25193.72413.94358.6615.04752.42435.43238.5587.85316503.73554.04502.07611.78603.66501.67540.16264.15427.11513.09478.21395.25650.14480.2485.79144.3839.1765.3271.8858.5754.33106.33257.21114.53127.49194.9331.09110.57194.42163.89389.9541.5573.03521.31645.21466.28558.83621.02515.02545.72786.75584.8957440.536.635.5337.4835.4333.8235.6336.5739.7736.0534.235.3835.623434.3834.7334.5837.235.1235.4234.7335.8533.7533.441.841.5840.1540.2640.724040.339.3738.8339.1538.9338.7838.4538.6338.2337.9237.236.5835.7335.7335.1535.5232.9534.0334.73534.2135.4336.141170.81707.21908.82072.42136.4930.82025.11397.81477.21517.214101886.619173471.42314.612501131.72726.71765245014951873.79701079177021591138.7152615911270.21177.41332.22311.81453.71482.71764.6227113671976.11530.83045.11255.61421.91346.6136016501014.31753.22810.22915.73362.71221.21111.7K阶自回归拟合及二阶自回归预测的Matlab实现%k阶自回归及二阶自回归预测%数据格式：单行%数据文件名：data_zhg.mat,变量名：testclear,clcloaddata_zhg.matM=length(test);k=input('请输入自回归阶数：');y=zeros(k+1,M-k);fori=1:k+1y(i,:)=test(i:i+M-k-1);endy=flipud(y);[m,N]=size(y);Y=y(1,:)';X=cat(2,ones(N,1),y(2:end,:)');A=X'*X;B=X'*Y;b=A\B;formatlonggstr=['首项为常数项，其余系数按距离预测值远近排列为：',num2str(b')];disp(str)%二元自回归预测ifk==2Yuce2=b(1)+b(2)*test(M)+b(3)*test(M-1);str=['二阶自回归预测下一个值为：',num2str(Yuce2)];disp(str)end附（直接粘贴到.mat文件中即可）：3149.443303.663010.33109.613639.213253.83466.53839.93894.664009.614253.254101.54119.884258.654401.79一次指数平滑预测的matlab实现%一次指数平滑及其预测%数据格式：单行(按时间序列排序)%数据文件名：data_ph.mat,变量名：testloaddata_ph.matM=length(test);disp('参考：时间序列稳定、数据波动较小时，a取(0.05,0.3)；否则取(0.7,0.95)')a=input('请输入平滑系数a：');y=zeros(1,M);fori=1:Mforj=0:i-1y(i)=y(i)+test(i-j)*a*(1-a)^j;endendyy=a*test(M)+(1-a)*y(M);formatcompactformatshortdata=test,result=yformatstr=['下一时段数值预测：',num2str(yy)];disp(str)附（直接复制到.mat文件中即可）：（某城市1999-2004年用水量数据）211.3260.18209.1248.79241250n次指数平滑及其预测%数据格式：单行(按时间序列排序)%数据文件名：data_ph.mat,变量名：testloaddata_ph.matM=length(test);disp('参考：时间序列稳定、数据波动较小时，a取(0.05,0.3)；否则取(0.7,0.95)')a=input('请输入平滑系数a：');k=input('请输入平滑次数：');y0=test;y=zeros(1,M);fort=1:kfori=1:Mforj=0:i-1y(i)=y(i)+y0(i-j)*a*(1-a)^j;endendift~=ky0=y;endendyy=a*y0(M)+(