发布时间 : 星期日 文章回归分析总结更新完毕开始阅读fe6e739b050876323012122b
表8-4 31个省市人口预期寿命与人均国内生产总值和体质得分数据 序号 预期寿命 体质得分 人均产值 序号 预期寿命 体质得分 人均产值 序号 预期寿命 体质得分 人均产值 1 2 3 4 5 6 7 8 9 11
71.54 73.92 73.27 71.20 73.91 72.54 70.66 71.85 71.08 74.70 66.165 71.25 70.135 65.125 69.99 65.765 67.29 67.71 66.525 67.13 69 .505 12857 12 24495 13 24250 14 10060 15 29931 16 18243 17 10763 18 9907 9088 19 21 65.49 56.775 68.95 73.34 65.96 72.37 70.07 66.01 67.97 62.9 66.1 64.51 8744 11494 20461 5382 19070 10935 22007 13594 11474 14335 7898 23 24 25 26 27 28 29 30 31 69.87 64.305 17717 67.41 60.485 15205 78.14 70.29 70622 76.10 69.345 47319 74.91 68.415 40643 72.91 66.495 11781 70.17 65.765 10658 66.03 64.37 63.28 62.84 11587 9725 72.55 68.385 71.65 66.205 65.77 73.10 67.065 67.47 63.605 13255 20 71.73, 33772 22 10 71.29, 模型的建立和求解 作表8-4数据(x1,y),(x2,y)的散点图如图8.3
图8.3 预期寿命与人均国内生产总值和体质得分的散点图
从图8.3可以看出人口预期寿命y与体质得分x2有较好的线性关系，y与人均国内生产总值x1的关系难以确定，我们建立二次函数的回归模型。
一般的多元二项式回归模型可表为 y??0??1x1???mxm?1?j,k?m??jkxjxk??
MATLAB统计工具箱提供了一个很方便的多元二项式回归命令：
Rstool(x,y, 'model',alpha)
输入x为自变量（n×m矩阵），y为因变量（n维向量），alpha为显著水平，model从下列4个模型中选择一个：
linear（只包含线性项）
purequadratic（包含线性项和纯二次项） interaction（包含线性项和纯交互项） quadratic（包含线性项和完全二次项）
输出一个交互式画面，对例3，编程如下：
y=[71.54 73.92 73.27 71.20 73.91 72.54 70.66 71.85 71.08 71.29,74.70 65.49 68.95 73.34 65.96 72.37 70.07 72.55 71.65 71.73,73.10 67.47 69.87 67.41 78.14 76.10 74.91 72.91 70.17 66.03 64.37];
x1=[12857 24495 24250 10060 29931 18243 10763 9907 13255 9088 33772 8744 11494 20461 5382 19070 10935 22007 13594 11474 14335 7898 17717 15205 70622 47319 40643 11781 10658 11587 9725];
x2=[66.165 71.25 70.135 65.125 69.99 65.765 67.29 67.71 66.525 67.13,69.505 56.775 66.01 67.97 62.9 66.1 64.51 68.385 66.205 65.77,67.065 63.605 64.305 60.485 70.29 69.345 68.415 66.495 65.765 63.28 62.84]; x=[x1',x2'];
rstool(x,y','purequadratic')
得到一个如图8.4的交互式画面
图8.4 预期寿命与人均国内生产总值和体质得分的一个交互式画面
左边一幅图形是x2固定时的曲线y(x1)及其置信区间，右边一幅图形是x1固定时的曲线
y(x2)及其置信区间。移动鼠标可改变x1，x2的值，同时图左边给出y的预测值及其置信
区间。如输入x1=128757，x2=66.165，则y=70.6948，其置信区间70.6948±1.1079。 图的左下方有两个下拉式菜单，上面的菜单Export用于输出数据（包括：回归系数parameters,残差residuals,剩余标准差RMSE等）, 在MATLAB工作空间中得到有关数据。通过下面的菜单在上述4个模型中变更选择，最后确定RMSE值较小的模型。例3则是包含线性项和完全二次项（quadratic）的模型最佳，即
2y??0??1x1??2x2??3x1x2??4x12??5x2??
剩余标准差为1.2622，因此，所得回归模型为：
2y?195.36?0.0045x1?5.5753x2?6.7338?10?5x1x2?3.3529?10?9x12?0.055556x2
利用此模型我们可以根据国内生产总值及体质得分，预测寿命
8.3 非线性回归分析
8.3.1 非线性最小二乘拟合
线性最小二乘拟合与线性回归中的“线性”并非指y与x的关系，而是指y是系数
?0,?1或??(?0,?1,,?m)的线性函数。拟合如y??0??1x2的函数仍然是最小二乘拟
?x合；如果拟合如y??0e1的曲线，y对?0,?1是非线性的，但取对数后lny对系数?0,?1是线性的，属于可化为线性回归的类型。下面讨论非线性拟合的情形。
非线性最小二乘拟合问题的提法是：已知模型
y?f(x,?),x?(x1,,xm),??(?0,?1,,?k)，
其中f对?是非线性的，为了估计参数?，收集n个独立观测数据
(xi,yi),xi?(xi1,差的平方和
xim)(i?1,,n),n?m。记拟合误差?i(?)?yi?f(xi,?)，求?使误
Q(?)???(?)??[yi?f(xi,?)]2
2ii?1i?1nn最小。
作为无约束非线性规划的特例，解非线性最小二乘拟合可用MATLAB优化工具箱命令lsqnonlin和lsqcurvefit。 8.3.2 非线性回归模型
非线性回归模型记作
y?f(x,?)??,x?(x1,2,xm),??(?0,?1,,?k)
其中f对回归系数?是非线性的，?~N(0,?)。求得回归系数?的最小二乘估计。 MATLAB统计工具箱中非线性回归的命令是： [b,R,J]=nlinfit(x,y, 'model',bo)
输入x是自变量数据矩阵，每列一个向量；y是因变量数据向量；model是模型的函数名（M文件），形式为y?f(b,x)，b为待估系数?；b0是回归系数?的初值。输出b是?的估计值，R是残差，J是用于估计预测误差的Jacobi矩阵。这个命令是依据高斯—牛顿法求解的。
将上面的输出作为命令 Bi=nlparci(b,R,J) 的输入，得到的bi是回归系数?的置信区间。用命令
nlintool(x,y, 'model',b)
可以得到一个交互式画面，其内容和用法与多项式回归的Polytool类似。 例4 酶促反应速度与底物浓度
酶促反应动力学简称酶动力学，主要研究酶促反应速度与底物（即反应物）浓度以及其它因素的关系。在底物浓度很低时酶促反应是一级反应；当底物浓度处于中间范围时，是混合级反应；当底物浓度增加时，向零级反应过渡。某生化系学生为了研究嘌呤霉素在某项酶促反应中对反应速度与底物浓度之间关系的影响，设计了两个实验，一个实验中所使用的酶
是经过嘌呤霉素处理的，而另一个实验所用的酶是未经嘌呤霉素处理的。所得实验数据见表8-5。试根据问题的背景和这些数据建立一个合适的数学模型，来反映这项酶促反应的速度与底物浓度以及嘌呤霉素处理与否之间的关系。
表8-5 嘌呤霉素实验中的反应速度与底物浓度数据
底物浓度（ppm） 0.02 0.06 86 107 0.11 98 123 0.22 0.56 1.10 未处理 67 51 84 反应 处理 76 47 97 速度
分析与假设
115 131 139 159 124 144 152 191 158 160 / 201 207 200 记酶促反应的速度为y,底物浓度为x，二者之间的关系写作y?f(x,?),其中β为参数(β可为一向量)。由酶促反应的基本性质可知，当底物浓度很低时酶促反应是一级反应，此时反应速度大致与底物浓度成正比；而当底物浓度很大，渐近饱和时，反应速度将趋于一个固定值（即零级反应）。下面的两个简单模型具有这种性质：
Michaelis-Menten 模型
y?f(x,?)?指数增长模型
???x??x
y?f(x,?)???(1?e???x)
非线性模型的求解
首先作出给出的经过嘌呤霉素处理和未经处理的反应速度与底物浓度的散点图，可以看出，上述两个模型与实际数据得到的散点图是大致符合的。
我们将主要对前一模型即Michaelis-Menten模型进行详细的分析。首先对经过嘌呤酶素处理的实验数据进行分析，在此基础上，再来讨论是否有更一般的模型来统一刻画处理前后的数据，进而揭示其中的联系。
我们用非线性回归的方法直接估计模型的参数?1，?2，模型的求解可利用MATLAB统计工具箱中的命令进行，使用格式为：
[beta,R,J]=nlinfit(x,y,'model',beta0)
其中输入x为自变量数据矩阵，每列一个变量；y为因变量数据向量；model为模型的M文件名，M函数形式为y=f (beta,x),beta为待估计参数；beta0为给定的参数初值。输出beta为参数估计值，R为残差，J为用于估计预测误差的Jacobi矩阵。参数beta的置信区间用命令 nlparci(beta,R,J)得到。
首先建立函数M文件huaxue.m，非线性模型参数估计的源程序如下： x=[0.02 0.02 0.06 0.06 0.11 0.11 0.22 0.22 0.56 0.56 1.10 1.10]; y=[76 47 97 107 123 139 159 152 191 201 207 200]; beta0=[195.8027 0.04841];
[beta,R,J]=nlinfit(x,y,'huaxue',beta0); betaci=nlparci(beta,R,J);