以酶促反应为例的统计回归问题及代码实现

背景和问题

研究嘌呤霉素在某项酶促反应中对反应速度与底物浓度之间关系的影响。设计了两个实验，一个实验中所使用的酶是经过嘌呤霉素处理的，而另一个实验中所使用的酶是未经过嘌呤霉素处理的。

分析与假设

记酶促反应的速度为y，底物浓度为x，二者之间的关系可写作 y=f(x,b),其中b为参数。由酶促反应的基本性质可知，当底物浓度较小时，反应速度大致与浓度成正比；而当底物浓度很大，渐进饱和时，反应速度将趋于固定值。
下面两个简单模型具有这种性质：

Michaelis-Menten模型

$y=f(x,b)= \dfrac{b1x}{b2+x}$

指数增长模型

$y=f(x,b)=b1(1-e^-(b2x))$

线性化模型

由于对参数b=(b1,b2)是[非线性的，但是可以通过下面的变量代换化为线性：
$\dfrac{1}{y}=\dfrac{1}{b1}+\dfrac{b2}{b1}\dfrac{1}{x}=k1+k2u$
此时因变量1/y对新的参数k=(k1,k2)是线性的。可以先单从线性回归模型计算，得到参数k1,k2.根据b与k的关系，可以得到b1和b2的估计值分别为195.802和0.0484.

非线性模型

由于线性模型中，底物浓度小时决定了参数变量，底物浓度大时会有很大偏差。
可利用MATLAB统计工具箱中的命令进行。格式为

[beta,R,J]=nlinfit(x,y,'model',beta0)

其中输入x为自变量的数据矩阵，y为因变量的数据矩阵，model为模型的M函数文件名，形式为 $y=f(beta,x)$,beta为待估参数，beta0为参数初值，R为残差，J为Jacobi矩阵，估计预测误差。

例：
“`beta0=[195.802 0.0484];
[beta,R,J]=nlinfit(x,y,’huaxue’,beta0);
betaci=nlparci(beta,R,J);
beta,betaci
yy=beta(1)*x./(beta(2)+x);
plot(x,y,’o’,x,yy,’+’),pause
nlintool(x,y,’huaxue’,beta)

function yhat=huaxue(beta,x)
yhat=beta(1)*x./(beta(2)+x);
“`

混合反应模型

可以引入示性变量，表示是否经过嘌呤霉素处理。
x1=1表示经过处理,x2=0表示未经过处理。
则
$y=f(x,b)= \dfrac{(b1+r1x2)x1}{(b2+r2x2)+x1}$
仍旧用nlinfit来计算求解。