[工学]Matlab与统计分析.doc

Matlab与统计分析一、 回归分析1、多元线性回归1.1 命令 regress( )， 实现多元线性回归，调用格式为b,bint,r,rint,stats=regress(y,x,alpha)其中因变量数据向量Y和自变量数据矩阵x按以下排列方式输人对一元线性回归，取k=1即可。alpha为显著性水平(缺省时设定为0.05)，输出向量b，bint为回归系数估计值和它们的置信区间，r，rint为残差及其置信区间，stats是用于检验回归模型的统计量，有三个数值，第一个是， 其中R是相关系数，第二个是F统计量值，第三个是与统计量F对应的概率P，当 时拒绝，回归模型成立.注：1、两组数据的相关系数在概率论的标准定义是：R= E(x - Ex) * (y - Ey) / (sqrt(x - Ex)2) * sqrt(y - Ey)2)E求取期望值。也就是两组数据协方差与两者标准差乘积的商。如果|R|=1说明两者相关，R=0说明两者不相关.1、 F是方差分析中的一个指标，一般方差分析是比较组间差异的。F值越大，P值越小，表示结果越可靠.1.2 命令 rcoplot(r，rint)，画出残差及其置信区间.1.3 实例 1已知某胡八年来湖水中COD浓度实测值(v)与影响因素湖区工业产值(x1)、总人口数(x2 )、捕鱼量(x3 )、降水量( x4)资料，建立污染物Y的水质分析模型.Step 1 输入数据x1=1.376, 1.375, 1.387, 1.401, 1.412, 1.428, 1.445, 1.477;x2=0.450,0.475,0.485,0.500,0.535,0.545,0.550,0.575;x3=2.170,2.554,2.676,2.713,2.823,3.088,3.122,3.262;x4=0.8922, 1.1610，0.5346，0.9589, 1.0239, 10499,1.1065, 1.1387;Y=5.19, 5.30，5.60，5.82，6.00，6.06，6.45，6.95;Step 2 保存数据（以数据文件.mat 形式保存，便于以后调用）save data x1 x2 x3 x4 yload data %取出数据Step 3 执行回归命令x=ones(8,1),x1,x2,x3,x4；b,bint,r,rint,stats=regress(y,x)得到结果：b=(-16.5283, 15.7206， 2.0327-0.2106，-0.1991)stats=(0.9908,80.9530，0.0022)即Y= -16.5283+15.7206x1+2.0327x2-0.2106xl+0.1991x4=0.9908, F=80.9530,P=0.00222、 非线性回归2.1 命令 nlinfit( ) 实现非线性回归，调用格式为beta,r,J=nlinfit(x,y,model,beta0)其中，输入数据x，y分别为nm矩阵和n维列向量，对一元非线性回归，x为n维列向量；model是事先用m-文件定义的非线性函数，beta0是回归系数的初值beta是估计出的回归系数，r是残差，J是Jacobian矩阵，它们是估计预测误差需要的数据.2.2 命令 nlpredci( ) 预测和预测误差的估计，调用格式为y,delta=npredci(model,x,beta,r,j)2.3 实例 2 对实例1中COD浓度实测值(Y)，建立时序预测模型，这里选用logistic模型，即Step 1 建立非线性函数 对所要拟合的非线性模型建立m-文件model.m如下function yhat=model(beta,t)yhat=beta(1).(1+beta(2)*exp(-beta(3)*t)Step 2 输入数据t= 1:8load data y(在datamat中取出数据y)beta0=50，10，1Step 3 求回归系数beta，r，J=nlinfit(t，Y，model, beta0)得结果：beta=(56.1157，10.4006，0.0445)即Step 4 预测及作图YY,delta=nlpredci(model,x,beta,r ,J)；plot(x,y,k+,x,YY,r)3、 逐步回归逐步回归的命令是stepwise， 它提供了一个交互式画面通过此工具可自由地选择变量，进行统计分析.调用格式为：stepwise（x，y，inmodel，alpha）其中x是自变量数据,是阶矩阵，y是因变量数据，阶矩阵，inmodel是矩阵的列数指标（给出初始模型中包括的子集（缺省时设定为全部自变量），alpha是显著性水平（缺省时为0.5）. 运行stepwise命令时产生三个图形窗口：Stepwise Plot，Stepwise Table，Stepwise History.在Stepwise Plot窗口，显示出各项的回归系数及其置信区间.Stepwise Table 窗口中列出了一个统计表，包括回归系数及其置信区间，以及模型的统计量剩余标准差（RMSE）、相关系数（R-square）、F值、与F对应的概率P.2、 主成分分析这里给出江苏省生态城市主成份分析实例。我们对江苏省十个城市的生态环境状况进行了调查，得到生态环境指标的指数值，见表1。现对生态环境水平进行分析和评价。我们利用Matlab65中的princomp命令实现。具体程序如下x=07883 07391 08111 06587 06543 08259 08486 06834 08495 0784607633 07287 07629 08552 07564 07455 07800 09490 08918 0895404745 05126 08810 08903 08288 07850 08032 08862 03987 0397008246 07603 06888 08977 07926 07856 06509 08902 06799 0987708791 08736 08183 09446 09202 09263 09185 09505 08620 0887309538 09257 09285 09434 09154 08871 09357 08760 09579 0974108785 08542 08537 09027 08729 08485 08473 09044 08866 0903506305 06187 06313 07415 06398 06142 05734 08980 06186 0738208928 07831 05608 08419 08464 07616 08234 06384 09604 08514x=x;stdr=std(x)； 求各变量标准差n,m=size(x)；sddata= x.stdr(ones(n,1),:)； 标准化变换p，princ，egenvalue=princomp(sddata) 调用主成分分析程序p3=p(：，1:3) 输出前三个主成分系数sc=princ(：，1:3) 输出前三个主成分得分egenvalue 输出特征根per=100*egenvaluesum(egenvalue) 输出各个主成分贡献率执行后得到所要结果，这里是前三个主成分、主成分得分、特征根。即egenvalue=38811，26407，10597 ， per=4312，2934，11971 这样，前三个主成分为Zl = -03677xl+ 03702x2+ 01364x3+ 04048x4+ 03355x5-01318x6+04236x7+ 04815x8-00643x9Z2 = 01442xl+ 02313x2-05299x3+ 01812x4-01601x5+ 05273x6+03116x7-00267x8+ 04589x9Z3 = -03282xl-03535x2+ 00498x3+ 00582x4+ 05664x5-00270x6-00958x7-02804x8+ 05933x9第一主成分贡献率为43.12，第二主成分贡献率为29.34，第三主成分贡献率为11.97，前三个主成分累计贡献率达8424。如果按80 以上的信息量选取新因子，则可以选取前三个新因子。第一新因子Z1包含的信息量最大为43.12，它的主要代表变量为x8(城市文明)、x7(生产效率)、 x4(城市绿化)，其权重系数分别为04815、04236、04048，反映了这三个变量与生态环境水平密切相关，第二新因子Z2包含的信息量次之为2934，它的主要代表变量为x3(地理结构)、x6(资源配置)、 x9(可持续性)，其权重系数分别为05299、05273、04589，第三新因子 Z3包含的信息量为1197，代表总量为 x9(可持续性)、 x5(物质还原)，权重系数分别为05933、05664。这些代表变量反映了各自对该新因子作用的大小，它们是生态环境系统中最重要的影响因素。根据前三个主成分得分，用其贡献率加权，即得十个城市各自的总得分F = 4312princ(：，1)+2934princ(：，2)+1197princ(：，3) =00970，-06069，-15170，11801，00640，-08178，-09562，11383，01107，13077根据总得分排序，结果见表1。3、 聚类分析我们对苏州所辖张家港市2003年七条河流中主要污染因子(指标)，即CODmn，BOD5，非离子氨，氨氮，挥发酚，石油类共6个变量(资料见表2，来源于张家港市2003年环境质量报告书)，进行聚类分析。我们利用Matlab65中的cluster命令实现。具体程序如下x= 3.14 8.41 23.78 25.79 4.17 6.47 5.47 9.57 26.48 23.79 6.42 6.58 3.1 4.31 21.2 22.48 5.34 6.54 5.67 9.54 10.23 20.87 4.2 6.8 6.81 9.05 16.18 24.56 5.2 5.45 6.21 7.08 21.05 31.56 6.15 8.214.87 8.97 26.54 34.56 5.58 8.07;n,m=size(x);stdr=sta(x);xx=x .stdr(ones(n,1),:); 标准化变换y=pdist(xx); 计算各样本间距离(这里为欧氏距离)z=linkage(y); 进行聚类(这里为最短距离法)h=dendrogram(z); 画聚类谱系图t=cluster(z,3) 将全部样本分为3类find(t=2); 找出属于第2类的样品编号执行后得到所要结果，聚类谱系图见图1