UA MATH571A 多元线性回归I 模型设定与推断

模型设定最小二乘法（Method of Least Square）系数Mean Response and Residual多元回归的ANOVAANOTA TableF检验与可决系数序贯方差分析增量平方和单个系数的F检验部分系数的F检验Lack-of-fit检验多元回归的统计推断回归模型标准化

模型设定

多元线性回归假定引起变量Y改变的因素有p-1个，即变量X∈Rp−1X\in R^{p-1}X∈Rp−1。变量X与Y有一组观测值(Xi,Yi)(X_i, Y_i)(Xi​,Yi​)，i=1,2,...,Ni = 1,2, ... , Ni=1,2,...,N，Xi∈Rp−1X_i \in R^{p-1}Xi​∈Rp−1，Yi∈RY_i \in RYi​∈R并且观测值满足如下关系：

Yi=β0+∑i=1p−1βiXi+ϵiY_i = \beta_0 + \sum_{i=1}^{p-1} \beta_i X_i + \epsilon_i Yi​=β0​+i=1∑p−1​βi​Xi​+ϵi​

可以将样本写成矩阵，Y=[Y1,Y2,...,YN]T∈RNY=[Y_1, Y_2, ... , Y_N]^T \in R^NY=[Y1​,Y2​,...,YN​]T∈RN, X=[1T,X1T,...,XNT]T∈RN×pX = [\textbf{1}^T, X_1^T, ... ,X_N^T]^T \in R^{N \times p}X=[1T,X1T​,...,XNT​]T∈RN×p, 其中1\textbf{1}1表示元素全为1的向量。矩阵X被称为设计矩阵（Design Matrix）。之所以称为设计矩阵是因为在回归模型中，解释变量并不是随机的，因此矩阵X的值完全取决于抽样方式的设计。残差项服从向量正态分布，

ϵ=[ϵ1,ϵ2,...,ϵN]T∼N(0,Σ)\epsilon = [\epsilon_1, \epsilon_2, ..., \epsilon_N]^T\sim N(0,\Sigma) ϵ=[ϵ1​,ϵ2​,...,ϵN​]T∼N(0,Σ)

一元线性回归中对残差项还有独立和同方差的假设，在多元回归中这个状态可以被一个表达式表示

Σ=σ2IN\Sigma=\sigma^2I_N Σ=σ2IN​

定义系数向量β=[β0,β1,...,βp−1]T∈Rp\beta = [\beta_0, \beta_1, ..., \beta_{p-1}]^T \in R^{p}β=[β0​,β1​,...,βp−1​]T∈Rp，单个回归系数，例如β1\beta_1β1​，其含义是在控制其他因素不变的情况下，变量X1X_1X1​增加一单位造成的变量Y的变化的均值。因为

β1=∂Y∂X1=lim⁡h→0Y(X1+h,X2,...,Xp)−Y(X1,X2,...,Xp)h\beta_1 = \frac{\partial Y}{\partial X_1} = {\lim_{h\to 0}} \frac{Y(X1+h,X_2,...,X_p)-Y(X1,X_2,...,X_p)}{h} β1​=∂X1​∂Y​=h→0lim​hY(X1+h,X2​,...,Xp​)−Y(X1,X2​,...,Xp​)​

从数学上看，β1\beta_1β1​衡量的就是在其他变量不发生改变的情况下，变量X1X_1X1​的微小变化对变量Y的影响。所以多元回归从数学上看就相当于控制实验。

回归方程可以写成矩阵形式，

Y=Xβ+ϵY = X\beta + \epsilon Y=Xβ+ϵ

事实上除了这种设定外还有其他模型设定方式。例如针对解释变量的取值可以将其分为实值变量和虚拟变量（dummy），或者称为数量型变量（quantitative variable）和质量型变量（qualitative variable）。从解释变量的阶可以将其分为线性回归（所有解释变量都是一阶的），多项式回归（polynomial regression，模型包含某些解释变量的二次方或者高次方项）。模型同样可以包含不同解释变量之间的交互效应，如

Y=β0+β1X1+β2X2+β3X1X2Y = \beta_0+ \beta_1 X_1 + \beta_2 X_2 + \beta_3 X_1 X_2 Y=β0​+β1​X1​+β2​X2​+β3​X1​X2​

就包含解释变量X1X_1X1​和X2X_2X2​的交互项X1X2X_1 X_2X1​X2​。之所以称为交互效应是因为解释变量X1X_1X1​对被解释变量Y的效应

∂Y∂X1=β1+β3X2\frac{\partial Y}{\partial X_1} = \beta_1 + \beta_3 X_2 ∂X1​∂Y​=β1​+β3​X2​

是与解释变量X2X_2X2​相关的。因此在设计回归模型时，可以根据需要加入不同的回归项，或者对数据做合理的变换。需要注意的是线性回归中线性的概念，将回归方程写成关于系数向量的函数Y=f(β)Y=f(\beta)Y=f(β)，则线性指的是fff是关于β\betaβ的一个线性映射。

最小二乘法（Method of Least Square）

多元回归模型的残差平方和为

Q=(Y−Xβ)2Q = (Y - X\beta)^2 Q=(Y−Xβ)2

所以最小二乘法的目标是

minβQ=(Y−Xβ)2=YTY+βTXTXβ−YTXβ−βTXTY{min_{\beta}} Q=(Y - X\beta)^2 = Y^TY + \beta^TX^TX\beta - Y^TX\beta - \beta^TX^TY minβ​Q=(Y−Xβ)2=YTY+βTXTXβ−YTXβ−βTXTY

系数

有两个比较有用的对向量求导数的结论，假设A∈Rn×nA \in R^{n \times n}A∈Rn×n是一个矩阵，v∈Rnv \in R^nv∈Rn和x∈Rnx \in R^nx∈Rn是列向量，则

∇x(Ax)=A∇x(vTx)=vT\nabla_x (Ax) = A \\ \nabla_x (v^Tx) = v^T ∇x​(Ax)=A∇x​(vTx)=vT

∇βQ=(XTX+(XTX)T)β^−(YTX)T−XTY=2XTXβ^−2XTY=0β^=(XTX)−1XTY\nabla_{\beta} Q = (X^TX + (X^TX)^T)\hat{\beta} - (Y^TX)^T - X^TY=2X^TX\hat{\beta}-2X^TY=0 \\ \hat{\beta} = (X^TX)^{-1}X^TY ∇β​Q=(XTX+(XTX)T)β^​−(YTX)T−XTY=2XTXβ^​−2XTY=0β^​=(XTX)−1XTY

上式是系数向量的最小二乘估计量。

Mean Response and Residual

Mean Response可以写成

EY=Xβ^=X(XTX)−1XTYEY = X\hat{\beta} = X(X^TX)^{-1}X^TY EY=Xβ^​=X(XTX)−1XTY

定义H=X(XTX)−1XTH = X(X^TX)^{-1}X^TH=X(XTX)−1XT，称为Hat matrix。其本质为投影矩阵，代表RNR^NRN上的一个线性变换，该变换将RNR^NRN中的任意向量投影到设计矩阵张成的线性空间span(X)中。Hat matrix具有对称性和幂等性：

HT=[X(XTX)−1XT]T=X(XTX)−1XT=HHH=X(XTX)−1XTX(XTX)−1XT=X(XTX)−1XT=HH^T = [X(X^TX)^{-1}X^T]^T = X(X^TX)^{-1}X^T = H \\ HH = X(X^TX)^{-1}X^TX(X^TX)^{-1}X^T=X(X^TX)^{-1}X^T=H HT=[X(XTX)−1XT]T=X(XTX)−1XT=HHH=X(XTX)−1XTX(XTX)−1XT=X(XTX)−1XT=H

残差可以写成

e=Y−Xβ^=Y−X(XTX)−1XTY=(IN−H)Ye = Y - X\hat{\beta} = Y - X(X^TX)^{-1}X^TY = (I_N - H)Y e=Y−Xβ^​=Y−X(XTX)−1XTY=(IN​−H)Y

定义M=IN−HM = I_N - HM=IN​−H, 则M同样是投影矩阵，它代表的线性变换将RNR^NRN中的任意向量投影到设计矩阵张成的线性空间span(X)的正交补空间中。该矩阵同样具有对称性和幂等性：

MT=(IN−H)T=INT−HT=IN−H=MMM=(IN−H)(IN−H)=IN−2H+H=IN−H=MM^T = (I_N - H)^T = I_N^T - H^T = I_N - H=M \\ MM = (I_N-H)(I_N-H)=I_N - 2H + H = I_N - H=M MT=(IN​−H)T=INT​−HT=IN​−H=MMM=(IN​−H)(IN​−H)=IN​−2H+H=IN​−H=M

关于残差与拟合值同样有如下三条性质：

性质一：残差和为零。（同样由Q关于β^0\hat{\beta}_0β^​0​的偏导为0导出）

∂Q∂β0=−2[1T(Y−Xβ^)]=01T(Y−Xβ^)=1Te=0\frac{\partial Q}{\partial \beta_0} = -2[\textbf{1}^T(Y-X\hat{\beta})]=0 \\ \textbf{1}^T(Y-X\hat{\beta}) = \textbf{1}^T e=0 ∂β0​∂Q​=−2[1T(Y−Xβ^​)]=01T(Y−Xβ^​)=1Te=0

性质二：残差与被解释变量正交（由Q关于β^\hat{\beta}β^​的偏导为0导出）

2XTXβ^−2XTY=0=2XT(Xβ^−Y)XTe=02X^TX\hat{\beta}-2X^TY=0=2X^T(X\hat{\beta}-Y) \\ X^T e = \textbf{0} 2XTXβ^​−2XTY=0=2XT(Xβ^​−Y)XTe=0

性质三：残差与拟合值正交（由性质一和性质二导出）

Y^Te=β^TXTe=β^T0=0\hat{Y}^Te = \hat{\beta}^TX^Te=\hat{\beta}^T \textbf{0} = 0 Y^Te=β^​TXTe=β^​T0=0

多元回归的ANOVA

根据总平方和的定义，

SST=∑i=1N(Yi−Yˉ)2=∑i=1NYi2−(∑i=1NYi)2NSST = \sum_{i=1}^N (Y_i-\bar{Y})^2 = \sum_{i=1}^N Y_i^2 - \frac{(\sum_{i=1}^{N} Y_i)^2}{N} SST=i=1∑N​(Yi​−Yˉ)2=i=1∑N​Yi2​−N(∑i=1N​Yi​)2​

其中∑i=1NYi2=YTY\sum_{i=1}^N Y_i^2 = Y^TY∑i=1N​Yi2​=YTY，假设JJJ是所有元素均为1的矩阵，则

(∑i=1NYi)2=YTJY(\sum_{i=1}^{N} Y_i)^2 = Y^TJY (i=1∑N​Yi​)2=YTJY

所以用矩阵表示的总平方和为

SST=YT(IN−JN)YSST = Y^T(I_N-\frac{J}{N})Y SST=YT(IN​−NJ​)Y

残差平方和为

SSE=eTe=(Y−Xβ^)T(Y−Xβ^)=YTY−YTXβ^−β^TXTY+β^TXTXβ^SSE = e^Te=(Y-X\hat{\beta})^T(Y-X\hat{\beta})=Y^TY-Y^TX\hat{\beta}-\hat{\beta}^TX^TY+\hat{\beta}^TX^TX\hat{\beta} SSE=eTe=(Y−Xβ^​)T(Y−Xβ^​)=YTY−YTXβ^​−β^​TXTY+β^​TXTXβ^​

其中

β^TXTXβ^=YTX(XTX)−1(XTX)(XTX)−1XTY=YTXβ^\hat{\beta}^TX^TX\hat{\beta} = Y^TX(X^TX)^{-1}(X^TX)(X^TX)^{-1}X^TY = Y^TX\hat{\beta} β^​TXTXβ^​=YTX(XTX)−1(XTX)(XTX)−1XTY=YTXβ^​

SSE=YTY−β^TXTYSSE = Y^TY - \hat{\beta}^TX^TY SSE=YTY−β^​TXTY

所以回归平方和为

SSR=SST−SSE=β^TXTY−1NYTJYSSR = SST - SSE = \hat{\beta}^TX^TY - \frac{1}{N}Y^TJY SSR=SST−SSE=β^​TXTY−N1​YTJY

ANOTA Table

根据上述推导，可以写出下列方差分析表

F检验与可决系数

考虑对系数向量的双边检验：

H0:β=0Ha:β≠0H_0: \beta = 0 \\ H_a: \beta \ne 0 H0​:β=0Ha​:β​=0

定义统计量

F∗=MSRMSEF^* = \frac{MSR}{MSE} F∗=MSEMSR​

在原假设下，F∗∼(p−1,N−p)F^* \sim (p-1,N-p)F∗∼(p−1,N−p)。假设检验水平为α\alphaα，若F∗≤F(1−α;p−1,N−p)F^*\le F(1-\alpha;p-1,N-p)F∗≤F(1−α;p−1,N−p)，接受原假设，若F∗>F(1−α;p−1,N−p)F^*>F(1-\alpha;p-1,N-p)F∗>F(1−α;p−1,N−p)，拒绝原假设。多元回归的可决系数也采用类似的定义

R2=SSRSST=1−SSESSTR^2 = \frac{SSR}{SST}=1-\frac{SSE}{SST} R2=SSTSSR​=1−SSTSSE​

但这个定义有一个固有的缺陷。当增加新的解释变量时，即使新的解释变量并没有解释力，R2R^2R2也是会增加的。为了弥补这个缺陷，多元回归通常使用调整的可决系数

Ra2=1−SSE/N−pSST/N−1R_a^2 =1-\frac{SSE/N-p}{SST/N-1} Ra2​=1−SST/N−1SSE/N−p​

在新加入的解释变量没有任何解释力时，调整可决系数Ra2R_a^2Ra2​会因为自由度调整而变小，因此多元回归一般看调整可决系数。

序贯方差分析

通常在用多元回归建模时，我们总是会考虑某个解释变量是不是应该加入到回归模型中。根据方差分析的思想，总平方和代表被解释变量的总信息，回归平方和代表能够被回归模型解释的信息，残差平方和代表不能被回归模型解释的那部分信息。以回归模型Y=β0+β1X1+ϵY = \beta_0 + \beta_1X_1 + \epsilonY=β0​+β1​X1​+ϵ为例，要考虑是否应该添加另一个解释变量X2X_2X2​，只需要考察增加X2X_2X2​后，回归平方和能够增加多少即可。如果回归平方和几乎没有变化，则加入X2X_2X2​不能给模型带来更多的解释力。

增量平方和

定义增量平方和（Extra Sum of Square）

SSR(X2∣X1)=SSR(X1,X2)−SSR(X1)SSR(X_2|X_1) = SSR(X_1,X_2)-SSR(X_1) SSR(X2​∣X1​)=SSR(X1​,X2​)−SSR(X1​)

其含义为在回归模型Y=β0+β1X1+ϵY = \beta_0 + \beta_1X_1 + \epsilonY=β0​+β1​X1​+ϵ中加入新的解释变量X2X_2X2​能够带来的回归平方和的增量。原来的一元回归的回归平方和为SSR(X1)SSR(X_1)SSR(X1​)，加入X2X_2X2​后新的二元回归额回归平方和为SSR(X1,X2)SSR(X_1,X_2)SSR(X1​,X2​)。增量平方和的自由度等于打算添加的解释变量的数量，比如SSR(X2∣X1)SSR(X_2|X_1)SSR(X2​∣X1​)考察的是要不要添加X2X_2X2​，因此增量平方和自由度为1，

MSR(X2∣X1)=SSR(X2∣X1)1MSR(X_2|X_1)=\frac{SSR(X_2|X_1)}{1} MSR(X2​∣X1​)=1SSR(X2​∣X1​)​

增量平方和有一个明显的性质：同时考虑添加多个解释变量时，增量平方和等于序贯添加这些解释变量带来的增量平方和

SSR(X2,X3∣X1)=SSR(X1,X2,X3)−SSR(X1)=SSR(X1,X2,X3)−SSR(X1,X2)+SSR(X1,X2)−SSR(X1)=SSR(X2∣X1)+SSR(X3∣X1,X2)SSR(X_2,X_3|X_1)=SSR(X_1,X_2,X_3)-SSR(X_1) \\ =SSR(X_1,X_2,X_3) - SSR(X_1,X_2) + SSR(X_1,X_2)-SSR(X_1)\\=SSR(X_2|X_1)+SSR(X_3|X_1,X_2) SSR(X2​,X3​∣X1​)=SSR(X1​,X2​,X3​)−SSR(X1​)=SSR(X1​,X2​,X3​)−SSR(X1​,X2​)+SSR(X1​,X2​)−SSR(X1​)=SSR(X2​∣X1​)+SSR(X3​∣X1​,X2​)

单个系数的F检验

在多元回归中，F检验也可以用来做单个系数的双边检验。考虑对系数βi\beta_iβi​的双边检验：

H0:βk=0Ha:βk≠0H_0: \beta_k = 0 \\ H_a: \beta_k \ne 0 H0​:βk​=0Ha​:βk​​=0

根据广义线性检验方法，多元回归的完整模型为

FM：Yi=β0+β1X1+...+βkXk+...+βp−1Xp−1+ϵiFM：Y_i = \beta_0 + \beta_1X_1+...+\beta_kX_k + ... +\beta_{p-1} X_{p-1}+ \epsilon_i FM：Yi​=β0​+β1​X1​+...+βk​Xk​+...+βp−1​Xp−1​+ϵi​

原假设下的简化模型为

RM：Yi=β0+β1X1+...+βk−1Xk−1+βk+1Xk+1+...+βp−1Xp−1+ϵiRM：Y_i = \beta_0 + \beta_1X_1+...+\beta_{k-1}X_{k-1} + \beta_{k+1}X_{k+1} + ... +\beta_{p-1} X_{p-1}+ \epsilon_i RM：Yi​=β0​+β1​X1​+...+βk−1​Xk−1​+βk+1​Xk+1​+...+βp−1​Xp−1​+ϵi​

定义统计量

F∗=SSE(RM)−SSE(FM)dfRM−dfFM/SSE(FM)dfFM∼F(dfRM−dfFM,dfFM)F^* = \frac{SSE(RM)-SSE(FM)}{df_{RM}-df_{FM}}/\frac{SSE(FM)}{df_{FM}} \sim F(df_{RM}-df_{FM},df_{FM}) F∗=dfRM​−dfFM​SSE(RM)−SSE(FM)​/dfFM​SSE(FM)​∼F(dfRM​−dfFM​,dfFM​)

其中SSE(RM)=SSE(X1,...,Xk−1,Xk+1,...Xp−1)SSE(RM)=SSE(X_1,...,X_{k-1},X_{k+1},...X_{p-1})SSE(RM)=SSE(X1​,...,Xk−1​,Xk+1​,...Xp−1​)，SSE(FM)=SSE(X1,...,Xp−1)SSE(FM)=SSE(X_1,...,X_{p-1})SSE(FM)=SSE(X1​,...,Xp−1​)。所以

SSE(RM)−SSE(FM)=SST(RM)−SSR(RM)−[SST(FM)−SSR(FM)]=SSR(FM)−SSR(RM)=SSR(FM∖RM∣RM)=SSR(Xk∣RM)SSE(RM)-SSE(FM)=SST(RM)-SSR(RM)-[SST(FM)-SSR(FM)]\\=SSR(FM)-SSR(RM) = SSR({FM\setminus RM}|RM) = SSR(X_k|RM) SSE(RM)−SSE(FM)=SST(RM)−SSR(RM)−[SST(FM)−SSR(FM)]=SSR(FM)−SSR(RM)=SSR(FM∖RM∣RM)=SSR(Xk​∣RM)

F∗=SSR(Xk∣RM)1/SSE(FM)dfFM=MSR(Xk∣X1,...,Xk−1,Xk+1,...Xp−1)MSE(X1,...,Xp−1)∼F(1,N−p)F^* = \frac{SSR(X_k|RM)}{1}/\frac{SSE(FM)}{df_{FM}}\\=\frac{MSR(X_k|X_1,...,X_{k-1},X_{k+1},...X_{p-1})}{MSE(X_1,...,X_{p-1})} \sim F(1,N-p) F∗=1SSR(Xk​∣RM)​/dfFM​SSE(FM)​=MSE(X1​,...,Xp−1​)MSR(Xk​∣X1​,...,Xk−1​,Xk+1​,...Xp−1​)​∼F(1,N−p)

假设检验水平为α\alphaα，若F∗≤F(1−α;1,N−p)F^*\le F(1-\alpha;1,N-p)F∗≤F(1−α;1,N−p)，接受原假设，若F∗>F(1−α;1,N−p)F^*>F(1-\alpha;1,N-p)F∗>F(1−α;1,N−p)，拒绝原假设。

部分系数的F检验

考虑如下双边检验：

H0:βk1=βk2=...=βkM=0Ha:notallofβk1,βk2,...βkMequaltozeroH_0: \beta_{k_1}=\beta_{k_2}=...=\beta_{k_M} = 0 \\ H_a: not\ all\ of\ \beta_{k_1}, \beta_{k_2},... \beta_{k_M}\ equal\ to\ zero H0​:βk1​​=βk2​​=...=βkM​​=0Ha​:notallofβk1​​,βk2​​,...βkM​​equaltozero

定义统计量

F∗=SSR(FM∖RM∣RM)M/SSE(FM)dfFM∼F(M,N−p)F^* = \frac{SSR({FM\setminus RM}|RM) }{M}/\frac{SSE(FM)}{df_{FM}} \sim F(M,N-p) F∗=MSSR(FM∖RM∣RM)​/dfFM​SSE(FM)​∼F(M,N−p)

假设检验水平为α\alphaα，若F∗≤F(1−α;M,N−p)F^*\le F(1-\alpha;M,N-p)F∗≤F(1−α;M,N−p)，接受原假设，若F∗>F(1−α;M,N−p)F^*>F(1-\alpha;M,N-p)F∗>F(1−α;M,N−p)，拒绝原假设。

Lack-of-fit检验

Lack-of-fit检验的目的是检验模型设定是否正确，

H0:E(Yij)=β0+∑i=1p−1βiXiHa:E(Yij)≠β0+∑i=1p−1βiXiH_0:E(Y_{ij})=\beta_0 + \sum_{i=1}^{p-1} \beta_i X_i \\ H_a:E(Y_{ij})\ne \beta_0 + \sum_{i=1}^{p-1} \beta_i X_i H0​:E(Yij​)=β0​+i=1∑p−1​βi​Xi​Ha​:E(Yij​)​=β0​+i=1∑p−1​βi​Xi​

定义F统计量

F∗=SSLFc−2/SSPEN−c=MSLFMSPEF^* = \frac{SSLF}{c-2}/\frac{SSPE}{N-c} = \frac{MSLF}{MSPE} F∗=c−2SSLF​/N−cSSPE​=MSPEMSLF​

假设检验水平为α\alphaα，若F∗≤F(1−α;c−2,N−c)F^*\le F(1-\alpha;c-2,N-c)F∗≤F(1−α;c−2,N−c)，接受原假设，若F∗>F(1−α;c−2,N−c)F^*>F(1-\alpha;c-2,N-c)F∗>F(1−α;c−2,N−c)，拒绝原假设。

多元回归的统计推断

多元回归做统计推断的原理与一元回归基本一致，下面用一个例子来说明。

回归模型标准化

为了避免样本数据单位对回归系数造成影响（事实上单位不会决定回归系数，决定的是回归系数的单位，单位不合适可能导致回归系数数值过大或者过小，进而导致计算过程中出现较明显舍入误差），可以将样本数据做标准化：

Yi∗=1N−1Yi−Yˉse(Y)Xji∗=1N−1Xji−Xˉjse(Xj)Y_i^* = \frac{1}{\sqrt{N-1}} \frac{Y_i - \bar{Y}}{se(Y)} \\ X_{ji}^* = \frac{1}{\sqrt{N-1}} \frac{X_{ji} - \bar{X}_j}{se(X_j)} Yi∗​=N−1​1​se(Y)Yi​−Yˉ​Xji∗​=N−1​1​se(Xj​)Xji​−Xˉj​​

按标准化后的样本数据做回归，则系数与不做标准化的回归系数的关系为

βj=se(Y)se(Xj)βj∗\beta_j= \frac{se(Y)}{se(X_j)} \beta_j^* βj​=se(Xj​)se(Y)​βj∗​