附录 B 线性回归与逻辑回归
B.1 总体普通最小二乘
假定 \((x_i,y_i)_{i=1}^{n}\stackrel{\mathrm{IID}}{\sim}(x,y)\)，其中 \(x\) 是 \(p\) 维随机向量（\(p=1\) 时就是标量），\(y\) 是随机标量。下面用 \((x,y)\) 表示一般观测，省掉下标 \(i\)。定义总体普通最小二乘（population ordinary least squares, OLS）系数为
\[ \beta=\arg\min_{b}\mathrm{E}\bigl\{(y-x^{\mathrm{T}}b)^2\bigr\}. \]
目标函数对 \(b\) 是二次的，于是可以证明：只要矩存在且 \(\mathrm{E}(xx^{\mathrm{T}})\) 可逆，最小化元就是
\[ \beta=\bigl\{\mathrm{E}(xx^{\mathrm{T}})\bigr\}^{-1}\mathrm{E}(xy). \]
有了 \(\beta\)，可以把 \(x^{\mathrm{T}}\beta\) 叫做 \(y\) 在 \(x\) 上的线性投影（linear projection），并把
\[ \varepsilon=y-x^{\mathrm{T}}\beta \tag{B.1} \]
叫做总体残差（population residual）。由 \(\beta\) 的定义，可以核实
\[ \mathrm{E}(x\varepsilon)=\mathrm{E}\bigl\{x(y-x^{\mathrm{T}}\beta)\bigr\}=\mathrm{E}(xy)-\mathrm{E}(xx^{\mathrm{T}})\beta=0. \]
例 B.1（带截距的总体 OLS） 若把 \(1\) 放进 \(x\) 的一个分量，则
\[ \mathrm{E}(\varepsilon)=\mathrm{E}(y-x^{\mathrm{T}}\beta)=0, \]
这进一步蕴含 \(\operatorname{cov}(x,\varepsilon)=0\)。所以 \(\beta\) 里有截距时，总体残差的均值必须是零，并且由构造就与其余协变量不相关。
例 B.2（带截距的一元总体 OLS） 一个重要特例是标量 \(x\) 与 \(y\)：可以定义
\[ (\alpha,\beta)=\arg\min_{a,b}\mathrm{E}\{(y-a-bx)^2\}, \]
它们有显式公式
\[ \beta=\frac{\operatorname{cov}(x,y)}{\operatorname{var}(x)}, \qquad \alpha=\mathrm{E}(y)-\beta\mathrm{E}(x). \]
例 B.3（不带截距的一元总体 OLS） 没有截距时，可以定义
\[ \gamma=\arg\min_{c}\mathrm{E}\{(y-cx)^2\}, \]
它等于
\[ \gamma=\frac{\mathrm{E}(xy)}{\mathrm{E}(x^2)}. \]
当 \(x\) 均值为零时，\(\gamma\) 等于例 B.2 里的 \(\beta\)。
也可以把 (B.1) 改写成
\[ y=x^{\mathrm{T}}\beta+\varepsilon, \tag{B.2} \]
它只来自总体 OLS 系数与残差的定义，不加任何建模假定。我们把 (B.2) 叫做总体 OLS 分解。
B.2 样本普通最小二乘
基于数据 \((x_i,y_i)_{i=1}^{n}\stackrel{\mathrm{IID}}{\sim}(x,y)\)，很容易得到总体 OLS 系数的矩估计
\[ \hat\beta = \left(n^{-1}\sum_{i=1}^{n}x_ix_i^{\mathrm{T}}\right)^{-1} \left(n^{-1}\sum_{i=1}^{n}x_iy_i\right), \]
以及残差 \(\hat\varepsilon_i=y_i-x_i^{\mathrm{T}}\hat\beta\)。这叫做样本 OLS，或干脆叫 OLS。OLS 系数 \(\hat\beta\) 最小化残差平方和
\[ \hat\beta=\arg\min_{b}\,n^{-1}\sum_{i=1}^{n}(y_i-x_i^{\mathrm{T}}b)^2, \]
因此必须满足
\[ \sum_{i=1}^{n}x_i(y_i-x_i^{\mathrm{T}}\hat\beta)=0, \]
有时叫做正规方程（Normal equation）。拟合值，也叫 \(y_i\) 在 \(x_i\) 上的线性投影，等于
\[ \hat y_i=x_i^{\mathrm{T}}\hat\beta\qquad(i=1,\ldots,n). \]
用矩阵记号
\[ X=\begin{pmatrix}x_1^{\mathrm{T}}\\ \vdots\\ x_n^{\mathrm{T}}\end{pmatrix}, \qquad Y=\begin{pmatrix}y_1\\ \vdots\\ y_n\end{pmatrix}, \]
可以把 OLS 系数写成
\[ \hat\beta=(X^{\mathrm{T}}X)^{-1}X^{\mathrm{T}}Y, \]
拟合向量写成
\[ \hat Y=X\hat\beta=X(X^{\mathrm{T}}X)^{-1}X^{\mathrm{T}}Y. \]
定义帽子矩阵（hat matrix）为
\[ H=X(X^{\mathrm{T}}X)^{-1}X^{\mathrm{T}}. \]
于是也有 \(\hat Y=HY\)，这正是「帽子矩阵」这个名字的来历。\(H\) 的对角元 \(h_{ii}\) 常常叫做杠杆值（leverage scores）。
假定 \((x,y)\) 有有限四阶矩，可以用大数定律与 CLT 证明
\[ \sqrt{n}(\hat\beta-\beta)\to\mathrm{N}(0,V) \]
依分布，其中 \(V=B^{-1}MB^{-1}\)，\(B=\mathrm{E}(xx^{\mathrm{T}})\)，\(M=\mathrm{E}(\varepsilon^2 xx^{\mathrm{T}})\)。于是 \(\hat\beta\) 渐近方差的一个矩估计是
\[ \hat V_{\mathrm{ehw}} = n^{-1} \left(n^{-1}\sum_{i=1}^{n}x_ix_i^{\mathrm{T}}\right)^{-1} \left(n^{-1}\sum_{i=1}^{n}\hat\varepsilon_i^2 x_ix_i^{\mathrm{T}}\right) \left(n^{-1}\sum_{i=1}^{n}x_ix_i^{\mathrm{T}}\right)^{-1}, \tag{B.3} \]
叫做 Eicker–Huber–White（EHW） 稳健协方差估计（Eicker, 1967; Huber, 1967; White, 1980）。可以证明 \(n\hat V_{\mathrm{ehw}}\) 依概率收敛到 \(V\)。有了 \(\hat\beta\) 与 \(\hat V_{\mathrm{ehw}}\)，就可以对总体 OLS 系数 \(\beta\) 做推断。
基于杠杆值，EHW 稳健协方差估计还有许多变体（Long and Ervin, 2000）。尤其是：HC1 变体把 \(\hat\varepsilon_i^2\) 改成 \(\hat\varepsilon_i^2/(n-p)\)，HC2 变体改成 \(\hat\varepsilon_i^2/(1-h_{ii})\)，HC3 变体改成 \(\hat\varepsilon_i^2/(1-h_{ii})^2\)，都代进 \(\hat V_{\mathrm{ehw}}\) 的定义。
B.3 Frisch–Waugh–Lovell 定理
Frisch–Waugh–Lovell（FWL）定理有两个版本：一个在总体层面，一个在样本层面。它把多元 OLS 收成一元 OLS，因而便于理解、也便于计算 OLS 系数。下面给本书够用的几个特例。
定理 B.1（总体 FWL） 把 \(y\) 对 \((x_1,x_2,\ldots,x_p)\) 做 OLS 时 \(x_1\) 的系数，等于把 \(y\) 或 \(\tilde y\) 对 \(\tilde x_1\) 做 OLS 时 \(\tilde x_1\) 的系数。其中 \(\tilde y\) 是把 \(y\) 对 \((x_2,\ldots,x_p)\) 做 OLS 的残差，\(\tilde x_1\) 是把 \(x_1\) 对 \((x_2,\ldots,x_p)\) 做 OLS 的残差。
定理 B.1 里，把 \(x_1\) 残差化是要紧的，把 \(y\) 残差化却不是。
定理 B.2（样本 FWL） 数据 \((Y,X_1,X_2,\ldots,X_p)\) 都是列向量。把 \(Y\) 对 \((X_1,X_2,\ldots,X_p)\) 做 OLS 时 \(X_1\) 的系数，等于把 \(Y\) 或 \(\tilde Y\) 对 \(\tilde X_1\) 做 OLS 时 \(\tilde X_1\) 的系数。其中 \(\tilde Y\) 是把 \(Y\) 对 \((X_2,\ldots,X_p)\) 做 OLS 的残差向量，\(\tilde X_1\) 是把 \(X_1\) 对 \((X_2,\ldots,X_p)\) 做 OLS 的残差向量。
同样，定理 B.2 里把 \(X_1\) 残差化是要紧的，把 \(Y\) 残差化却不是。Ding (2021) 给过与样本 FWL 定理有关的更多数值性质。
B.4 线性模型
有时我们加上更强的模型假定：给定 \(x\) 时 \(y\) 的条件均值是线性的，
\[ \mathrm{E}(y\mid x)=x^{\mathrm{T}}\beta, \]
或者等价地，
\[ y=x^{\mathrm{T}}\beta+\varepsilon \qquad\text{且}\qquad \mathrm{E}(\varepsilon\mid x)=0. \]
这叫做受限均值模型（restricted mean model）。在这个模型下，总体 OLS 系数就是真正关心的参数：
\[\begin{align*} \bigl\{\mathrm{E}(xx^{\mathrm{T}})\bigr\}^{-1}\mathrm{E}(xy) &= \bigl\{\mathrm{E}(xx^{\mathrm{T}})\bigr\}^{-1}\mathrm{E}\bigl\{x\mathrm{E}(y\mid x)\bigr\} \\ &= \bigl\{\mathrm{E}(xx^{\mathrm{T}})\bigr\}^{-1}\mathrm{E}(xx^{\mathrm{T}}\beta) \\ &= \beta. \end{align*}\]
再者，总体 OLS 系数并不依赖 \(x\) 的分布。第 B.2 节的渐近推断对这个模型同样适用。
特例 \(\operatorname{var}(\varepsilon\mid x)=\sigma^2\) 时，OLS 系数的渐近方差收成
\[ V=\sigma^2\bigl\{\mathrm{E}(xx^{\mathrm{T}})\bigr\}^{-1}, \]
于是 \(\hat\beta\) 渐近方差有一个更简单的矩估计
\[ \hat V_{\mathrm{ols}} = \hat\sigma^2 \left(\sum_{i=1}^{n}x_ix_i^{\mathrm{T}}\right)^{-1}, \tag{B.4} \]
其中 \(\hat\sigma^2=(n-p)^{-1}\sum_{i=1}^{n}\hat\varepsilon_i^2\) 是 \(\sigma^2\) 的无偏估计；\(n\) 是样本量，\(p\) 是 \(x\) 的维数。这就是 lm 函数给出的标准协方差估计。
基于 BostonHousing 数据，先看 lm 的标准输出。
> library("mlbench")
> data(BostonHousing)
> ols.fit = lm(medv ~ ., data = BostonHousing)
> summary(ols.fit)残差：Min \(-15.595\)，\(1\)Q \(-2.730\)，中位数 \(-0.518\)，\(3\)Q \(1.777\)，Max \(26.199\)。系数如下（原书完整 summary 输出）。
| Estimate | Std. Error | \(t\) | \(\operatorname{Pr}(>\lvert t\rvert)\) | |
|---|---|---|---|---|
| (Intercept) | \(3.646\times 10^{1}\) | \(5.103\) | \(7.144\) | \(3.28\times 10^{-12}\) |
| crim | \(-1.080\times 10^{-1}\) | \(3.286\times 10^{-2}\) | \(-3.287\) | \(0.001087\) |
| zn | \(4.642\times 10^{-2}\) | \(1.373\times 10^{-2}\) | \(3.382\) | \(0.000778\) |
| indus | \(2.056\times 10^{-2}\) | \(6.150\times 10^{-2}\) | \(0.334\) | \(0.738288\) |
| chas1 | \(2.687\) | \(8.616\times 10^{-1}\) | \(3.118\) | \(0.001925\) |
| nox | \(-1.777\times 10^{1}\) | \(3.820\) | \(-4.651\) | \(4.25\times 10^{-6}\) |
| rm | \(3.810\) | \(4.179\times 10^{-1}\) | \(9.116\) | \(<2\times 10^{-16}\) |
| age | \(6.922\times 10^{-4}\) | \(1.321\times 10^{-2}\) | \(0.052\) | \(0.958229\) |
| dis | \(-1.476\) | \(1.995\times 10^{-1}\) | \(-7.398\) | \(6.01\times 10^{-13}\) |
| rad | \(3.060\times 10^{-1}\) | \(6.635\times 10^{-2}\) | \(4.613\) | \(5.07\times 10^{-6}\) |
| tax | \(-1.233\times 10^{-2}\) | \(3.760\times 10^{-3}\) | \(-3.280\) | \(0.001112\) |
| ptratio | \(-9.527\times 10^{-1}\) | \(1.308\times 10^{-1}\) | \(-7.283\) | \(1.31\times 10^{-12}\) |
| b | \(9.312\times 10^{-3}\) | \(2.686\times 10^{-3}\) | \(3.467\) | \(0.000573\) |
| lstat | \(-5.248\times 10^{-1}\) | \(5.072\times 10^{-2}\) | \(-10.347\) | \(<2\times 10^{-16}\) |
在 R 里，lm 可以算 \(\hat\beta\)，car 包的 hccm 可以算 \(\hat V_{\mathrm{ehw}}\) 及其变体。下面比较不同标准误给出的 \(t\) 统计量。这个例子里，有些回归系数的 EHW 标准误差得很远。
> library("car")
> ols.fit.hc0 = sqrt(diag(hccm(ols.fit, type = "hc0")))
> ols.fit.hc1 = sqrt(diag(hccm(ols.fit, type = "hc1")))
> ols.fit.hc2 = sqrt(diag(hccm(ols.fit, type = "hc2")))
> ols.fit.hc3 = sqrt(diag(hccm(ols.fit, type = "hc3")))
> tvalues = summary(ols.fit)$coef[,1] /
+ cbind(summary(ols.fit)$coef[,2],
+ ols.fit.hc0,
+ ols.fit.hc1,
+ ols.fit.hc2,
+ ols.fit.hc3)
> colnames(tvalues) = c("ols", "hc0", "hc1", "hc2", "hc3")
> round(tvalues, 2)| ols | hc0 | hc1 | hc2 | hc3 | |
|---|---|---|---|---|---|
| (Intercept) | \(7.14\) | \(4.62\) | \(4.56\) | \(4.48\) | \(4.33\) |
| crim | \(-3.29\) | \(-3.78\) | \(-3.73\) | \(-3.48\) | \(-3.17\) |
| zn | \(3.38\) | \(3.42\) | \(3.37\) | \(3.35\) | \(3.27\) |
| indus | \(0.33\) | \(0.41\) | \(0.41\) | \(0.41\) | \(0.40\) |
| chas1 | \(3.12\) | \(2.11\) | \(2.08\) | \(2.05\) | \(2.00\) |
| nox | \(-4.65\) | \(-4.76\) | \(-4.69\) | \(-4.64\) | \(-4.53\) |
| rm | \(9.12\) | \(4.57\) | \(4.51\) | \(4.43\) | \(4.28\) |
| age | \(0.05\) | \(0.04\) | \(0.04\) | \(0.04\) | \(0.04\) |
| dis | \(-7.40\) | \(-6.97\) | \(-6.87\) | \(-6.81\) | \(-6.66\) |
| rad | \(4.61\) | \(5.05\) | \(4.98\) | \(4.91\) | \(4.76\) |
| tax | \(-3.28\) | \(-4.65\) | \(-4.58\) | \(-4.54\) | \(-4.43\) |
| ptratio | \(-7.28\) | \(-8.23\) | \(-8.11\) | \(-8.06\) | \(-7.89\) |
| b | \(3.47\) | \(3.53\) | \(3.48\) | \(3.44\) | \(3.34\) |
| lstat | \(-10.35\) | \(-5.34\) | \(-5.27\) | \(-5.18\) | \(-5.01\) |
B.5 加权最小二乘
假定 \((w_i,x_i,y_i)\stackrel{\mathrm{IID}}{\sim}(w,x,y)\)，且 \(w\ne 0\)。在总体层面，可以定义加权最小二乘（weighted least squares, WLS）系数为
\[ \beta_w=\arg\min_{b}\mathrm{E}\bigl\{w(y-x^{\mathrm{T}}b)^2\bigr\}, \]
它满足
\[ \mathrm{E}\bigl\{wx(y-x^{\mathrm{T}}\beta_w)\bigr\}=0, \]
从而等于
\[ \beta_w=\bigl\{\mathrm{E}(wxx^{\mathrm{T}})\bigr\}^{-1}\mathrm{E}(wxy), \]
只要 \(\mathrm{E}(wxx^{\mathrm{T}})\) 可逆。
在样本层面，可以定义 WLS 系数为
\[ \hat\beta_w=\arg\min_{b}\sum_{i=1}^{n}w_i(y_i-x_i^{\mathrm{T}}b)^2, \]
它满足
\[ \sum_{i=1}^{n}w_ix_i(y_i-x_i^{\mathrm{T}}\hat\beta_w)=0, \]
从而等于
\[ \hat\beta_w = \left(n^{-1}\sum_{i=1}^{n}w_ix_ix_i^{\mathrm{T}}\right)^{-1} \left(n^{-1}\sum_{i=1}^{n}w_ix_iy_i\right), \]
只要 \(\sum_{i=1}^{n}w_ix_ix_i^{\mathrm{T}}\) 可逆。
在 R 里，可以在 lm 里指定 weights 来实现 WLS。
B.6 逻辑回归
B.6.1 模型
技术上，即便结果 \(y\) 是二值的，也可以用 OLS。可预测概率跑到 \([0,1]\) 外面，终究有点别扭。这动机了下面的模型：
\[ \operatorname{pr}(y_i=1\mid x_i)=g(x_i^{\mathrm{T}}\beta), \]
其中 \(g(\cdot):\mathbb{R}\to[0,1]\) 是单调函数，它的逆常常叫做连接函数（link function）。\(g(\cdot)\) 可以是任何一个随机变量的分布函数，但我们盯逻辑形式：
\[ g(z)=\frac{e^z}{1+e^z}=(1+e^{-z})^{-1}. \]
也可以把逻辑模型写成
\[ \operatorname{pr}(y_i=1\mid x_i)=\frac{e^{x_i^{\mathrm{T}}\beta}}{1+e^{x_i^{\mathrm{T}}\beta}}, \]
或者等价地，用定义 \(\operatorname{logit}(z)=\log\{z/(1-z)\}\)，有
\[ \operatorname{logit}\{\operatorname{pr}(y_i=1\mid x_i)\}=x_i^{\mathrm{T}}\beta. \]
假定 \(x_{i1}\) 是二值的。在逻辑模型下，
\[\begin{align*} \beta_1 &= \operatorname{logit}\{\operatorname{pr}(y_i=1\mid x_{i1}=1,\ldots)\} - \operatorname{logit}\{\operatorname{pr}(y_i=1\mid x_{i1}=0,\ldots)\} \\ &= \log \frac {\operatorname{pr}(y_i=1\mid x_{i1}=1,\ldots)/\operatorname{pr}(y_i=0\mid x_{i1}=1,\ldots)} {\operatorname{pr}(y_i=1\mid x_{i1}=0,\ldots)/\operatorname{pr}(y_i=0\mid x_{i1}=0,\ldots)}, \end{align*}\]
其中 \(\ldots\) 装着其余回归元 \(x_{i2},\ldots,x_{ip}\)。因此，系数 \(\beta_1\) 等于：条件于其他回归元时，\(x_{i1}\) 对 \(y_i\) 的优势比的对数。
B.6.2 最大似然估计
令 \(\operatorname{pr}(y_i=1\mid x_i)=\pi(x_i,\beta)\)。要估参数 \(\beta\)，可以最大化下面的似然函数：
\[\begin{align*} L(\beta) &= \prod_{i=1}^{n} \bigl\{\pi(x_i,\beta)\bigr\}^{y_i} \bigl\{1-\pi(x_i,\beta)\bigr\}^{1-y_i} \\ &= \prod_{i=1}^{n} \left\{\frac{\pi(x_i,\beta)}{1-\pi(x_i,\beta)}\right\}^{y_i} \bigl\{1-\pi(x_i,\beta)\bigr\} \\ &= \prod_{i=1}^{n} \bigl(e^{x_i^{\mathrm{T}}\beta}\bigr)^{y_i} \frac{1}{1+e^{x_i^{\mathrm{T}}\beta}} \\ &= \prod_{i=1}^{n} \frac{e^{y_i x_i^{\mathrm{T}}\beta}}{1+e^{x_i^{\mathrm{T}}\beta}}. \end{align*}\]
把最大化元记作 \(\hat\beta\)，叫做最大似然估计（maximum likelihood estimate, MLE）。对 \(L(\beta)\) 取对数再对 \(\beta\) 求导，可以证明 MLE 必须满足一阶条件：
\[ \sum_{i=1}^{n}x_i\bigl\{y_i-\pi(x_i,\hat\beta)\bigr\}=0. \]
若 \(x_i\) 含截距，MLE 还必须满足
\[ \sum_{i=1}^{n}\bigl\{y_i-\pi(x_i,\hat\beta)\bigr\}=0, \]
也就是说，观测到的 \(y_i\) 的平均，必须等于拟合概率 \(\pi(x_i,\hat\beta)\) 的平均。
用 MLE 的一般理论，可以证明它对真参数 \(\beta\) 相合，并且渐近正态：
\[ \sqrt{n}(\hat\beta-\beta)\to\mathrm{N}(0,V) \]
依分布，其中
\[ V = \Bigl(\mathrm{E}\bigl[\pi(x_i,\beta)\bigl\{1-\pi(x_i,\beta)\bigr\}xx^{\mathrm{T}}\bigr]\Bigr)^{-1}. \]
于是可以用
\[ \left( \sum_{i=1}^{n} \pi(x_i,\hat\beta)\bigl\{1-\pi(x_i,\hat\beta)\bigr\}x_ix_i^{\mathrm{T}} \right)^{-1} \]
去近似 \(\hat\beta\) 的协方差阵。在 R 里，glm 可以找到 MLE，并报告估计协方差阵。我们用 lalonde 数据说明逻辑回归：二值结果是 1978 年真实收入是否为正。
> library(Matching)
> data(lalonde)
> logit.re78 = glm(I(re78 > 0) ~ ., family = binomial,
+ data = lalonde)
> summary(logit.re78)偏差残差：Min \(-2.1789\)，\(1\)Q \(-1.3170\)，中位数 \(0.7568\)，\(3\)Q \(0.9413\)，Max \(1.0882\)。系数如下（原书完整 summary 输出）。
| Estimate | Std. Error | \(z\) | \(\operatorname{Pr}(>\lvert z\rvert)\) | |
|---|---|---|---|---|
| (Intercept) | \(1.910\) | \(1.241\) | \(1.539\) | \(0.1238\) |
| age | \(-2.812\times 10^{-3}\) | \(1.533\times 10^{-2}\) | \(-0.183\) | \(0.8545\) |
| educ | \(-2.179\times 10^{-2}\) | \(7.831\times 10^{-2}\) | \(-0.278\) | \(0.7808\) |
| black | \(-1.060\) | \(5.041\times 10^{-1}\) | \(-2.103\) | \(0.0354\) |
| hisp | \(2.741\times 10^{-1}\) | \(6.967\times 10^{-1}\) | \(0.393\) | \(0.6940\) |
| married | \(7.577\times 10^{-2}\) | \(3.057\times 10^{-1}\) | \(0.248\) | \(0.8042\) |
| nodegr | \(-1.984\times 10^{-1}\) | \(3.460\times 10^{-1}\) | \(-0.573\) | \(0.5664\) |
| re74 | \(7.857\times 10^{-6}\) | \(3.173\times 10^{-5}\) | \(0.248\) | \(0.8044\) |
| re75 | \(4.016\times 10^{-5}\) | \(6.058\times 10^{-5}\) | \(0.663\) | \(0.5074\) |
| u74 | \(-6.177\times 10^{-2}\) | \(4.095\times 10^{-1}\) | \(-0.151\) | \(0.8801\) |
| u75 | \(1.505\times 10^{-2}\) | \(3.518\times 10^{-1}\) | \(0.043\) | \(0.9659\) |
| treat | \(5.412\times 10^{-1}\) | \(2.222\times 10^{-1}\) | \(2.435\) | \(0.0149\) |
B.6.3 伸到病例对照研究
病例对照研究（case-control study）里，抽样是条件于二值结果的：结果 \(y_i=1\) 与 \(y_i=0\) 的单元，被抽中的概率不同。令 \(s_i\) 为抽样指示。病例对照研究里，
\[ \operatorname{pr}(s_i=1\mid x_i,y_i)=\operatorname{pr}(s_i=1\mid y_i) \]
是 \(y_i\) 的函数，我们只观测到 \(s_i=1\) 的那些单元。
就模型和病例对照的抽样机制而言，还能不能用逻辑回归去估系数，看起来并不显然。可 Prentice and Pyke (1979) 证明了一个正面结果。病例对照研究里，逻辑回归仍然能一致地估出除截距以外的所有系数。
B.6.4 带权重的逻辑回归
有时单元 \(i\) 有权重 \(w_i\)。于是可以通过解
\[ \sum_{i=1}^{n}w_ix_i\bigl\{y_i-\pi(x_i,\hat\beta)\bigr\}=0 \]
来拟合加权逻辑回归。在 R 里，可以在 glm 里指定 weights 来实现加权逻辑回归。
B.7 习题
B.1 带截距的样本 WLS
假定回归元 \(x_i\) 含截距。证明
\[ \bar y_w=\bar x_w^{\mathrm{T}}\hat\beta_w \tag{B.5} \]
其中 \(\bar x_w=\sum_{i=1}^{n}w_ix_i/\sum_{i=1}^{n}w_i\)、\(\bar y_w=\sum_{i=1}^{n}w_iy_i/\sum_{i=1}^{n}w_i\) 是 \(x_i\) 与 \(y_i\) 的加权平均。
B.2 二值回归元的总体 OLS
假定 \(x\) 是二值的。定义总体 OLS：
\[ (\alpha,\beta)=\arg\min_{(a,b)}\mathrm{E}\{(y-a-bx)^2\}. \]
证明 \(\beta=\mathrm{E}(y\mid x=1)-\mathrm{E}(y\mid x=0)\)，以及 \(\alpha=\mathrm{E}(y\mid x=0)\)。
B.3 一元 WLS
作为 WLS 的特例，定义
\[ (\hat\alpha_w,\hat\beta_w) = \arg\min_{(a,b)} \sum_{i=1}^{n}w_i(y_i-a-bx_i)^2, \]
其中 \(w_i\ge 0\)。证明
\[ \hat\beta_w = \frac{\sum_{i=1}^{n}w_i(x_i-\bar x_w)(y_i-\bar y_w)}{\sum_{i=1}^{n}w_i(x_i-\bar x_w)^2} \tag{B.6} \]
以及
\[ \hat\alpha_w=\bar y_w-\hat\beta_w\bar x_w, \tag{B.7} \]
其中 \(\bar x_w=\sum_{i=1}^{n}w_ix_i/\sum_{i=1}^{n}w_i\)、\(\bar y_w=\sum_{i=1}^{n}w_iy_i/\sum_{i=1}^{n}w_i\) 是 \(x_i\) 与 \(y_i\) 的加权平均。
进一步假定 \(x_i\) 是二值的。证明
\[ \hat\beta_w = \frac{\sum_{i=1}^{n}w_ix_iy_i}{\sum_{i=1}^{n}w_ix_i} - \frac{\sum_{i=1}^{n}w_i(1-x_i)y_i}{\sum_{i=1}^{n}w_i(1-x_i)}. \tag{B.8} \]
也就是说：一元 WLS 里回归元是二值时，回归元的系数等于加权均值之差。
注： 证明 (B.8) 时，对 WLS 问题做一次合适的再参数化。否则推导会很烦。
B.4 正交回归元的 OLS
考虑把 \(n\) 维向量 \(Y\) 对 \(n\times p\) 矩阵 \(X\) 做样本 OLS，系数是 \(\hat\beta\)。把 \(X\) 拆成 \(X=(X_1,X_2)\)，其中 \(X_1\) 是 \(n\times k\) 矩阵，\(X_2\) 是 \(n\times l\) 矩阵，\(p=k+l\)。相应地，把 \(\hat\beta\) 拆成
\[ \hat\beta=\begin{pmatrix}\hat\beta_1\\ \hat\beta_2\end{pmatrix}. \]
假定 \(X_1\) 与 \(X_2\) 正交，也就是 \(X_1^{\mathrm{T}}X_2=0\)。证明：\(\hat\beta_1\) 等于把 \(Y\) 对 \(X_1\) 做 OLS 的系数，\(\hat\beta_2\) 等于把 \(Y\) 对 \(X_2\) 做 OLS 的系数。
B.5 回归元做非退化变换后的 OLS
令 \(\hat\beta\) 为把 \(n\) 维向量 \(Y\) 对 \(n\times p\) 矩阵 \(X\) 做样本 OLS 的系数。令 \(\Gamma\) 为 \(p\times p\) 非退化矩阵，定义 \(X'=X\Gamma\)。令 \(\hat\beta'\) 为把 \(Y\) 对 \(X'\) 做样本 OLS 的系数。
证明
\[ \hat\beta=\Gamma\hat\beta'. \]
B.6 OLS 估计量的方差
假定 \((x_i,y_i)_{i=1}^{n}\stackrel{\mathrm{IID}}{\sim}(x,y)\)，并且
\[ \mathrm{E}(y\mid x)=x^{\mathrm{T}}\beta, \qquad \operatorname{var}(y\mid x)=\sigma^2(x). \]
证明 OLS 估计量 \(\hat\beta\) 的条件均值是
\[ \mathrm{E}(\hat\beta\mid x_1,\ldots,x_n)=\beta, \]
条件方差是
\[ \operatorname{var}(\hat\beta\mid x_1,\ldots,x_n) = \left(\sum_{i=1}^{n}x_ix_i^{\mathrm{T}}\right)^{-1} \left(\sum_{i=1}^{n}\sigma^2(x_i)x_ix_i^{\mathrm{T}}\right) \left(\sum_{i=1}^{n}x_ix_i^{\mathrm{T}}\right)^{-1}. \]
注： 这道题也许能给方差估计 (B.3) 与 (B.4) 一点直觉。