MATLAB在数据统计中的应用培训资料

1、Good is good, but better carries it.精益求精，善益求善。MATLAB在数据统计中的应用-MATLAB在数据统计中的应用_目录：HYPERLINKlm11、一元线性回归的matlab实现（含检验）【更新】HYPERLINKlm22、一维数据滑动平均的matlab实现HYPERLINKlm33、多元线性回归的matlab实现HYPERLINKlm44、K阶自回归拟合及二阶自回归预测的Matlab实现HYPERLINKlm55、一次指数平滑预测的matlab实现HYPERLINKlm66、n次指数平滑及其预测HYPERLINKlm77、一维数据移动平滑的matla

2、b实现HYPERLINKlm88、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,:);%注：也可

3、以用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

4、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(R2拟合优度检验：)str

5、in=R2=,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.

6、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+1Fa，则通过检验。)disp(2.方法一求回归平方和时用到了因变量的拟合值，方法二用的是偏回归系数。)附（数据较多，显示不便，直接复制到.

7、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.7

8、5584.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.113

9、97.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_z

10、hg.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=AB;formatlonggstr=首项为常数项，其余系数按距离预测值远近排列为：,num2str(b);disp(str)%二元自回归预测ifk=2Yuce2=b(1)+b(2)*test(M

11、)+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(参考：时间序列稳定、数据波

12、动较小时，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次指数平滑及其预测%

13、数据格式：单行(按时间序列排序)%数据文件名：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)+(1-a)*y(M);formatcompactformatshor

14、tdata=test,result=yformat附（直接复制到.mat文件中即可）：（某城市1999-2004年用水量数据）211.3260.18209.1248.79241250一维数据移动平滑的matlab实现%平滑是为了消除偶然因素对地理变量的影响%以使地理变量随时间发展变化的趋势和方向明显化%常用的有三种：移动平均、滑动平均、指数平滑%数据格式：单行(按时间序列排序)%数据文件名：data_ph.mat,变量名：testloaddata_ph.matM=length(test);k=input(请输入移动点数（时距）：);y=zeros(1,M);ifkMfori=1:M-kforj

15、=1:ky(i+k)=y(i+k)+test(i+k-j);endy(i+k)=y(i+k)/k;endy(1:k)=NaN;yy=y(M)+(y(M)-y(M-k)/k;str=int2str(k),点移动平均结果如下：;disp(str)formatcompactformatshortdata=test,result=yifisnan(y(M-k)stri=预测下一时段数值为：,num2str(yy);disp(stri)elsedisp(移动点数太大，无法预测下一时段数值！)endformatelsedisp(Error:数据个数不足！)end附（直接复制到.mat文件中即可）：（某城市

16、1999-2004年用水量数据）211.3260.18209.1248.79241250K阶自相关系数的matlab实现（含置信度检验）%k阶自相关系数的计算%数据格式：单行%数据文件名：data_zxgk.mat，变量名：testclear,clcloaddata_zxgk.matM=length(test);k=input(请输入自相关阶数K：);ifkra则通过检验)elsedisp(Error:K值太大！)end附（直接复制到.mat文件中即可）：3149.443303.663010.33109.613639.213253.83466.53839.93894.664009.614253

17、.254101.54119.884258.654401.79虽然matlab里面有这些函数，但是要求自己编写，计算机视觉上有这个实验，是别人编写的。别人到网上找了半天才零散的找到一些碎片，整理以后发上来的！MatLab自编的均值滤波、中值滤波、高斯滤波图像处理函数。%自编的均值滤波函数。x是需要滤波的图像,n是模板大小(即nn)functiond=avefilt(x,n)a(1:n,1:n)=1;%a即nn模板,元素全是1p=size(x);%输入图像是pq的,且pn,qnx1=double(x);x2=x1;%A(a:b,c:d)表示A矩阵的第a到b行,第c到d列的所有元素fori=1:p(

18、1)-n+1forj=1:p(2)-n+1c=x1(i:i+(n-1),j:j+(n-1).*a;%取出x1中从(i,j)开始的n行n列元素与模板相乘s=sum(sum(c);%求c矩阵(即模板)中各元素之和x2(i+(n-1)/2,j+(n-1)/2)=s/(n*n);%将模板各元素的均值赋给模板中心位置的元素endend%未被赋值的元素取原值d=uint8(x2);%自编的中值滤波函数。x是需要滤波的图像,n是模板大小(即nn)functiond=midfilt(x,n)p=size(x);%输入图像是pq的,且pn,qnx1=double(x);x2=x1;fori=1:p(1)-n+1

19、forj=1:p(2)-n+1c=x1(i:i+(n-1),j:j+(n-1);%取出x1中从(i,j)开始的n行n列元素,即模板(nn的)e=c(1,:);%是c矩阵的第一行foru=2:ne=e,c(u,:);%将c矩阵变为一个行矩阵endmm=median(e);%mm是中值x2(i+(n-1)/2,j+(n-1)/2)=mm;%将模板各元素的中值赋给模板中心位置的元素endend%未被赋值的元素取原值d=uint8(x2);%自编的高斯滤波函数，S是需要滤波的图象，n是均值，k是方差functiond=gaussfilt(k,n,s)Img=double(s);n1=floor(n+1

20、)/2);%计算图象中心fori=1:nforj=1:nb(i,j)=exp(-(i-n1)2+(j-n1)2)/(4*k)/(4*pi*k);endend%生成高斯序列b。Img1=conv2(Img,b,same);%用生成的高斯序列卷积运算，进行高斯滤波d=uint8(Img1);%此为程序主文件，包含主要功能单元，以及对子函数进行调用try%实验步骤一：彩色、灰度变换h=imread(photo.jpg);%读入彩色图片c=rgb2gray(h);%把彩色图片转化成灰度图片，256级figure,imshow(c),title(原始图象);%显示原始图象g=imnoise(c,gaussian,0.1,0.002);%加入高斯噪声figure,imshow(g),title(加入高斯噪声之后的图象);%显示加入高斯噪声之后的图象%实验步骤二：用系统预定义滤波器进行均值滤波n=input(请输入均值滤波器模板大小n);A=fspecial(average,n);%生成系统预定义的3X3滤波器Y=filter2(A,g)/255;%用生成的滤波器进行滤波,并归一化figure,imshow(Y),title(用系统函数进行均值滤波后的结果);%显示滤波后的图象%实验步骤三:用自己的编写的函数