。
当然,也可以利用直线拟合得到同一方程。
只不过不能得到参数置信区间和对模型进行检验。
拟合程序如下:
y=[143145146147149150153154155156157158159160162164];
x=[8885889192939395969897969899100102];
a=polyfit(x,y,1)
temp=polyval(a,x);
plot(x,y,'r*',x,temp)
注意:
函数相同,但输出一次函数参数顺序与回归分析(升幂排列)中不同。
另一个差别是拟合不能发现奇异数据。
8.2多元线性回归分析
8.2.1多元线性回归模型的建模步骤及其MATLAB实现
如果根据经验和有关知识认为与因变量有关联的自变量不止一个,那么就应该考虑用最小二乘准则建立多元线性回归模型。
设影响因变量
的主要因素(自变量)有m个,记
,假设它们有如下的线性关系式:
,
如果对变量
与自变量
同时作n次观察(n>m)得n组观察值,采用最小二乘估计求得回归方程
.
建立回归模型是一个相当复杂的过程,概括起来主要有以下几个方面工作
(1)根据研究目的收集数据和预分析;
(2)根据散点图是否具有线性关系建立基本回归模型;(3)模型的精细分析;(4)模型的确认与应用等。
收集数据的一个经验准则是收集的数据量(样本容量)至少应为可能的自变量数目的6~10倍。
在建模过程中首先要根据所研究问题的目的设置因变量,然后再选取与该因变量有统计关系的一些变量作为自变量。
我们当然希望选择与问题关系密切的变量,同时这些变量之间相关性不太强,这可以在得到初步的模型后利用MATLAB软件进行相关性检验。
下面通过一个案例探讨MATLAB软件在回归分析建模各个环节中如何应用。
多元线性回归的MATLAB实现
仍然用命令regress(y,X),只是要注意矩阵X的形式,将通过如下例子说明其用法。
8.2.2某类研究学者的年薪
1.问题
例2工薪阶层关心年薪与哪些因素有关,以此可制定出它们自己的奋斗目标。
某科学基金会希望估计从事某研究的学者的年薪Y与他们的研究成果(论文、着作等)的质量指标X1、从事研究工作的时间X2、能成功获得资助的指标X3之间的关系,为此按一定的实验设计方法调查了24位研究学者,得到如下数据(i为学者序号):
表8-2从事某种研究的学者的相关指标数据
i
1
2
3
4
5
6
7
8
9
10
11
12
3.5
5.3
5.1
5.8
4.2
6.0
6.8
5.5
3.1
7.2
4.5
4.9
9
20
18
33
31
13
25
30
5
47
25
11
6.1
6.4
7.4
6.7
7.5
5.9
6.0
4.0
5.8
8.3
5.0
6.4
33.2
40.3
38.7
46.8
41.4
37.5
39.0
40.7
30.1
52.9
38.2
31.8
i
13
14
15
16
17
18
19
20
21
22
23
24
8.0
6.5
6.6
3.7
6.2
7.0
4.0
4.5
5.9
5.6
4.8
3.9
23
35
39
21
7
40
35
23
33
27
34
15
7.6
7.0
5.0
4.4
5.5
7.0
6.0
3.5
4.9
4.3
8.0
5.8
43.3
44.1
42.5
33.6
34.2
48.0
38.0
35.9
40.4
36.8
45.2
35.1
试建立Y与
之间关系的数学模型,并得出有关结论和作统计分析。
2.作出因变量Y与各自变量的样本散点图
作散点图的目的主要是观察因变量Y与各自变量间是否有比较好的线性关系,以便选择恰当的数学模型形式。
下图分别为年薪Y与成果质量指标
、研究工作时间
、获得资助的指标
之间的散点图,
subplot(1,3,1),plot(x1,Y,'g*'),
subplot(1,3,2),plot(x2,Y,'k+'),
subplot(1,3,3),plot(x3,Y,'ro'),
从图可以看出这些点大致分布在一条直线旁边,因此,有比较好的线性关系,可以采用线性回归。
Y与x1的散点图Y与x2的散点图Y与x3的散点图
图8.1因变量Y与各自变量的样本散点图
3.利用MATLAB统计工具箱得到初步的回归方程
设回归方程为:
.
建立m-文件输入如下程序数据:
x1=[3.55.35.15.84.26.06.85.53.17.24.54.98.06.56.53.76.27.04.04.55.95.64.83.9];
x2=[920183331132530547251123353921740352333273415];
x3=[6.16.47.46.77.55.96.04.05.88.35.06.47.67.05.04.05.57.06.03.54.94.38.05.0];
Y=[33.240.338.746.841.437.539.040.730.152.938.231.843.344.142.533.634.248.038.035.940.436.845.235.1];
n=24;m=3;
X=[ones(n,1),x1',x2',x3'];
[b,bint,r,rint,s]=regress(Y',X,0.05);
b,bint,r,rint,s,
运行后即得到结果如表8-3所示。
表8-3对初步回归模型的计算结果
回归系数
回归系数的估计值
回归系数的置信区间
18.0157
[13.905222.1262]
1.0817
[0.39001.7733]
0.3212
[0.24400.3984]
1.2835
[0.66911.8979]
=0.9106F=67.9195p<0.0001
=3.0719
计算结果包括回归系数b=(
)=(18.0157,1.0817,0.3212,1.2835),且置信区间均不包含零点,;残差及其置信区间;统计变量stats,它包含四个检验统计量:
相关系数的平方
,假设检验统计量F,与F对应的概率p,
的值(7.0以前版本
也可由程序sum(r.^2)/(n-m-1)计算)。
因此我们得到初步的回归方程为:
由结果对模型的判断:
回归系数置信区间不包含零点表示模型较好,残差在零点附近也表示模型较好,接着就是利用检验统计量R,F,p 的值判断该模型是否可用。
(1)相关系数R的评价:
一般地,相关系数绝对值在0.8~1范围内,可判断回归自变量与因变量具有较强的线性相关性。
本例R的绝对值为0.9542,表明线性相关性较强。
(2)F检验法:
当
,即认为因变量
与自变量
之间显着地有线性相关关系;否则认为因变量
与自变量
之间线性相关关系不显着。
本例 F=67.919>
=3.10(查F分布表或输入命令finv(0.95,3,20)计算)。
(3)p值检验:
若
(
为预定显着水平),则说明因变量
与自变量
之间显着地有线性相关关系。
本例输出结果,p<0.0001,显然满足P<
=0.05。
以上三种统计推断方法推断的结果是一致的,说明因变量
与自变量之间显着地有线性相关关系,所得线性回归模型可用。
当然越小越好,这主要在模型改进时作为参考。
4.模型的精细分析和改进
(1)残差分析
残差
,是各观测值
与回归方程所对应得到的拟合值
之差,实际上,它是线性回归模型中误差
的估计值。
即有零均值和常值方差,利用残差的这种特性反过来考察原模型的合理性就是残差分析的基本思想。
利用MATLAB进行残差分析则是通过残差图或时序残差图。
残差图是指以残差为纵坐标,以其他指定的量为横坐标的散点图。
主要包括:
(1)横坐标为观测时间或观测值序号;
(2)横坐标为某个自变量的观测值;(3)横坐标为因变量的拟合值。
通过观察残差图,可以对奇异点进行分析,还可以对误差的等方差性以及对回归函数中是否包含其他自变量、自变量的高次项及交叉项等问题给出直观的检验。
以观测值序号为横坐标,残差为纵坐标所得到的散点图称为时序残差图,画出时序残差图的MATLAB语句为rcoplot(r,rint)(图8.2)。
可以清楚看到残差大都分布在零的附近,因此还是比较好的,不过第4、12、19这三个样本点的残差偏离原点较远,如果作为奇异点看待,去掉后重新拟合,则得回归模型为:
且回归系数的置信区间更小均不包含原点,统计变量stats包含的三个检验统计量:
相关系数的平方
,假设检验统计量F,概率P ,分别为:
0.9533;115.5586;0.0000,比较可知R,F均增加模型得到改进。
图8.2时序残差图
(2)变量间的交互作用讨论
变量间的交互作用包括:
不同自变量之间的交互作用以及同一变量的自相关性。
不同自变量之间的交互作用:
有时,在实验中不仅单因素对指标有影响,而且因素间还会联合起来对指标产生影响,常称这种联合作用为交互作用。
处理两个因素间交互作用的一个简单办法是加入这两个自变量的乘积项。
本文案例如果加入交互项则为:
用表8.2的数据,利用MATLAB统计工具箱得到回归系数分别为:
27.0727,1.1147,-0.0215,-0.1843,0.0033,-0.0054,0.0511。
但它们的置信区间均包含原点,其他指标也不理想,因此,本例中其交互作用并不显着,该模型不如前面两个模型好。
自相关性的诊断和处理:
若数据是以时间为序的,称为时间序列数据。
在时间序列数据中,同一变量的顺序观测值之间出现的相关现象称为自相关。
一旦数据中存在这种自相关序列,如果仍采用普通的回归模型直接处理,将产生不良后果,使预测失去意义。
自相关的诊断主要有图示检验法、相关系数法和DW检验法。
图示检验法是通过绘制残差
散点图观察,如果散布点
大部分点落在第Ⅰ,Ⅲ象限,表明存在着正的序列相关;如果大部分点落在第Ⅱ,Ⅳ象限,表明存在着负的序列相关。
对DW检验法可以利用MATLAB软件编程计算统计量:
,
然后查阅DW检验上下界表,以决定模型的自相关状态。
当一个回归模型存在序列相关性时,首先要查明序列相关产生的原因。
如果是回归模型选用不当,则应改用适当的回归模型;如果是缺少重要的自变量,则应增加自变量;如果以上方法都不能消除序列相关性,则需要采用差分法、迭代法等处理,更详细内容参见相关概率统计参考文献。
8.2.3逐步回归方法建模
逐步回归就是一种从众多自变量中有效地选择重要变量的方法。
逐步回归的基本思路是,先确定一个包含若干自变量的初始集合,然后每次从集合外的变量中引入一个对因变量影响最大的,再对集合中的变量进行检验,从变得不显着的变量中移出一个影响最小的,依此进行,直到不能引入和移出为止。
引入和移出都以给定的显着性水平为标准。
MATLAB统计工具箱中逐步回归的命令是stepwise,它提供了一个人机交互式画面,通过此工具可以自由地选择变量进行统计分析。
该命令的用法是:
stepwise(X,Y,inmodel,alpha)
其中X是自变量数据,排成
矩阵(m为自变量个数,n为每个变量的数据量),Y是因变量数据,排成
向量,inmodel是自变量初始集合的指标,缺省时为全部自变量,alpha为显着水平,缺省时为0.05。
运行stepwise命令时产生图形窗口:
StepwisePlot,StepwiseTable,StepwiseHistory.当鼠标移到图形某个区域时,鼠标点击后产生交互作用。
StepwisePlot窗口中的虚线表示回归系数的置信区间包含零点,即该回归系数与零无显着差异,一般应将该变量移去;实线则表明该回归系数与零有显着差异,应保留在模型中(蓝色表示该变量已进入模型,红色表示该变量已移出模型)。
引入和移出变量还可参考StepwiseHistory窗口中剩余标准差RMSE是否在下降,剩余标准差RMSE最小的就是最好的模型。
StepwiseTable窗口中列出了一个统计表,包括回归系数及其置信区间,以及模型的统计量剩余标准差RMSE、相关系数R-square、F值、与F对应的概率。
关于本节案例2,如果引入新的自变量
.也可以采用逐步回归法解决,源程序如下:
A=[3.55.35.15.84.26.06.85.53.17.24.54.98.06.56.53.76.27.04.04.55.95.64.83.9;920183331132530547251123353921740352333273415;6.16.47.46.77.55.96.04.05.88.35.06.47.67.05.04.05.57.06.03.54.94.38.05.0]';
Y=[33.240.338.746.841.437.539.040.730.152.938.231.843.344.142.533.634.248.038.035.940.436.845.235.1]';
x1=A(:
1);
x2=A(:
2);
x3=A(:
3);
x4=x1.*x2;
x5=x1.*x3;
x6=x2.*x3;
X=[A,x4,x5,x6];
stepwise(X,Y)
运行并按上述步骤操作后可以得到本文前面线性回归相同的结论,即不含交互项的模型是最好的。
在此只介绍操作过程,其交互界面,只要在MATLAB软件上一试便知。
8.2.4多项式回归
多项式回归仍然属于多元线性回归,可以是一元多项式回归或多元多项式回归。
一元多项式回归模型的一般形式为
用MATLAB求解一元多项式回归,除了使用命令polyfit(x,y,m)外,还可以使用如下命令:
Polytool(x,y,m,alpha)
输入x,y,m同命令polyfit,alpha是显着性水平(默认0.05),则输出一个交互式画面,画面显示回归曲线及其置信区间,通过图左下方的export下拉式菜单,还可以输出回归系数估计值及其置信区间、残差等。
下面通过一个用多元多项式回归的实例说明什么时候用多项式回归以及如何通过MATLAB软件进行处理。
例3为了了解人口平均预期寿命与人均国内生产总值和体质得分的关系,我们查阅了国家统计局资料,北京体育大学出版社出版的《2000国民体质监测报告》,表8-4是我国大陆31个省市的有关数据。
我们希望通过这几组数据考察它们是否具有良好的相关关系,并通过它们的关系从人均国内生产总值(可以看作反映生活水平的一个指标)、体质得分预测其寿命可能的变化范围。
体质是指人体的质量,是遗传性和获得性的基础上表现出来的人体形态结构,生理机能和心理因素综合的、相对稳定的特征。
体质是人的生命活动和工作能力的物质基础。
它在形成、发展和消亡过程中,具有明显的个体差异和阶段性。
中国体育科学学会体质研究会研究表明,体质应包括身体形态发育水平、生理功能水平、身体素质和运动能力发展水平、心理发育水平和适应能力等五个方面。
目前,体质的综合评价主要是形态、机能和身体素质三类指标按一定的权重进行换算而得。
表8-431个省市人口预期寿命与人均国内生产总值和体质得分数据
序号
预期寿命
体质得分
人均产值
序号
预期寿命
体质得分
人均产值
序号
预期寿命
体质得分
人均产值
1
71.54
66.165
12857
12
65.49
56.775
8744
2