回歸分析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 迭代算法，求得數值解。