统计回归模型举例

1、多元线性回归命令：[b,bint,r,rint,stats]=regress(y,x,alpha)2、一元多项式回归命令：[p,s]=polyfit(x,y,m)3、多元二项式回归命令：rstool(x,y,’model’,alpha)线性（linear），完全二次（quadratic),纯二次（purequadratic),交叉（interaction)4、非线性回归命令：[beta,r,j]=nlinfit(x,y,’model’,beta0)几个常见回归命令例1牙膏的销售量

问题建立牙膏销售量与价格、广告投入之间的模型

预测在不同价格和广告费用下的牙膏销售量

收集了30个销售周期本公司牙膏销售量、价格、广告费用，及同期其它厂家同类牙膏的平均售价

9.260.556.804.253.70307.930.055.803.853.8029

8.510.256.754.003.7527.38-0.055.503.803.851销售量(百万支)价格差（元）广告费用(百万元)其它厂家价格(元)本公司价格(元)销售周期

令y表示公司牙膏的销售量，x1表示其它厂家与本公司价格差，

x2表示公司广告费用，则数据如下：>>x1=[-0.050.250.600.250.20.150.05-0.150.150.20.10.40.450.350.30.50.50.4-0.05-0.05-0.10.20.10.50.6-0.0500.050.55];x2=[5.56.757.255.576.56.755.255.2566.56.2576.96.86.87.176.86.56.2566.576.86.86.55.755.86.8];>>y=[7.388.519.527.59.338.288.757.877.187.898.159.18.868.98.879.2698.757.957.657.2788.58.759.218.277.677.939.26];

下面探讨y与x1、x2的关系：用matlab软件作图：plot(x1,y,’*’);plot(x2,y,’*’)运行得如下图形：x1y从右图看出，y与x1成线性关系，y与x2成二次曲线关系。x2y

>>x3=x2.^2;>>x=[ones(30,1)x1'x2'x3'];>>[b,bint,r,rint,stats]=regress(y',x)运行结果：b=17.3244，1.3070，-3.6956，0.3486bint=5.728228.92060.68291.9311-7.49890.10770.03790.6594stats=0.9054，82.9409，0.0000，0.0490模型求解MATLAB统计工具箱

结果分析y的90.54%可由模型确定参数参数估计值置信区间17.3244[5.728228.9206]1.3070[0.68291.9311]-3.6956[-7.49890.1077]0.3486[0.03790.6594]R2=0.9054F=82.9409p=0.0000

0

1

2

3F远超过F检验的临界值P<<

=0.05

2的置信区间包含零点(右端点距零点很近)x2对因变量y的影响不太显著由于x22项显著可将x2保留在模型中模型从整体上看成立销售量预测价格差x1=其它厂家价格x3-本公司价格x4估计x3调整x4控制价格差x1=0.2元，投入广告费x2=650万元控制x1通过x1,x2预测yx1=0.2;x2=6.5;Y=b(1)+b(2)*x1+b(3)*x2+b(4)*(x2.^2)运行结果：Y=8.2933即预测牙膏销售量为8.2933百万支。

上述模型中的回归变量x1,x2对因变量y的影响是相互独立的。即牙膏销售量y的均值与广告费x2的二次关系由回归系数β2和β3确定，而不必依赖于差价x1,同样y的均值与x1的线性关系仅由回归系数β1确定，不依赖于x2.根据直觉和经验可以猜想，x1和x2之间的交互作用也会对y有影响，不妨简单地用x1,x2的乘积来表示他们的相互作用，于是上述模型中增加一项，得到：模型改进

>>x=[ones(30,1)x1',x2'(x2.^2)'(x1.*x2)'];>>[b,bint,r,rint,stats]=regress(y',x)b=29.113311.1342-7.60800.6712-1.4777bint=3.701344.52521.977820.2906-12.6932-2.52280.25381.0887-2.8518-0.1037stats=0.9209，72.7771，0.0000，0.0426模型比较x1和x2对y的影响独立

参数参数估计值置信区间17.3244[5.728228.9206]1.3070[0.68291.9311]-3.6956[-7.49890.1077]0.3486[0.03790.6594]R2=0.9054F=82.9409p=0.0000

0

1

2

3参数参数估计值置信区间29.1133[13.701344.5252]11.1342[1.977820.2906]-7.6080[-12.6932-2.5228]0.6712[0.25381.0887]-1.4777[-2.8518-0.1037]R2=0.9209F=72.7771p=0.0000

3

0

1

2

4x1和x2对y的影响有交互作用

由于R2有所提高，所以模型(**)比模型（*）有所改进，并且参数的置信区间不再包含0点，所以有理由认为模型（**）比模型（*）更符合实际。预测比较：x1=0.2;x2=6.5;Y=b(1)+b(2)*x1+b(3)*x2+b(4)*(x2.^2)+b(5)*(x1.*x2)Y=8.3272两模型销售量预测比较(百万支)区间[7.8230，8.7636]区间[7.8953，8.7592](百万支)控制价格差x1=0.2元，投入广告费x2=6.5百万元预测区间长度更短略有增加完全二次多项式模型>>

x=[x1'x2'];>>rstool(x,y','quadratic')运行结果：beta=2.098414.7436-8.6367-2.10381.10740.7594rmse=0.2083剩余标准差为0.2.83较小，说明回归模型的显著性比较好。问题：一家高技术公司人事部门为研究软件开发人员的薪金与他们的资历、管理责任、教育程度等因素之间的关系，要建立一个数学模型，以便分析公司人事策略的合理性，并作为新聘人员的薪金的参考。他们认为目前公司人员的薪金总体上是合理的，可以作为建模的依据。于是调查了46名软件开发人员的档案资料，如下表，其中资历一列指从事专业工作的年数，管理一列中：1表示管理人员，0表示非管理人员，教育一列中：1表示中学程度，2表示大学程度，3表示更高程度（研究生）。例2软件开发人员的薪金

编号薪金资历管理教育编号薪金资历管理教育113876111131980031321160810314114174013187011131520263413411283102161323140351176710317128844026208722121813245502711772202191367750381053520120159655119121952032112366601101231330222213526131114975311231383960212213713122422884612编号薪金资历管理教育编号薪金资历管理教育251697871136168821202261480380237241701213271740481138159901301282218481339263301312291354880140179491402301446710014125685151331159421002422783716123223174101343188381602332378010124417483160134254101112451920717023514861110146193462001

分析与假设——按照常识，薪金自然随着资历（年）的增长而增加，管理人员的薪金应高于非管理人员，教育程度越高薪金也越高。令y表示薪金，x1表示资历，x2表示是否管理人员，x3表示学历

基本模型——假设薪金y与资历x1、管理x2、学历x3成线性关系:

y=[13876116081870111283117672087211772105351219512313149752137119800114172026313231128841324513677159651236621352138392288416978148031740422184135481446715942231742378025410148611688224170159902633017949256852783718838174831920719346];x1=[111112222333344445556666788881010101011111212131314151616161720];

x2=[1010010000111010000101011011000111001010110000];x3=[1332322132123133223113221213112323123122322121];x=[ones(46,1),x1',x2',x3'];>>[b,bint,r,rint,stats]=regress(y',x)b=1.0e+003*6.9333,0.5659,6.5936,1.6134bint=1.0e+003*5.66128.20540.49110.64065.81347.37371.11112.1156stats=0.9327194.016901603719.76601

由于R2=0.9327接近于1，F=194.0169大于临界值，p<<0.05所以模型的显著性较好。回归模型为：残差分析：Rcoplot(r,rint)

模型修正在上述模型中，资历、管理、学历对薪金的影响都是独立的。事实上，管理与学历对薪金应具有交叉影响，为此增加交叉项x2x3,得模型：模型求解——

x=[ones(46,1),x1',x2',x3',(x2.*x3)'];>>[b,bint,r,rint,stats]=regress(y',x)b=8135.915538.3664525.2491077.1391019.748bint=6410.8499860.982461.000615.7322311.7486738.750351.3631802.914-6.3142045.810stats=0.9387157.01201495857.511

R2=0.9387>0.9327,所以，该模型较好。

为了表示三种教育程度，也可引进两个0——1变量来表示：

y=[13876116081870111283117672087211772105351219512313149752137119800114172026313231128841324513677159651236621352138392288416978148031740422184135481446715942231742378025410148611688224170159902633017949256852783718838174831920719346];x1=[111112222333344445556666788881010101011111212131314151616161720];x2=[1010010000111010000101011011000111001010110000];

X3=[1000000100100100000110001010110000100100000101];X4=[0001011001010000110000110100001011010011011010];x=[ones(46,1),x1',x2',x3‘,x4‘];[b,bint,r,rint,stats]=regress(y',x)

b=11032.7343011281546.12764929776882.53291698754-2994.17834433349147.737980069428stats=0.956691811962102226.42579883577701057144.84841479

R2=0.956691811962102F=226.425798835777p<<0.05所以模型的显著性较好。残差分析：rcoplot(r,rint)

模型修正在上述模型中，资历、管理、学历对薪金的影响都是独立的。事实上，管理与学历对薪金应具有交叉影响，为此增加交叉项x2x3,x2x4,得模型模型求解：X=[x(x2.*x3)'(x2.*x4)'];>>[b,bint,r,rint,stats]=regress(y',X)

b=11203.753782227896.8639299123927047.99973466834-1726.5041924628-348.392543178968-3070.596188012791835.9676370463stats=0.998829102890402,5544.79903960134,030047.093445917

R2=0.998829102890402>0.956691811962102

所以，该模型较好。例3投资额与国民生产总值和物价指数

问题建立投资额模型，研究某地区实际投资额与国民生产总值(GNP)及物价指数(PI)的关系2.06883073.0424.5201.00001185.9195.0101.95142954.7474.9190.96011077.6166.491.78422631.7401.9180.9145992.7144.281.63422417.8423.0170.8679944.0149.371.50422163.9386.6160.8254873.4133.361.40051918.3324.1150.7906799.0122.851.32341718.0257.9140.7676756.0125.741.25791549.2206.1130.7436691.1113.531.15081434.2228.7120.7277637.797.421.05751326.4229.8110.7167596.790.91物价指数国民生产总值投资额年份序号物价指数国民生产总值投资额年份序号根据对未来GNP及PI的估计，预测未来投资额该地区连续20年的统计数据

时间序列中同一变量的顺序观测值之间存在自相关以时间为序的数据，称为时间序列

分析许多经济数据在时间上有一定的滞后性

需要诊断并消除数据的自相关性，建立新的模型若采用普通回归模型直接处理，将会出现不良后果

投资额与国民生产总值和物价指数

……………………1.32341718.0257.9140.7676756.0125.741.25791549.2206.1130.7436691.1113.531.15081434.2228.7120.7277637.797.421.05751326.4229.8110.7167596.790.91物价指数国民生产总值投资额年份序号物价指数国民生产总值投资额年份序号

y=[90.997.4113.5125.7122.8133.3149.3144.2166.4195.0229.8228.7206.1257.9324.1386.6423.0401.9474.9424.5];>>x1=[596.7637.7691.1756.0799.0873.4944.0992.71077.61185.91326.41434.21549.21718.01918.32163.92417.82631.72954.73073.0];>>x2=[0.71670.72770.74360.76760.79060.82540.86790.91450.96011.01.05751.15081.25791.32341.40051.50421.63421.78421.95142.0688];t~年份，yt~投资额，x1t~GNP,x2t~物价指数画出散点图投资额与GNP及物价指数间均有很强的线性关系

0,1,2~回归系数x1tytx2tyt

t~对t相互独立的零均值正态随机变量Plot(x1,y,’*’)Plot(x2,y,’*’)

x=[ones(20,1)x1‘x2’]；[b,bint,r,rint,stats]=regress(y',x)b=322.724963028216；0.618456651396034；-859.478998026578bint=224.338557066255；421.1113689901780.477272347881124；0.759640954910943-1121.47567088142；-597.482325171736stats=0.990843999027999919.8528940192380161.707321609199基本回归模型的结果与分析

MATLAB统计工具箱

参数参数估计值置信区间

0322.7250[224.3386421.1114]

10.6185[0.47730.7596]

2-859.4790[-1121.48,-597.48]R2=0.9908F=919.8529p=0.0000剩余标准差

s=12.7164

没有考虑时间序列数据的滞后性影响R2＝0.9908，拟合度高模型优点模型缺点可能忽视了随机误差存在自相关；如果存在自相关性，用此模型会有不良后果

例4、教学评估问题：为了考评教师的教学质量，教学研究部门设计了一个教学评估表，对学生进行一次问卷调查，要求学生对12名教师的15门课程（其中3位教师有2门课）按以下7项内容打分，分值为1—5分（5分最好，1分最差）。X1—课程内容的合理性；x2—主要问题展开的逻辑性；X3—回答学生问题的有效性；x4—课下交流的有助性；X5—教科书的帮助性；x6—考试平分的公平性；y—对教师的总体评价。

收回问卷调查后，得到了学生对12位教师15门课的各项评分的平均值，见下表：教师编号课程编号x1x2x3x4x5x6y12014.464.424.234.14.564.374.1122244.113.823.293.63.993.823.3833013.583.313.243.764.393.753.1743014.424.374.344.43.634.274.3953014.624.474.534.674.634.574.6963093.183.823.923.623.54.143.2573112.472.793.583.52.843.842.8483114.293.924.053.762.764.113.9593124.414.364.274.754.594.114.18103124.594.344.244.392.644.384.44113334.554.454.434.574.454.44.47124244.674.644.524.393.484.214.6133513.713.413.394.184.084.063.1744114.284.454.14.073.764.434.1594244.244.384.354.484.154.54.33

教学研究部门认为，所列各项具体内容x1——x6不一定每项都对教师总体评价y有显著影响，并且各项内容之间也可能存在很强的相关性，他们希望得到一个总体评价与各项具体内容之间的模型，这个模型应尽量简单和有效，并且由由此能给教师一个合理的建议，以提高总体评价。逐步回归的基本思想——先确定一个包含若干自变量的初始集合，然后每次从集合外的变量中引入一个对因变量影响最大的，再对集合中的变量进行检验，从变得不显著的变量中移出一个影响最小的。依次进行，直到不能引入和移出为止。引入和移出都以给定的显著性水平为标准。虽然给出了6个变量，但是我们希望从中挑选出对因变量y有显著影响的哪些来建立回归模型。为此我们采用逐步回归方法。

MATLAB统计工具箱中逐步回归命令为：stepwise通常的用法为：Stepwise(x,y,inmdel,penter,premove)x：自变量数据矩阵；y：因变量数据；Inmodel：自变量初始集合的指标（即矩阵x中哪些列进入初始集合），缺省时设定为没有选取任何x的列向量；Penter：引入变量时设定的最大p值，缺省时为0.05；Premove：移出变量时设定的最小p值，缺省时为0.10。注意：Premove的值不能小于Penter的值。

x1=[4.464.113.584.424.623.182.474.294.414.594.554.673.714.284.24];x2=[4.423.823.314.374.473.822.793.924.364.344.454.643.414.454.38];x3=[4.233.293.244.344.533.923.584.054.274.244.434.523.394.104.35];x4=[4.103.603.764.404.673.623.503.764.754.394.574.394.184.074.48];x5=[4.563.994.393.634.633.502.842.764.592.644.453.484.063.764.15];x6=[4.373.823.754.274.574.143.844.114.114.384.404.214.064.434.50];y=[4.113.383.174.394.693.252.843.954.184.444.474.613.174.154.33];x=[x1'x2'x3'x4'x5'x6'];

模型解释：在最终模型里回归变量只有x1,x2，是一个简单易用的模型。据此可把课程内容组织的合理性（x1）和回答学生问题的有效性(x3)，列入考评的重点，模型（*）表明，x1的分值每增加一分，对教师的总体评价就增加0.5分；x3的每增加1分，对教师的总体评价就增加0.77分，应建议教师注重这两方面的工作。为了分析其他自变量没有最终进入模型的原因，可以计算x1~x6,y的相关系数。

>>

A=[xy'];>>corrcoef(A)ans=1.0000，0.9008，0.6752，0.7361，0.2910，0.6471，0.89730.9008，1.0000，0.8504，0.7399，0.2775，0.8026，0.93630.6752，0.8504，1.0000，0.7499，0.0808，0.8490，0.91160.7361，0.7399，0.7499，1.0000，0.4370，0.7041，0.82190.2910，0.2775，0.0808，0.4370，1.0000，0.1872，0.17830.6471，0.8026，0.8490，0.7041，0.1872，1.0000，0.82460.8973，0.9363，0.9116，0.8219，0.1783，0.8246，1.0000一般认为，两个变量的相关系数超过0.85时才具有显著的相关性。由上面结果知道，与y相关性显著的只有x1,x2,x3，而X2未进入最终模型，是由于它与x1,x3的相关性显著（r12=0.9008，r23=0.8504），可以说，模型中有了x1,x3之后，变量X2是多余的，应该去掉。

例6冠心病与年龄问题：冠心病简称CHD，是一种常见的心脏疾病，严重地危害着人类的健康。到目前为止，其疾病尚未完全研究清楚，医学界普遍认同的、重要的易患因素是高领、高血压、糖尿病、动脉粥样硬化及家族史等。多项研究表明，冠心病发病率随着年龄的增加而上升，在冠心病的流行病学研究中，年龄也最常见的混杂因素之一。为了更好地说明冠心病发病率与年龄的关系，医学界对100名不同年龄的人进行观察，表1给出了这100名被观察者的年龄及是否患冠心病的数据。

表1100名被观察者的年龄与是否患冠心病的观察数据序号年龄冠心病序号年龄冠心病序号年龄冠心病序号年龄冠心病120026350514417655122302735052441775613240283605345078561425029361544517956152513036055460805706260313705646181570726032371574708257182803337058470835719280343805947184571102903538060480855711130036390614818658012300373916248187581133003840063490885811430039401644908959115300404106549190591163014141066500916001732042420675019260118320434206851093611193304442069520946212033045421705219562121340464307153196631223404743072531976402334148431735419864124340494407455099651253405044075551100691

表1冠心病一栏中，1代表患冠心病，0表示不患冠心病。试根据这些数据建立数学模型，来分析冠心病发病率与年龄的关系，并进行统计预测。分析与假设——假设这100名被观察者是独立选取的，记x被观察者的年龄，Y为观察者患冠心病的情况（Y=1表示患冠心病，Y=0表示未患）x=[20232425252626282829303030303030323233333434343434353536363637373738383939404041414242424243434344444444454546464747474848484949495050515252535354555555565656575757575757585858595960606162626364646569;

Y=[0000100000000001000000100000100100001010000010010011010100101100101001111011111001111011110111110111];作出Y对x的散点图Plot(x,Y.’*’)从右图可以看出，直接对上述数据建立回归模型是行不通的，需要对数据进行处理。

数据处理的一种常见方法是将被观察者按年龄进行分组，并统计各年龄段中患冠心病的人数，及患病人数占该组人数的比例（以下简称患病比例）为方便起见，将年龄分成8个年龄段，分组数据如下表：表2各年龄段的冠心病患病人数及比例年龄段组中值人数患病人数患病比例20-2924.51010.130-34321520.1335-39371230.2540-44421550.3345-49471360.4650-5452850.6355-5957171