回归分析11：含定性因变量的回归模型

目录

• Chapter 11：含定型因变量的回归模型
  • 6.1 含定型因变量的线性回归模型
    • 6.1.1 含定性因变量的回归函数的意义
    • 6.1.2 含定性因变量的回归模型的特殊性
    • 6.1.3 Logistic 回归模型
  • 6.2 Logistic 回归模型的参数估计
    • 6.2.1 分组数据情形
    • 6.2.2 未分组数据情形

Chapter 11：含定型因变量的回归模型

6.1 含定型因变量的线性回归模型

6.1.1 含定性因变量的回归函数的意义

这里我们研究因变量只有两个可能结果的情况，这样的因变量可以用取值为 \(0\) 或 \(1\) 的虚拟变量来表示。

假设因变量 \(y\) 为只取 \(0\) 和 \(1\) 两个值的定性变量，考虑如下的简单线性回归模型

\[y_i=\beta_0+\beta_1x_i+e_i \ . \]

我们通常假设 \({\rm E}(e_i)=0\) ，在因变量只取 \(0\) 和 \(1\) 两个值时，假设

\[P(y_i=1)=\pi_i \ , \quad P(y_i=0)=1-\pi_i \ , \]

于是 \({\rm E}(y_i|x_i)=\pi_i\) ，所以

\[{\rm E}(y_i|x_i)=\pi_i=\beta_0+\beta_1x_i \ . \]

这表明回归函数 \({\rm E}(y_i|x_i)=\beta_0+\beta_1x_i\) 是给定自变量取值为 \(x_i\) 时 \(y_i=1\) 的概率。

6.1.2 含定性因变量的回归模型的特殊性

(1) 离散非正态误差项。对只取 \(0\) 和 \(1\) 的定性因变量 \(y\) ，若它关于自变量 \(x\) 满足简单回归模型

\[y_i=\beta_0+\beta_1x_i+e_i \ , \]

则其误差项 \(e_i\) 也是只能取两个值的离散变量。即

• 当 \(y_i=1\) 时，\(e_i=1-\beta_0-\beta_1x_i=1-\pi_i\) ；
• 当 \(y_i=0\) 时，\(e_i=0-\beta_0-\beta_1x_i=-\pi_i\) ；

所以，误差项为两点分布的随机变量，故正态误差回归模型的假定不在适用。

(2) 误差项仍保持零均值性质但异方差。利用离散型随机变量的性质计算 \(e_i\) 的数学期望，

\[{\rm E}(e_i)=(1-\pi_i)\pi_i-\pi_i(1-\pi_i)=0 \ , \]

说明误差项 \(e_i\) 仍然具有零均值性质，但是

\[{\rm Var}(e_i)={\rm Var}(y_i)=\pi_i(1-\pi_i)=(\beta_0+\beta_1x_i)(1-\beta_0-\beta_1x_i) \ . \]

说明误差项 \(e_i\) 的方差不相等。因此，模型误差为异方差，不满足线性回归模型的基本假定。这表明，对因变量为定性变量的线性回归模型，最小二乘估计的效果不会很好。

(3) 回归函数的限制。当因变量 \(y\) 为只取 \(0\) 和 \(1\) 两个值的定性变量时，由于 \({\rm E}(y_i|x_i)\) 表示给定自变量取值为 \(x_i\) 时 \(y_i=1\) 的概率，所以 \({\rm E}(y_i|x_i)\) 受如下限制：

\[0\leq{\rm E}(y_i|x_i)=\pi_i=\beta_0+\beta_1x_i\leq 1 \ . \]

然而，一般回归函数并不具有这种限制。也就是说，对定性因变量直接建立回归模型，等号右边得到的拟合值是可以超出此范围限制的，因此这是不可取的而且得不到合理的解释的。

6.1.3 Logistic 回归模型

当因变量 \(y\) 为只取 \(0\) 和 \(1\) 两个值的定性变量时，如果我们对影响 \(y\) 的因素 \(x_1,x_2,\cdots,x_p\) 建立线性回归模型，则会遇到以下两个问题：

1. 因变量 \(y\) 本身为只取 \(0\) 和 \(1\) 的离散值，而等号右端的取值可在某个范围内连续变化；
2. 因变量 \(y\) 的取值最大为 \(1\) 最小为 \(0\) ，而等号右端的取值可超出 \([0,1]\) ，甚至可能在 \((-\infty,\infty)\) 上取值。

对于上述的第一个问题，我们可以考虑因变量的均值。对于上述的第二个问题，可以考虑因变量均值的某个函数，使得该函数的取值范围是 \((-\infty,\infty)\) 。其中，最常用的就是 Logit 函数：

\[{\rm Logit}(z)=\ln\frac{z}{1-z} \ , \quad z\in[0,1] \ . \]

因此，我们可以考虑如下的模型：

\[{\rm Logit}(\pi_i)=\ln\frac{\pi_i}{1-\pi_i}=x_i'\beta \ , \quad i=1,2,\cdots,n \ , \]

其中 \(x_i'=(1,x_{i1},\cdots,x_{ip}),\,\beta=\left(\beta_0,\beta_1,\cdots,\beta_p\right)'\) 。该模型可以等价地写为

\[{\rm E}(y_i|X)=\pi_i=\frac{\exp\left\{\beta_0+\sum_{j=1}^n\beta_jx_{ij}\right\}}{1+\exp\left\{\beta_0+\sum_{j=1}^n\beta_jx_{ij}\right\}} \ , \quad i=1,2,\cdots,n \ . \]

这个模型被称为 Logistic 回归模型。

我们将 \(\pi_i/(1-\pi_i)\) 称为“事件发生”与“事件没有发生”的优势比，它是 \(\pi_i\) 的严格增函数。

6.2 Logistic 回归模型的参数估计

6.2.1 分组数据情形

对于 Logistic 回归模型的参数估计问题，首先考虑分组数据情形。假设某一事件 \(A\) 发生的概率 \(\pi\) 依赖于自变量 \(x_1,x_2,\cdots,x_p\) ，我们对事件 \(A\) 在 \(m\) 个不同的自变量组合条件下作了 \(n\) 次观测，其中对应于一个组合 \(x_i=(x_{i1},x_{i2},\cdots,x_{ip})'\) 观测了 \(n_i\) 个结果，\(i=1,2,\cdots,m\) ，且有 \(\sum_{i=1}^mn_i=n\) 。

假设在这 \(n_i\) 个观测中，事件 \(A\) 发生了 \(r_i\) 次，于是事件 \(A\) 发生的概率可以用 \(\hat\pi_i=r_i/n_i\) 来估计。我们把这种结构的数据称为分组数据。用 \(\pi_i\) 的估计值 \(\hat\pi_i\) 作为替代，可以得到

\[y_i^*\xlongequal{def}\ln\frac{\hat\pi_i}{1-\hat\pi_i}=\ln\frac{\pi_i}{1-\pi_i}+e_i=x_i'\beta+e_i \ , \quad i=1,2,\cdots,m \ . \]

这是我们常见的线性回归模型。若假设 \(e_1,e_2,\cdots,e_m\) 互不相关，且 \({\rm E}(e_i)=0\) 和 \({\rm Var}(e_i)=v_i\) ，则参数 \(\beta\) 的广义最小二乘估计为

\[\hat\beta=\left(X'V^{-1}X\right)^{-1}X'V^{-1}Y^* \ , \]

其中

\[Y^*=\begin{bmatrix} y_1^* \\ y_2^* \\ \vdots \\ y_m^* \end{bmatrix} \ , \quad X=\begin{bmatrix} 1 & x_{11} & x_{12} & \cdots & x_{1p} \\ 1 & x_{21} & x_{22} & \cdots & x_{2p} \\ \vdots & \vdots & \vdots & & \vdots \\ 1 & x_{m1} & x_{m2} & \cdots & x_{mp} \\ \end{bmatrix} \ , \quad V=\begin{bmatrix} v_1 & 0 & \cdots & 0 \\ 0 & v_2 & \cdots & 0 \\ \vdots & \vdots & & \vdots \\ 0 & 0 & \cdots & v_m \end{bmatrix} \ . \]

要考察 \(x_j\) 对事件 \(A\) 发生的概率是否有影响，也需要检验 \(x_j\) 的回归系数 \(\beta_j=0\) 这一假设是否成立。但是这需要一个前提，即 \(e_i\) 满足正态性假设。在定性因变量的情况下，下面我们证明这一假设在大样本下仍然成立。

**引理 (Delta Method) **：令 \(\{Y_n\}\) 是一列随机变量满足

\[\sqrt{n}(Y_n-\theta)\xrightarrow{d}N(0,\sigma^2) \ . \]

对于给定函数的函数 \(g(\cdot)\) 和给定的特殊值 \(\theta\) ，假设 \(g'(\theta)\) 存在且非零，则有

\[\sqrt{n}\left(g(Y_n)-g(\theta)\right)\xrightarrow{d} N\left(0,\sigma^2[g'(\theta)]^2\right) \ . \]

由于 \(\hat\pi_i=r_i/n_i\) 是样本的频率，因此由大数定律和中心极限定理可知，当 \(n_i\to\infty\) 时，\(\hat\pi_i\) 以概率 \(1\) 收敛到 \(\pi_i\) ，且有

\[\sqrt{n_i}\left(\hat\pi_i-\pi_i\right)\xrightarrow{d} N(0,\pi_i(1-\pi_i)) \ . \]

下面我们继续推导 \(y_i^*\) 的极限分布。由 \(f(z)=\ln\dfrac{z}{1-z}\) 可得

\[f'(z)=\frac{1}{z(1-z)} \ , \quad f'(z)\bigg|_{z=\pi_i}=\frac{1}{\pi_i(1-\pi_i)} \ . \]

由 Delta 方法可知，当 \(n_i\to\infty\) 时，有

\[\sqrt{n_i}\left(\ln\frac{\hat\pi_i}{1-\hat\pi_i}-\ln\frac{\pi_i}{1-\pi_i}\right)\xrightarrow{d} N\left(0,\frac{1}{\pi_i(1-\pi_i)}\right) \ . \]

这表明，当 \(\min\{n_1.n_2,\cdots,n_m\}\) 充分大时，我们可以认为 \(y_i^*\) 服从正态分布 \(N\left(x_i'\beta,v_i\right)\) ，其中

\[v_i=\frac{1}{n_i\pi_i(1-\pi_i)} \ . \]

由于 \(\pi_i\) 是未知的，我们用 \(\hat\pi_i\) 代替 \(\pi_i\) ，得到

\[ \hat v_i=\frac{1}{n_i\hat\pi_i(1-\hat\pi_i)} \ , \]

用 \(\hat v_i\) 代替 \(V\) 中的 \(v_i\) ，从而可以得到广义最小二乘估计，以及回归系数的检验统计量。

6.2.2 未分组数据情形

假设 \(y_i\sim B(1,\pi_i)\) ，而 \(x_1,x_2,\cdots,x_p\) 是对 \(y\) 有影响的 \(p\) 个自变量。在 \(( x_1,x_2,\cdots,x_p)\) 的 \(n\) 个不同的样本点，对 \(y\) 进行了 \(n\) 次独立观测得到 \(y_1,y_2,\cdots,y_n\) ，显然 \(y_1,y_2,\cdots,y_n\) 是相互独立的两点分布随机变量，\(y_i\) 的概率分布为

\[\pi_i^{y_i}(1-\pi_i)^{1-y_i} \ , \quad y_i=0,1 \ . \]

于是 \(y_1,y_2,\cdots,y_n\) 的似然函数为

\[L(\pi_1,\pi_2,\cdots,\pi_n)=\prod_{i=1}^n\pi_i^{y_i}(1-\pi_i)^{1-y_i} \ . \]

其对数似然函数为

\[l(\pi_1,\pi_2,\cdots,\pi_n)=\sum_{i=1}^ny_i\ln\pi_i+(1-y_i)\ln(1-\pi_i) \ . \]

将 \(\ln\dfrac{\pi_i}{1-\pi_i}=x_i'\beta\) 代入上式得

\[l(\beta)=\sum_{i=1}^n\left[y_ix_i'\beta-\ln(1+\exp\left\{x_i'\beta\right\})\right] \ . \]

求 \(\beta\) 的极大似然估计，就是寻找 \(\beta\) 使得 \(l(\beta)\) 达到最大，为此对其一阶求导

\[\frac{\partial l(\beta)}{\partial\beta}=\sum_{i=1}^n\left(y_i-\frac{\exp\left\{x_i'\beta\right\}}{1+\exp\left\{x_i'\beta\right\}}\right)x_i=X'\varepsilon \ . \]

其中

\[\begin{aligned} &X=(x_1,x_2,\cdots,x_n)'=\begin{bmatrix} 1 & x_{11} & x_{12} & \cdots & x_{1p} \\ 1 & x_{21} & x_{22} & \cdots & x_{2p} \\ \vdots & \vdots & \vdots & & \vdots \\ 1 & x_{m1} & x_{m2} & \cdots & x_{mp} \\ \end{bmatrix} \ , \\ \\ &\varepsilon=(\varepsilon_1,\varepsilon_2,\cdots,\varepsilon_n)' \ , \quad \varepsilon_i=y_i-\frac{\exp\left\{x_i'\beta\right\}}{1+\exp\left\{x_i'\beta\right\}} \ . \end{aligned} \]

令 \(X'\varepsilon=0\) ，求解方程

\[X'\varepsilon=\sum_{i=1}^n\left(y_i-\frac{\exp\left\{x_i'\beta\right\}}{1+\exp\left\{x_i'\beta\right\}}\right)x_i=0 \ , \]

即可得到 \(\beta\) 的极大似然估计。

但上述方程是关于参数 \(\beta\) 的一个较复杂的非线性函数，要获得 \(\beta\) 的极大似然估计 \(\hat\beta\) ，一般需要使用迭代算法，如 Newton-Raphson 迭代算法，求得数值解。