多元回归分析matlab剖析

回归分析MATLAB工具箱

一、多元线性回归

多元线性回归：y??0??1x1?...??pxp 1、确定回归系数的点估计值： 命令为：b=regress(Y, X ) ????0????1?①b表示b???...?

??????p???Y1??Y?2②Y表示Y???

?...????Yn??1x11?1x21③X表示X???......???1xn1x12x22...xn2...x1p?...x2p?? ......??...xnp??2、求回归系数的点估计和区间估计、并检验回归模型：

命令为：[b, bint,r,rint,stats]=regress(Y,X,alpha) ①bint表示回归系数的区间估计. ②r表示残差.

③rint表示置信区间.

④stats表示用于检验回归模型的统计量,有三个数值：相关系数r2、F值、与F对应的概率p.

说明：相关系数r越接近1,说明回归方程越显著；F?F1??(k,n?k?1)时拒绝H0,F越大,说明回归方程越显著；与F对应的概率p??时拒绝H0,回归模型成立. ⑤alpha表示显著性水平(缺省时为0.05) 3、画出残差及其置信区间. 命令为：rcoplot(r,rint) 例1.如下程序. 解：(1)输入数据.

x=[143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164]'; X=[ones(16,1) x];

Y=[88 85 88 91 92 93 93 95 96 98 97 96 98 99 100 102]'; (2)回归分析及检验.

[b,bint,r,rint,stats]=regress(Y,X) b,bint,stats

得结果：b = bint =

2 -16.0730 -33.7071 1.5612 0.7194 0.6047 0.8340 stats =

0.9282 180.9531 0.0000

???16.073??0.7194；?的置信区间为[-33.7017,1.5612], ??的置信区间为[0.6047,0.834]; 即?,??0101r2=0.9282, F=180.9531, p=0.0000,我们知道p<0.05就符合条件, 可知回归模型 y=-16.073+0.7194x成立. (3)残差分析,作残差图. rcoplot(r,rint)

Residual Case Order Plot4321Residuals0-1-2-3-4-5246810Case Number121416 从残差图可以看出,除第二个数据外,其余数据的残差离零点均较近,且残差的置信区间均包含零点,这说明回归模型 y=-16.073+0.7194x能较好的符合原始数据,而第二个数据可视为异常点.

(4)预测及作图.

z=b(1)+b(2)*x plot(x,Y,'k+',x,z,'r')

二、多项式回归 (一)一元多项式回归.

1、一元多项式回归：y?a1xm?a2xm?1?...?amx?am?1

(1)确定多项式系数的命令：[p,S]=polyfit(x,y,m)

说明：x=(x1,x2,…,xn),y=(y1,y2,…,yn)；p=(a1,a2,…,am+1)是多项式y=a1xm+a2xm-1+…+amx+am+1的系数；S是一个矩阵,用来估计预测误差. (2)一元多项式回归命令：polytool(x,y,m) 2、预测和预测误差估计.

(1)Y=polyval(p,x)求polyfit所得的回归多项式在x处的预测值Y；

(2)[Y,DELTA]=polyconf(p,x,S,alpha)求polyfit所得的回归多项式在x处的预测值Y及预测值的显著性为1-alpha的置信区间Y±DELTA；alpha缺省时为0.5.

例1. 观测物体降落的距离s与时间t的关系,得到数据如下表,求s. (关于t的回归方程

??a?bt?ct2) st (s) s (cm) t (s) s (cm) 1/30 11.86 8/30 61.49 2/30 15.67 9/30 72.90 3/30 20.60 4/30 26.69 11/30 99.08 5/30 33.71 12/30 113.77 6/30 41.93 13/30 129.54 7/30 51.13 14/30 146.48 10/30 85.44 解法一：直接作二次多项式回归. t=1/30:1/30:14/30;

s=[11.86 15.67 20.60 26.69 33.71 41.93 51.13 61.49 72.90 85.44 99.08 113.77 129.54 146.48]; [p,S]=polyfit(t,s,2) 得回归模型为：

??489.2946st2?65.8896t?9.1329

解法二：化为多元线性回归. t=1/30:1/30:14/30;

s=[11.86 15.67 20.60 26.69 33.71 41.93 51.13 61.49 72.90 85.44 99.08 113.77 129.54 146.48]; T=[ones(14,1) t' (t.^2)'];

[b,bint,r,rint,stats]=regress(s',T); b,stats

得回归模型为：

??9.1329?65.8896st?489.2946t2

预测及作图： Y=polyconf(p,t,S) plot(t,s,'k+',t,Y,'r')

(二)多元二项式回归

多元二项式回归命令：rstool(x,y,’model’, alpha)

说明：x表示n?m矩阵；Y表示n维列向量；alpha：显著性水平(缺省时为0.05)；model表示由下列4个模型中选择1个(用字符串输入,缺省时为线性模型)：

linear(线性)：y??0??1x1????mxm

purequadratic(纯二次)：y??0??1x1????mxm?interaction(交叉)：y??0??1x1????mxm?n??j?1jjx2j

1?j?k?m??jkxjxk

jkquadratic(完全二次)：y??0??1x1????mxm?1?j,k?m??90

xjxk

例1. 设某商品的需求量与消费者的平均收入、商品价格的统计数据如下,建立回归模型,预测

平均收入为1000、价格为6时的商品需求量. 需求量 100 收入 价格 75

80 1200 6 70 500 6 50 300 8 65 400 7 100 1100 4 110 1300 3 60 300 9 1000 600 5 7 1300

5 22解法一：选择纯二次模型,即y??0??1x1??2x2??11x1. ??22x2直接用多元二项式回归：

x1=[1000 600 1200 500 300 400 1300 1100 1300 300]; x2=[5 7 6 6 8 7 5 4 3 9];

y=[100 75 80 70 50 65 90 100 110 60]'; x=[x1' x2'];

rstool(x,y,'purequadratic')

在左边图形下方的方框中输入1000,右边图形下方的方框中输入6，则画面左边的“Predicted Y”下方的数据变为88.47981,即预测出平均收入为1000、价格为6时的商品需求量为88.4791. 在画面左下方的下拉式菜单中选”all”, 则beta、rmse和residuals都传送到Matlab工作区中. 在Matlab工作区中输入命令：beta, rmse 得结果：beta =

110.5313 0.1464 -26.5709 -0.0001 1.8475 rmse =

4.5362

22故回归模型为：y?110.5313?0.1464x1?26.5709x2?0.0001 x1?1.8475x2剩余标准差为4.5362, 说明此回归模型的显著性较好.

22解法二：将y??0??1x1??2x2??11x1化为多元线性回归： ??22x2X=[ones(10,1) x1' x2' (x1.^2)' (x2.^2)']; [b,bint,r,rint,stats]=regress(y,X); b,stats

结果为: b =

110.5313 0.1464 -26.5709 -0.0001 1.8475 stats =

0.9702 40.6656 0.0005

三、非线性回归

1、非线性回归：

(1)确定回归系数的命令：[beta,r,J]=nlinfit(x,y,’model’, beta0)

说明：beta表示估计出的回归系数；r表示残差；J表示Jacobian矩阵；x,y表示输入数据x、y分别为矩阵和n维列向量,对一元非线性回归,x为n维列向量；model表示是事先用m-文件定义的非线性函数；beta0表示回归系数的初值. (2)非线性回归命令：nlintool(x,y,’model’, beta0,alpha) 2、预测和预测误差估计：

[Y,DELTA]=nlpredci(’model’, x,beta,r,J)

表示nlinfit 或nlintool所得的回归函数在x处的预测值Y及预测值的显著性为1-alpha的置信区间Y±DELTA. 例1. 如下程序.

解：(1)对将要拟合的非线性模型y=aeb/x,建立m-文件volum.m如下：

function yhat=volum(beta,x) yhat=beta(1)*exp(beta(2)./x); (2)输入数据： x=2:16;

y=[6.42 8.20 9.58 9.5 9.7 10 9.93 9.99 10.49 10.59 10.60 10.80 10.60 10.90 10.76]; beta0=[8 2]'; (3)求回归系数：

[beta,r ,J]=nlinfit(x',y','volum',beta0)； beta (4)运行结果：

beta =

11.6036 -1.0641 即得回归模型为：

y?11.6036e?1.10641x

(5)预测及作图：

[YY,delta]=nlpredci('volum',x',beta,r ,J)； plot(x,y,'k+',x,YY,'r')

四、逐步回归

1、逐步回归的命令：stepwise(x,y,inmodel,alpha)

说明：x表示自变量数据,n?m阶矩阵；y表示因变量数据,n?1阶矩阵；inmodel表示矩阵的列数的指标,给出初始模型中包括的子集(缺省时设定为全部自变量)；alpha表示显著性水平(缺省时为0.5).

2、运行stepwise命令时产生三个图形窗口：Stepwise Plot,Stepwise Table,Stepwise History. 在Stepwise Plot窗口,显示出各项的回归系数及其置信区间.

(1)Stepwise Table窗口中列出了一个统计表,包括回归系数及其置信区间,以及模型的统计量剩余标准差(RMSE)、相关系数(R-square)、F值、与F对应的概率P.

例1. 水泥凝固时放出的热量y与水泥中4种化学成分x1、x2、x3、 x4有关,今测得一组数据如下,试用逐步回归法确定一个线性模型. 1 2 3 4 5 6 7 8 9 10 11 12 13