统计Ⅱ-习题记录

以下是对于一些习题的记录。

一些常用结论

这里记录的是在解题过程中的一些关键步骤，或是一些有用的结论。

\(SSR=\hat\beta^2S_{XX}={S^2_{XY}\over S_{XX}}\)

\(\beta={S_{XY}\over S_{XX}}\)

\(Corr^2(Y,\hat Y)=R^2\)

\((X'X+\lambda I)^{-1}X=X(X'X+\lambda I)^{-1}\)

一个常用的定理

定理 1：对于 n 个独立同分布的正态观测 \(X_1,...,X_n\overset{iid}{\sim}N(\mu, \sigma^2)\) ，我们有：

• \(\bar{X}\sim N(\mu, {\sigma^2\over n})\)
• \(s^2={1\over n-1}\sum(X_i-\bar X)^2\sim {1\over n-1}\sigma^2\chi^2_{n-1}\)
• \(s^2\bot \bar X\)

---

e.g. 假设检验

Consider a multiple linear regression model $Y = β_0 1_n + Xβ + \epsilon, E(\epsilon) = 0, cov(\epsilon) = σ^2 I $, where \(Y\) is the \(n × 1\) observation vector and \(X\) is the \(n × k\) design matrix with full rank. Find the F-statistics for the following two hypothesis testing problems:

(1) \(H_{01} : β_1 = ... = β_k = c\)

(2) \(H_{02} : β_1 = ... = β_k\)

Here, c is some given (known) constant.

一道经典的假设检验问题，我会先给出解答，之后谈谈 Reduced Model 与假设检验之间的关系。

• 在（1）下，\(β_1 = ... = β_k\) 给定了，这时，我们不能把\(X_i\beta_i\) 看成是一个随机变量，我们尝试把这些项排除出去。先来考虑\(c=0\) 情况：

这时，显然问题就变成了验证线性模型是否存在的情况，可以使用 F 检验；为了和后面的内容呼应，在这里先加以补充说明：来看\(H_0\)，若将其写成矩阵形式 \(C\beta=0\)，则这里的 C 是一个 \(k\times k\)的矩阵。这时，模型退化为

\[Y=\beta_01_n+\epsilon \]

再来看 F 检验的具体形式：

\[F={MSR\over MSE}={SSR/k\over SSE/n-k-1}\overset{H_0}{\sim}F_{k,n-k-1} \]

（可以看到，这个 F 分布的自由度和在 Reduced Model 中的情况是一致的；此前并没有将两者联系起来，现在看来是有关系的。）

再来看 \(c\ne0\) 的情况，这时的模型为

\[Y=\beta_01_n+X1_nc+\epsilon\\\tilde Y=Y-c X1_n=\beta_01_n+\epsilon\\ \]

我们做了简单的变形，退化到了\(c=0\) 时的情况。

• 在（2）中，没有给定具体的值，也就是说还带有随机性，不能像（1）中那样来处理了。我们令\(\beta_i=r\)，则模型可表为

\[Y=\beta_01_n+r(X1_n)+\epsilon\\=\beta_01_n+rZ+\epsilon \]

我们可以把 \(X'1_n\) 看成单独的一个新的变量\(Z\)。在新的变量下，我们的 design 变为 \(\tilde X=(1_n,X1_n)\)。在这种 Reduced Model 下的 SSE 为

\[SSE_R=Y'(\tilde H-{1\over n}11')Y, \tilde H=\tilde X(\tilde X'\tilde X)^{-1}\tilde X'\tag{1} \]

并且，我们不加证明地指出：

\[SSE_R-SSR_F\sim \chi^2_{k-1} \]

这里用到了一个重要的结论，即对于 Reduced Model 来说，假设约束 \(C\beta=0\) 中约束的个数（C 的行数）为 d，则 \(SSE_R-SSR_F\) 服从自由度为 d 的卡方分布；并且其与 \(SSR_F\) 是独立的。在（2）中，零假设可表为 \(k-1\) 个方程，因此得到上式。这样的话，我们容易得出检验统计量

\[F={SSE_R-SSR_F/k-1\over SSE_F/n-k-1}\overset{H_0}{\sim}F_{k-1,n-k-1} \]

至此，解题结束。

老师课上只是讲了模型的假设检验，之后讲了 Reduced Model，也介绍了如何对 Reduced Model 进行检验；然而由于这门课没有教材，一直没搞清楚两者的关系，基于这道题似乎看到了两者的关系：对于模型的检验实际上可以看做是对 Reduced Model 检验的一种特殊情况。

一直以为， Reduced Model 是在一定的约束之下去求估计，这样的理解当然没有问题；那么我们再往前想一步呢？这种线性约束是哪里来的？或者说，对于它的检验是怎么做的？

想到这一点，就容易理解了：我们的假设检验的过程，实际上就是在一个 Reduced Model 下进行的；只是，这时我们的约束比较特殊（对于模型的检验假定了所有的系数均为零，对于单个系数的检验假定了某个系数为零），在这样的约束之下（Reduced Model 之下）我们推导其分布并给出检验统计量。这种情况下，我们感兴趣的量是 \(SSR\)（因为它描述了模型能假设的变异的程度）并根据其分布给出了 F 统计量。我们进一步来看，

\[SSR=SST-SSE\\=SSE_R-SSE_F \]

第一行是我们所熟知的分解；然而，我们来看 \(SST\)，在零假设之下（所有变量均与因变量无关），我们用样本均值来估计所有的点，这不正是 \(SST\) ？上面已经提到了，\(H_0\) 可以表示成一个 \(k\times k\)的矩阵，所以 \(SSR\) ，或者说 \(SSE\) 之差服从自由度为 k 的卡方分布。

从而说明了，对于模型的检验可以理解成对一个 Reduced Model 检验的过程；进一步，我们对任意的一个线性的假设检验，均可以从 Reduced Model 出发，逐步构建统计量。1. 一种方法是将约束表示成矩阵的形式，使用已有的公式直接暴力求解；2. 实际上，问题的难点就在于给出\(SSE_R\) 的表达形式，对于一些特殊的约束，我们可以通过对模型的变换或化简更为容易得得到其具体的形式，如在问题（1）中，我们变换成为对于模型的 F 检验，而在问题（2）中，我们通过变量代换最终给出了式（1），计算要比直接套 Reduced Model 的长串公式简单很多。

e.g. 误差相关下的 BLUE

Consider the linear regression model, \(Y = Xβ + ε\), where \(Eε = 0\) and \(cov(ε) = σ^2 Σ\), in which \(Σ\) is a known positive definite matrix. Find the BLUE for \(β\) and derive its variance-covariance matrix.

这是对于原来假定误差之间独立的线性回归的一个拓展。在这种情况下，我们可以对 Y 作变换，我们知道 \(Y\sim N(X\beta,\sigma^2\Sigma)\)，则有

\[\Sigma^{-1/2} Y\sim N(\Sigma^{-1/2}X\beta,\sigma^2I)\\ \widetilde Y\sim N(\tilde X\beta,\sigma^2I) \]

于是，我们代入原有 LR 中的结论

\[\hat\beta=(\tilde X'\tilde X)^{-1}\tilde X\tilde Y=(X'\Sigma^{-1}X)^{-1}X'\Sigma^{-1}XY \]

并不加证明地指出，此 LSE 即为 BLUE（事实上，可用在误差不相关的情况下同样的思路进行证明，参考 wiki）。并有

\[Cov(\hat\beta)=\sigma^2(X'\Sigma X)^{-1} \]

e.g. 构建 CI

Let \(Y_1^* , · · · , Y_m^*\) be independent responses of m future (not yet observed) items such that

\[Y^*= (Y_1^* , · · · , Y_m^*)'\sim N_m(X^*\beta, \sigma^2I_m) \]

for a given \(m × d\) design matrix \(X^∗\) . Let \(a = (a_1 ,... , a_m )'\) be a vector of \(m\) known constants. Suppose we are interested in predicting

\[l=\sum_{j=1}^m a_jY_j^*=a'Y^* \]

Give the \(100(1 − α)%\) prediction interval for \(l\).

先来看看我们要求的量的分布如何：

\[l=a'Y^*\sim N(a'X^*\beta, \sigma^2a'a) \]

服从一个正态分布，那么很显然，我们对于\(l\) 的估计可采用其均值；注意，我们需要估计其中的参数，因此

\[\hat l=a'X^*\hat\beta \]

另外，此题是对于一个未知量进行 predict，其本身就带有随机性（区别于对 mean response 进行 predict），因此，我们需要考虑的变量是 \(\hat l-l\)，显然有

\[E[\hat l-l]=E[a'X^*\hat\beta-a'X^*\beta-a'\epsilon]=0\\ Var(\hat l-l)=Var(a'X^*\hat\beta-a'X^*\beta-a'\epsilon)=Var(a'X^*(\hat\beta-\beta))+Var(a'\epsilon)\\ =a'X^*Var(\hat\beta)X^{*'}a+a'Var(\epsilon)a=\sigma^2[a'X^*(X^{*'}X^*)^{-1}X^{*'}a+a'a] \]

这样，我们就可以构建检验统计量

\[{(\hat l-l)/\sigma\sqrt{a'X^*(X^{*'}X^*)^{-1}X^{*'}a+a'a}\over s/\sigma} \]

分子服从正态分布，分母服从 \(\sqrt{\chi^2_{m-d}/m-d}\)，因此，总体服从 \(t_{m-d}\) 分布。这里要注意自由度，从 design 来看，列数为 \(d\)，可知自由参数为\(d-1\)个。为说明清楚，这里再赘述一下 \(s^2=\hat\sigma^2={1\over m-d}SSE\) 。再基于此构建 CI。

总结：这道题的表述看似有点复杂，实际的意思就是说，我们计划了 m 个点 \(x_i\)，然后要从这些点去 predict 各 \(y_i\)，而最后要得到的结果是这整个预测向量Y 的一个线性组合。可以看到，这里的流程其实和求只给定了一个点 x，去求y 的 CI 的流程是一样的：1. 先分析待求量 \(G\) 的分布（一般为正态）；2. 构建预测表达式\(\hat G\)（包含参数 \(\hat\beta\)）；3. 求出 \(\hat G-G\) 的分布（因为 \(G\) 本身可能是随机变量；这里就是求 mean 和 var）；4. 使用 t-test，注意把其中的方差替换为其预测量 \(s^2\)（注意其自由度）。

e.g. 证明：约定系数等于相关系数的平方

For a multiple linear regression model, prove that the coefficient of multiple determination is the square of the sample correlation between \(y_i\) and \(\hat y_i\) .

一个重要的结论：**多元回归的决定系数 \(R^2\) 等于观测值和预测值之间的样本相关系的平方。

我们来看 corr 是什么形式，因为这两个 sample 的均值均为\(\bar Y\)，因此该 corr 可表为

\[Corr(Y,\hat Y)={\sum(y_i-\bar y)(\hat y-\bar y)\over \sqrt{\sum(y_i-\bar y)^2\sum(\hat y-\bar y)^2}}={\sum(y_i-\bar y)(\hat y-\bar y)\over \sqrt{SST SSR}} \]

可以看到已经比较相似了，相较于 \(R^2={SSR\over SST}\)，为证

\[Corr^2(Y,\hat Y)=R^2 \]

只需要证明

\[\sum(y_i-\bar y)(\hat y-\bar y)=SSR=\sum(y_i-\bar y)^2\tag{1}\\ \]

因此我们可以对左边进行添项展开

\[\sum(y_i-\bar y)(\hat y-\bar y)=\sum(y_i-\hat y+\hat y-\bar y)(\hat y-\bar y)=SSR+\sum(y_i-\hat y)(\hat y-\bar y) \]

问题即变为证明 \(\sum(y_i-\hat y)(\hat y-\bar y)=0\)，为此，我们可以考虑其矩阵表示

\[\sum(y_i-\hat y)(\hat y-\bar y)=(Y-HY)'(HY-{1\over n}11'Y)=Y'(I-H)(H-{1\over n}11')Y\\ =Y'(I-H)HY+Y'(I-H){1\over n}11'Y\tag{2} \]

对于前项显然为 0；对于后项来说，注意 H 为 X 所对应的表示矩阵，而 X 的第一列即为 1，因此\(H1=1\)，可见后项也为 0。证毕。

注：

1. 对于（1）式有更简单的方法：

\[\sum(y_i-\bar y)(\hat y-\bar y)=(Y-HY)'(HY-{1\over n}11'Y)=Y'(I-H)(H-{1\over n}11')Y\\ =Y'(I-H)(-{1\over n}11')Y=Y'(H-{1\over n}11')Y \]

第二行中两个等式成立的原因和（2）式中一样。

2. 另，对于 \(\bar{\hat Y}=\bar{Y}\) 也简要说明：表示成矩阵形式

\[\overline{\hat Y}=1'HY=1'Y \]

这里再次用到了这条等式：

\[H1=1\tag{3} \]

e.g. \(AIC\) 和\(C_p\) 的等价性

证明：变量选择中，若已知方差\(\sigma^2\) ，证明 max AIC 和 min \(C_p\) 是等价的

我们先来回忆两者的定义

\[AIC(p)=\log(L)-p\\ C_p={SSE\over \sigma^2}-[n-2(p+1)] \]

AIC 实际上是对数似然加了一个惩罚项 \(p\)，其越大越好；而\(C_p\) 表征的则是「错误率」，越小越好。

我们把 AIC 进行变换

\[AIC(p)=-{n\over 2}\log(2\pi\sigma^2)-{1\over 2\sigma^2}\sum(y_i-x_i'\hat\beta)-p\\ =-{1\over 2\sigma^2}(Y-X\hat\beta)'(Y-X\hat\beta)-p+C=-{1\over 2}{1\over \sigma^2}SSE-p+C\\ =-{1\over 2}[{SSE\over \sigma^2}-[n-2(p+1)]]-{1\over2}[n-2(p+1)]-p+C=-{1\over 2}C_p-{1\over n}+1+C \]

可见

\[\arg\max_p AIC(p)=\arg\max_p-{1\over 2}C_p-{1\over n}+1+C=\arg\min_p{1\over 2}C_p=\arg\min_pC_p \]

e.g. 一个技巧

Show that the matrix \(H(λ) = X(X^T X + λI)^{−1} X^T\) , associated with ridge regression, is not a projection matrix (for any \(λ > 0\))

这里并不是想要解这道题，只是这里有个很妙的想法：因为有

\[X(X'X+\lambda I)=(XX'+\lambda I)X \]

两边各乘一个逆，可得

\[(X'X+\lambda I)^{-1}X=X(X'X+\lambda I)^{-1}\tag{1} \]

得到了一个非常优美的恒等式。

由此，我们可以得到

\[H(λ) = X(X' X + λI)^{−1} X'=(X' X + λI)^{−1}XX' \]

\[H^2(λ) =(X' X + λI)^{−1}XX'(X' X + λI)^{−1}XX' \]

由于 \(λ > 0\)，因此 \(X' X + λI\ne X'X\) ，所以 \(H^2(λ)\ne H(λ)\)，即 \(H(λ)\) 不是投影阵。

e.g. ridge

Recall that there exists \(λ >\) 0 such that \(MSE(\hat β) > MSE(\hat β(λ))\). Verify that this carries over to the linear predictor. That is, there exists a \(λ > 0\) such that

\[MSE(\hat Y ) = MSE(X \hat β) > MSE(X \hat β(λ)). \]

这里需要事先说明几点：这题是针对 ridge 的，其中

\[\hat \beta(\lambda)=(X'X+\lambda I)^{-1}X'Y\\ =(X'X+\lambda I)^{-1}X'X\hat\beta\tag{1} \]

第一行为 LSE，第二行为其与无 ridge 的估计量的关系。注意该估计有偏

\[bias(\hat \beta(\lambda))=E[\hat \beta(\lambda)-\beta]=-\lambda(X'X+\lambda I)^{-1}\beta\tag{2} \]

代入（1）式第二行，并利用\(\hat\beta\)的无偏性即可。

另外，我们在这里定义的 \(MSE\) 也要说明一下，注意它是对于一个向量而言的，我们定义为

\[MSE(\beta)=E||\hat\beta-\beta||^2=E||\hat\beta-E\hat\beta+E\hat\beta-\beta||^2\\ =E||\hat\beta-E\hat\beta||^2+||E\hat\beta-\beta||^2\\ =Var(\hat\beta)+bias^2(\hat\beta) \]

注意，这里的 \(Var(\hat\beta)\) 不同于我们日常用的 \(Cov(\hat\beta)\)（尽管我日常也会把后者写成 var 的形式），我们在这里将两者进行区别，var 表示一个随机向量与其均值的距离平方的期望（是一个数字）；而将 cor 表为随机向量的协方差矩阵。我们也用平方来简化表示 bias 的内积。

\[Var(\hat\beta)=E(\hat\beta-E\hat\beta)'(\hat\beta-E\hat\beta)=Tr E(\hat\beta-E\hat\beta)(\hat\beta-E\hat\beta)'=Tr Cor(\hat\beta)\tag{3} \]

对于\(MSE(\hat β) > MSE(\hat β(λ))\) 的证明是课上所讲，在此不给出了，其过程和下面是完全一样的。说明完毕，下面开始证明

\[MSE(\hat Y ) = MSE(X \hat β)=Var(X\hat\beta)+bias^2(X\hat\beta)\\ =Tr(Cor(X\hat\beta))=Tr(X\sigma^2(X'X)^{-1}X')=p\sigma^2 \]

第二行第一个等式用到了\(\hat\beta\) 的无偏性，以及关系（3）。

\[MSE(X \hat β(λ))=Var(X\hat\beta(λ))+bias^2(X\hat\beta(λ))\\ =Tr(XCor(\hat\beta(\lambda))X')+(\lambda X(X'X+\lambda I)^{-1}\beta)'(\lambda X(X'X+\lambda I)^{-1}\beta)\\ =\sigma^2 Tr(X'X+\lambda I)^{-1}X'X(X'X+\lambda I)^{-1}X'X+ \lambda^2 \beta'(X'X+\lambda I)^{-1}X'X(X'X+\lambda I)^{-1}\beta\\ =Tr(X'X+\lambda I)^{-1}X'X(X'X+\lambda I)^{-1}[\sigma^2X'X+\lambda^2\beta\beta']\overset{\triangle}{=}f(\lambda) \]

在第二行中，我们将 var 转为 cor 并把式（2） bias 代入，第三行中我们将 \(Cor(\hat\beta(\lambda))\) （容易得到）代入；最终得到了一个\(\lambda\) 的函数，我们对其求导

\[f'(\lambda)=-Tr(X'X+\lambda I)^{-2}X'X(X'X+\lambda I)^{-1}[\sigma^2X'X+\lambda^2\beta\beta']\\-Tr(X'X+\lambda I)^{-1}X'X(X'X+\lambda I)^{-2}[\sigma^2X'X+\lambda^2\beta\beta']\\+2\lambda Tr(X'X+\lambda I)^{-1}X'X(X'X+\lambda I)^{-1}\beta\beta'\\ \]

矩阵求导不太会，不知道有没有算错 Orz。Anyway，我们容易注意到

\[f(0)=MSE(\hat\beta) \]

而在 0 点，\(f(\lambda)\) 的导数值为

\[f'(\lambda)|_{\lambda=0}=-2Tr(X'X+\lambda I)^{-2}(\sigma^2X'X)<0 \]

可见，总存在一个大于 0 的值，使得

\[MSE(X \hat β) > MSE(X \hat β(λ)) \]

小结：这题有点偏，结合了很多课上的结论，考试应该不会涉及；不过这里的思路还是很清楚的，里面的技巧也很有用，我把它写下来算是对矩阵的相关运算进行复习。

e.g. lasso 的唯一性

Lasso regression fits the linear regression model \(Y = Xβ + \epsilon\) using the penalty with 1−norm, that is, the lasso regression loss function is given by

\[L(β; λ_1 ) = ||Y − Xβ||^2_2 + λ_1 ||β||_1 \]

and

\[\hat β(λ_1 ) = \arg \min L(β; λ_1 ) \]

Show that its linear predictor \(X β(λ_1 )\) is unique.

这里提到了 lasso，要注意 lasso 可能没有唯一解；而这里要求证明虽然其解可能不唯一，但对于不同的估计来说其预测 \(X β(λ_1 )\) 却是唯一的。

我们先假定问题有两个不同的最优解 \(\beta_1, \beta_2\)

\[L(\beta_1,\lambda)=L(\beta_2,\lambda)=c^* \]

反证法，若假设有

\[X\beta_1\ne X\beta_2 \]

则考虑到二次函数\(||Y-a||^2\)的强凸性，以及 1 范数的凸性，则对于线性组合 \(\tilde\beta=\alpha\beta_1+(1-\alpha)\beta_2\)

\[L(\tilde\beta\lambda)=L(\alpha\beta_1+(1-\alpha)\beta_2;\lambda)<\alpha L(\beta_1;\lambda)+(1-\alpha) L(\beta_1;\lambda) \]

也就是说， \(\beta_1, \beta_2\) 并非最优解，矛盾。

注：对于 lasso 来说，我们说明其唯一解的条件。

我们记

\[sign(\beta)=(sing(\beta_0),...)' \]

于是有

\[L(β; λ_1 ) = ||Y − Xβ||^2_2 + λ_1 sign(\beta)'β\\ \frac{\part L(β; λ_1 )}{\part \beta}=2X'(X\beta-Y)+λ_1 sign(\beta)\\ \frac{\part L(β; λ_1 )}{\part \beta\part \beta'}=2X'X \]

可见，只有当\(X'X\) 正定，也即 X 满秩 \(rank(X)=p\) 的情况下，loss 是强凸的，有唯一解。

e.g. RL 视角下的 ANOVA

Consider a 1-way ANOVA model. By rewriting the model using the regression approach, prove that \(S_e /σ^2\) (or \(SSE/σ^2\) )\(\sim χ_{n−r}^2\) and is independent of \(\bar{ Y_{i.}}\) .

写出 ANOVA 表示成线性回归的形式，然后直接利用 LR 中的方式来证明卡方分布。即证明 \(I-H\) 对称幂等，因此特征值均为1/0，并有\(Tr(I-H)=n-r\)，最终利用特征分解的形式得出结论。

为证明独立性，可将两者表示为：

\[SSE=\sum_{ij}(Y_{ij}-\overline Y_{i.})^2=\sum_{ij}(\epsilon_{ij}-\bar \epsilon_{i.})^2\\ \overline Y_{i.}=\mu_i+\bar\epsilon_{i.} \]

即要证明\(\bar\epsilon_{i.}\)与\(\sum_{j}(\epsilon_{ij}-\bar \epsilon_{i.})^2\) 之间独立（\(\bar\epsilon_{i.}\) 与其他项显然独立）。两者分别为\(\{\epsilon_{i1},...,\epsilon_{in_i}\}\) 这些独立正态观测的方差和均值，由定理1 可知其独立。