多元统计分析大纲（适用于一轮复习）

本文基于《应用多元统计分析》教材整理，适用于备考浙江大学多元统计分析（3学分，课程号06123040），并且跳过其中回归分析的部分。基于不同教材的多元统计侧重点不同，考点也不同，请勿完全参考。整理过程多有疏漏，如果错误，欢迎在评论区指出。

大纲如下：

1. 多元正态分布
   • 多元正态分布基础
   • 参数估计与常用统计量
   • 常用分布
   • 假设检验
2. 归类
   • 判别分析（Discriminatory Analysis）
   • 聚类分析（Cluster Analysis）
3. 降维
   • 主成分分析（Principal Component Analysis）
   • 因子分析（Factor Analysis）
   • 典型相关分析（Canonical Correlation Analysis）
4. 典型计算（另开一篇）

由于是复习笔记，所以内容较为简略。

目录

• Part 1：多元正态分布
  • 1-1 多元正态分布基础
  • 1-2 参数估计与常用统计量
  • 1-3 常用分布
  • 1-4 均值假设检验
  • 1-5 其他假设检验
• Part 2：归类
  • 2-1 距离判别
  • 2-2 贝叶斯判别
  • 2-3 费希尔判别
  • 2-4 系统聚类
• Part 3：降维
  • 3-1 主成分分析
  • 3-2 主成分分析实操
  • 3-3 因子模型
  • 3-4 因子模型实操
  • 3-5 典型相关变量求解
  • 3-6 典型相关实操

Part 1：多元正态分布

1-1 多元正态分布基础

样本数据阵的排布，通常是将每一个样本视为一行，\(n\)个样本得到的样本数据阵应当是\(n\times p\)的，\(p\)是所感兴趣的维度。

接下来假设随机向量是\(p\)维的，总体是\(X=(X_1,X_2,\cdots,X_p)\)。随机向量具有如下的基本概念：

1. 联合分布函数：\(p\)元函数\(F\)，定义为

   \[F(x_1,\cdots,x_p)=\mathbb{P}(X_1\le x_1,\cdots,X_p\le x_p). \]

2. 联合密度函数：\(p\)元函数\(f\)，如果某个联合分布可以表示为

   \[F(x_1,\cdots,x_p)=\int_{-\infty}^{x_1}\cdots\int_{-\infty}^{x_p} f(u_1,\cdots,u_p){\rm d}u_1\cdots{\rm d}u_p. \]

   满足非负性与规范性。

3. 边缘分布：随机向量\(X\)的部分维度\((X_{i_1},\cdots,X_{i_m})\)的分布函数。类似定义边缘密度。

4. 条件分布：将\(X\)分为\(r\)维分量\(X^{(1)}\)与\(p-r\)维分量\(X^{(2)}\)，在给定另一个分量的取值时，分量的分布。类似定义条件密度。

5. 特征函数：是\(p\)元函数，定义为

   \[\Phi(t)=\mathbb{E}(e^{{\rm i}t'X}). \]

   类比一元情况下的特征函数。

在边缘分布下可定义分量的独立性：联合分布函数等于边缘分布函数的乘积。相互独立的分量之间，条件分布等于边缘分布。

设\(X=(X_1,\cdots,X_p)\)是\(p\)维随机向量，\(Y=(Y_1,\cdots,Y_q)\)是\(q\)维随机向量。随机向量具有如下常用的数字特征：

1. 均值向量，需要每一个\(X_i\)都有\(\mathbb{E}(X_i)=\mu_i\)存在。

   \[\mathbb{E}(X)=\begin{bmatrix} \mathbb{E}(X_1) \\ \vdots \\ \mathbb{E}(X_p) \end{bmatrix}=\begin{bmatrix} \mu_1 \\ \vdots \\ \mu_p \end{bmatrix}\xlongequal{def}\mu. \]

2. 自协方差阵，需要每一个\(X_i, X_j\)的协方差\(\mathbb{Cov}(X_i, X_j)\)存在。

   \[\mathbb{D}(X)=\mathbb{E}[(X-\mathbb{E}(X))(X-\mathbb{E}(X))']\xlongequal{def}(\sigma_{ij})_{p\times p}, \\ \mathbb{D}(X)=\begin{bmatrix} \mathbb{Cov}(X_1, X_1) & \mathbb{Cov}(X_1, X_2) & \cdots & \mathbb{Cov}(X_1, X_p) \\ \mathbb{Cov}(X_2, X_1) & \mathbb{Cov}(X_2, X_2) & \cdots & \mathbb{Cov}(X_2, X_p) \\ \vdots & \vdots & \ddots & \vdots \\ \mathbb{Cov}(X_p, X_1) & \mathbb{Cov}(X_p, X_2) & \cdots & \mathbb{Cov}(X_p, X_p) \end{bmatrix}\xlongequal{def}\Sigma. \]

3. 互协方差阵，需要每一个\(X_i, Y_j\)的协方差\(\mathbb{Cov}(X_i, Y_j)\)存在。

   \[\mathbb{COV}(X, Y)=\begin{bmatrix} \mathbb{Cov}(X_1, Y_1) & \mathbb{Cov}(X_2, Y_1) & \cdots & \mathbb{Cov}(X_1, Y_q) \\ \mathbb{Cov}(X_2, Y_1) & \mathbb{Cov}(X_2, Y_2) & \cdots & \mathbb{Cov}(X_2, Y_q) \\ \vdots & \vdots & \ddots & \vdots \\ \mathbb{Cov}(X_p, Y_1) & \mathbb{Cov}(X_p, Y_2) & \cdots & \mathbb{Cov}(X_p, Y_q) \end{bmatrix}_{p\times q}. \]

4. 相关系数阵，需要自协方差阵存在，即相关系数构成的矩阵。

   \[r_{ij}=\frac{\sigma_{ij}}{\sqrt{\sigma_{ii}\sqrt{\sigma_{jj}}}},R\xlongequal{def}(r_{ij})_{p\times p}. \]

5. 标准差阵，每一个随机变量的标准差构成的对角阵，极少用到。

   \[V^{1/2} ={\rm diag}(\sqrt{\sigma_{11}}, \cdots, \sqrt{\sigma_{pp}}). \]

   有\(\Sigma=V^{1/2}RV^{1/2}\)和\(R=V^{-1/2}\Sigma V^{-1/2}\)。

数字特征的相关性质：

1. 计算性质，有

   \[\mathbb{E}(AXB)=A\mathbb{E}(X)B,\\ \mathbb{D}(AX)=A\mathbb{D}(X)A',\\ \mathbb{COV}(AX, BY)=A\mathbb{COV}(X, Y)B'. \]

2. \(\Sigma=\mathbb{D}(X)\)是对称非负定的，且当\(X\)线性无关时是正定的。对称非负定矩阵具有平方根，即\(\exists L\)使得\(\Sigma=L^2\)。

   \[\Sigma=U\Lambda U',\quad L=U\Lambda^{1/2}U'. \]

3. 如果\(\Sigma\)追加是正定的，则有Cholesky分解：\(\Sigma=LL'\)，这里\(L\)是下三角矩阵。

多元正态分布的定义方式：

1. 独立标准正态分布的线性组合：设\(U=(U_1,\cdots,U_q)'\)，每一个\(U_q\)独立同分布服从\(N(0, 1)\)。任意\(p\times q\)矩阵\(A\)和\(p\)维常数列向量\(\mu\)构造的随机向量：

   \[AU+\mu\xlongequal{def}X\sim N_p(\mu, AA'). \]

2. 特征函数法（几乎不用）：对于\(\Sigma\ge 0\)和常数向量\(\mu\)，多元正态分布随机向量\(X\)具有如下的特征函数：

   \[\Phi_X(t)=\exp\left[{\rm i}t'\mu-\frac{1}{2}t'\Sigma t \right],\\ X\sim N_p(\mu,\Sigma). \]

3. 分量的线性组合（可理解为一种性质，用于证否）：若\(p\)维随机向量\(X\)的任意线性组合均服从正态分布，则\(X\)服从多元正态分布。

4. 联合密度法：对于\(\Sigma>0\)和常数向量\(\mu\)，多元正态分布随机向量\(X\)具有如下的联合密度：

   \[f(x)=\frac{1}{(2\pi)^{p/2}|\Sigma|^{1/2}}\exp\left[-\frac{1}{2}(x-\mu)'\Sigma^{-1}(x-\mu) \right],\\ X\sim N_p(\mu,\Sigma). \]

   这给出非退化正态随机向量的定义方式。

多元正态分布的常用性质：

1. 如果\(X\sim N_p(\mu,\Sigma)\)，常数阵\(C_{q\times p}\)，则

   \[CX\sim N_q(C\mu, C\Sigma C'). \]

2. 如果\(X\sim N_p(\mu,\Sigma\)，将\(X\)分成\(r\)维的\(X^{(1)}\)与\(p-r\)维的\(X^{(2)}\)，则

   \[X^{(1)}\sim N_r(\mu^{(1)}, \Sigma_{11}), \\ X^{(2)}\sim N_{p-r}(\mu^{(2)}, \Sigma_{22}), \\ \mathbb{COV}(X^{(1)}, X^{(2)})=\Sigma_{12}=\Sigma_{21}'. \]

3. 多元正态分布分量间不相关性与独立性等价。

4. 独立性投影：\(X^{(1)}\)在\(X^{(2)}\)方向上的投影是

   \[\Sigma_{12}\Sigma_{22}^{-1}X^{(2)}, \]

   这相当于\(X^{(1)}-\Sigma_{12}\Sigma_{22}^{-1}X^{(2)}\)与\(X^{(2)}\)不相关（即独立）。由此，可以作可逆线性变换如下：

   \[B=\begin{bmatrix} I_r & -\Sigma_{12}\Sigma_{22}^{-1} \\ O & I_{p-r} \end{bmatrix},\\ Z=BX=\begin{bmatrix} X^{(1)}-\Sigma_{12}\Sigma_{22}^{-1}X^{(2)} \\ X^{(2)} \end{bmatrix}\xlongequal{def}\begin{bmatrix} Z^{(1)} \\ Z^{(2)} \end{bmatrix}. \]

   由此，再结合随机向量函数密度公式，可以计算得到

   \[(X^{(1)}|X^{(2)})\sim N_r(\mu_{1\cdot 2}, \Sigma_{11\cdot2}),\\ \mu=\mu^{(1)}+\Sigma_{12}\Sigma_{22}^{-1}(X^{(2)}-\mu^{(2)}),\\ \Sigma_{11\cdot 2}=\Sigma_{11}-\Sigma_{12}\Sigma_{22}^{-1}\Sigma_{21}. \]

1-2 参数估计与常用统计量

以下假设\(X\)是\(n\times p\)数据矩阵，第\(i\)行代表第\(i\)个样品的观测值\(X_{(i)}\)；\(\boldsymbol{1}_n\)代表\(n\)维纯\(1\)列向量。常用统计量：

1. 样本均值向量\(\bar X\)：常用作样本均值的估计。

   \[\bar X=\frac{1}{n}\sum_{i=1}^nX_{(i)}=\frac{1}{n}X'\boldsymbol{1}_n. \]

2. 样本离差阵\(A\)：类比一元情况下的离差平方和，是一个\(p\times p\)矩阵。

   \[A=\sum_{i=1}^n(X_{(i)}-\bar X)(X_{(i)}-\bar X)'=X'X-\frac{1}{n}\bar X\bar X'=X'(I_n-\frac{1}{n}\boldsymbol{1}_n\boldsymbol{1}'_n)X. \]

3. 样本协方差阵\(S\)：类比一元情况下的样本方差，常用作自协方差矩阵的估计。

   \[S=\frac{1}{n-1}A\xlongequal{def}(s_{ij})_{p\times p}. \]

4. 样本相关阵\(R\)：常用作相关阵的估计。

   \[R\xlongequal{def}(r_{ij})_{p\times p},\quad r_{ij}=\frac{s_{ij}}{\sqrt{s_{ii}}\sqrt{s_{jj}}}. \]

   在实际应用中，要求相关阵，一般先对数据矩阵进行标准化，再求协方差阵即可。标准化指的是扣掉均值、除以标准差。

正态总体\(N_p(\mu,\Sigma)\)中，常用统计量的性质：

1. \(\bar X\sim N_p(\mu,\Sigma/n)\)。

2. 设\(Z_1,\cdots,Z_{n-1}\)独立同分布于\(N_p(0,\Sigma)\)，则

   \[A\xlongequal{d}\sum_{i=1}^{n-1}Z_iZ_i'. \]

3. \(\bar{X}\)和\(A\)相互独立。

4. \(\mathbb{P}(A>0)=1\Leftrightarrow n>p\)。

5. 前三点通过构造正交矩阵\(\Gamma\)（最后一行为\(\frac{1}{\sqrt{n}}\)）来证明，第四点令\(B=(Z_1,\cdots,Z_{n-1})\)，有\(A=BB'\)，由\(r(A)=r(B)\)只需证明\(B\)的前\(p\)列线性相关概率为0即可。

多元正态分布参数\((\mu,\Sigma)\)的极大似然估计：

• 如果\(\mu,\Sigma\)均未知，则

  \[\hat{\mu}=\bar{X}, \hat{\Sigma}=\frac{1}{n}A. \]

  注意极大似然估计不是样本协方差阵，经过无偏调整后才是样本协方差阵。

• 如果\(\mu\)已知，则

  \[\hat{\Sigma}=\frac{1}{n}\sum_{i=1}^n(X_{(i)}-\bar{X})(X_{(i)}-\bar{X})'. \]

• 用到的矩阵求导公式：\(A\)为实对称矩阵时，有

  \[\frac{{\rm d}\ln|A|}{{\rm d}A}=A^{-1},\\ \frac{{\rm d}x'Ax}{{\rm d}A}=xx'. \]

  似然函数为

  \[L(\mu,\Sigma)=\frac{1}{(2\pi)^{np/2}|\Sigma|^{n/2}}\exp\left[-\frac{1}{2}\sum_{i=1}^n(x_{(i)}-\mu)'\Sigma^{-1}(x_{(i)}-\mu) \right], \\ l(\mu,\Sigma)=C+\frac{n}{2}\ln|\Sigma^{-1}|-\frac{1}{2}\sum_{i=1}^n(x_{(i)}-\mu)\Sigma^{-1}(x_{(i)}-\mu). \]

• \(\hat{\Sigma}\)不是无偏的而\(S\)是无偏的；\(\bar{X}, S\)是\(\mu, \Sigma\)的UMVUE，而且是相合的。

• MLE具有可映射性，对于Borel可测映射\(g(\theta)\)，如果\(\hat\theta\)是\(\theta\)的MLE，则\(g(\hat{\theta})\)是\(g(\theta)\)的MLE。

1-3 常用分布

如同一元数理统计一般，为了给出区间估计，需要构造一些常用分布。本节过于理论，建议结合1-4理解。

设\(X=(X_1,\cdots,X_p)\)且各分量相互独立，卡方分布与二次型切实相关，有以下相关结论，但卡方分布相关的结论在多元统计分析中并不常用：

1. 如果对所有分量都有\(\mu_i=0\)，\(\sigma^2_i=\sigma^2\)，则

   \[\frac{1}{\sigma^2}X'X\sim \chi^2(n). \]

2. 如果\(\mu\ne\boldsymbol{0}\)，令\(\delta\)为非中心参数，则

   \[\frac{1}{\sigma^2}X'X\sim\chi^2_n(\delta),\quad \delta\xlongequal{def}\frac{\mu'\mu}{\sigma^2}. \]

3. 若\(X\sim N_n(\boldsymbol{0}_n,\sigma^2I_n)\)，\(A\)为对称矩阵且\(r(A)=r\)，则

   \[\frac{X'AX}{\sigma^2}\sim\chi^2(r)\Leftrightarrow A^2=A. \]

4. 若\(X\sim N_n(\boldsymbol{0}_n,I_n)\)，\(A\)为\(n\)阶对称阵，\(B\)为\(m\times n\)矩阵，则\(BA=O\Leftrightarrow X'AX\)与\(BX\)相互独立。

5. 若\(X\sim N_n(\boldsymbol{0}_n,I_n)\)，\(A\)和\(B\)都是\(n\)阶对称阵，则\(AB=O\Leftrightarrow X'AX\)和\(X'BX\)相互独立。

6. 均值检验 若\(X\sim N_p(\boldsymbol{0}_p,\Sigma)\)，且\(\Sigma>0\)，则

   \[X'\Sigma^{-1}X\sim \chi^2(p). \]

   如果\(X\sim N_p(\mu,\Sigma)\)，则

   \[X'\Sigma^{-1}X\sim \chi^2_p(\delta),\quad \delta=\mu'\Sigma^{-1}\mu. \]

7. 若\(X\sim N_p(\mu,\Sigma)\)，\(A\)为对称矩阵且\(r(A)=r\)，则

   \[(X-\mu)'\Sigma^{-1}(X-\mu)\sim \chi^2(r)\Leftrightarrow \Sigma A\Sigma=\Sigma A\Sigma A\Sigma. \]

8. 若\(X\sim N_p(\mu, \Sigma)\)，\(\Sigma>0\)且\(A\)和\(B\)都是\(p\)阶对称阵，则\((X-\mu)'A(X-\mu)\)与\((X-\mu)'B(X-\mu)\)独立等价于

   \[\Sigma A\Sigma B\Sigma = O_{p\times p}. \]

Wishart分布常用于表述样本离差阵、样本协方差阵的分布：设\(X_{(\alpha)}\sim N_p(0,\Sigma)\)相互独立，\(X=(X_{(1)}, \cdots, X_{(n)})'\)，则

\[W_{}=\sum_{\alpha=1}^n X_{(\alpha)}X_{(\alpha)}'=X'X\sim W_p(n,\Sigma). \]

Wishart分布的性质（个人认为这里不是很重要）：

1. Wishart分布可以看成卡方分布的推广，其中\(p\)代表矩阵阶数，除此之外的两个参数是\(n\)（代表独立同分布多元正态分布的个数）和\(\Sigma\)（零均值正态分布的协方差阵）。

2. 设\(X_{(\alpha)}\sim N_p(\mu,\Sigma)\)相互独立，则样本离差阵\(A\)服从Wishart分布，即

   \[A=\sum_{\alpha=1}^{n}(X_{(\alpha)}-\bar X)(X_{(\alpha)}-\bar X)'\sim W_p(n-1,\Sigma). \]

   注意：自由度为\(n-1\)，与正态分布样本方差类似。

3. Wishart分布关于自由度具有可加性，即\(W_i\sim W_p(n_i,\Sigma)\)相互独立，则

   \[W\xlongequal{def}\sum_{i=1}^n W_i\sim W_p\left(\sum_{i=1}^nn_i,\Sigma\right). \]

4. 对称线性变换：若\(W\sim W_p(n,\Sigma)\)，则对于\(m\times p\)常数阵\(C\)，有

   \[CWC'\sim W_m(n, C\Sigma C'). \]

5. 分块Wishart分布：独立同分布的\(X_{(\alpha)}\sim N_p(0, \Sigma)\)，如果可以对\(\Sigma\)分块成\(\Sigma_{11},\Sigma_{12},\Sigma_{21},\Sigma_{22}\)，则\(W=\sum_{\alpha=1}^n X_{(\alpha)}X_{(\alpha)}'\)也可以分块成\(W_{11},W_{12},W_{21},W_{22}\)，并且

   \[W_{11}\sim N_r(n, \Sigma_{11}),W_{22}\sim N_{p-r}(n, \Sigma_{22}),\\ \]

   \(\Sigma_{12}=O\)时，\(W_{11}\)和\(W_{22}\)相互独立。

6. Wishart分布的期望：\(W\sim W_p(n,\Sigma)\)，则

   \[\mathbb{E}(W)=n\Sigma. \]

霍特林\(T^2\)分布：设\(X\sim N_p(0,\Sigma)\)与\(W_p(n,\Sigma)\)相互独立，则

\[T^2=X'(\frac{W}{n})^{-1}X\sim T^2(p,n). \]

霍特林\(T^2\)分布的性质（比较重要）：

1. 均值检验 设\(X_{(\alpha)}\)独立同分布于\(N_p(\mu,\Sigma)\)，则由于

   \[X_{(\alpha)}-\mu \sim N_p(0,\Sigma),\\ \sqrt{n}(\bar X-\mu) \sim N_p(0,\Sigma) \\ A\sim W_p(n-1,\Sigma), \]

   所以

   \[n(n-1)(\bar X-\mu)'A^{-1}(\bar X-\mu)=n(\bar X-\mu)'S^{-1}(\bar X-\mu)\sim T^2(p,n-1). \]

   这一性质常常用于协方差阵未知时的均值检验。

2. \(T^2\)分布与\(F\)分布的关系：如果\(T^2\sim T^2(p, n)\)，则

   \[\frac{n-p+1}{np}T^2\sim F(p,n-p+1). \]

   由此，可以将\(T^2\)分布统计量转化为\(F\)统计量进行假设检验。

3. \(T^2\)分布与\(\Sigma\)无关，因为

   \[nU'W_I^{-1}U\sim T^2(p, n), \]

   这里\(U\sim N_p(0,I_p)\)，\(W_I\sim W_p(n, I_p)\)，而

   \[X\xlongequal{d}\Sigma^{1/2}U, \\ W\xlongequal{d}\Sigma^{1/2}W_I\Sigma^{1/2}, \\ nU'W_I^{-1}U\xlongequal{d}nX'\Sigma^{-1/2}\Sigma^{-1/2}W_I^{-1}\Sigma^{1/2}\Sigma^{-1/2}X =nX'W^{-1}X. \]

4. 在假设检验中，\(T^2\)统计量对非退化变换保持不变，即对\(X_{(\alpha)}\)对\(\mu\)作假设检验的效果，与对\(CX_{(\alpha)}+d\)对\(C\mu+d\)作假设检验的效果一致。

威尔克斯\(\Lambda\)分布：设\(A_1\sim W_p(n_1,\Sigma)\)，\(A_2\sim W_p(n_2,\Sigma)\)，这里\(\Sigma>0\)，\(n_1\ge p\)，则

\[\Lambda=\frac{|A_1|}{|A_1+A_2|}\sim \Lambda(p, n_1, n_2). \]

威尔克斯\(\Lambda\)分布的性质：

1. 若\(\Lambda\sim \Lambda(p, n_1, n_2)\)，则

   \[\Lambda\xlongequal{d}B_1B_2\cdots B_p, \]

   这里\(B_k\)相互独立，服从不同的分布为

   \[B_k\sim \beta\left(\frac{n_1-p+k}{2}, \frac{n_2}{2} \right). \]

2. 若\(n_2<p\)，则

   \[\Lambda(p, n_1, n_2)\xlongequal{d}\Lambda(n_2, p, n_1+n_2-p). \]

1-4 均值假设检验

方差已知时，均值向量的假设检验：

\[H_0:\mu=\mu_0\Leftarrow H_1:\mu\ne\mu_0. \]

利用

\[\bar X\sim N_p(\mu,\frac{\Sigma}n), \]

可以得到

\[\sqrt{n}(\bar X-\mu)\sim N_p(\mu,\Sigma). \]

利用正态分布的二次型，可以构造

\[n(\bar X-\mu)'\Sigma^{-1}(\bar X-\mu)\stackrel{H_0}\sim \chi^2(p). \]

方差未知时，均值向量的假设检验，利用

\[\sqrt{n}(\bar X-\mu)\sim N_p(0,\Sigma),\\ A=\sum_{\alpha=1}^{n}(X_{(\alpha)}-\bar X)(X_{(\alpha)}-\bar X)'\sim W_p(p, n-1), \]

可以得到

\[n(\bar X-\mu)'\left(\frac{A}{n-1} \right)^{-1}(\bar X-\mu)\stackrel{H_0}{\sim }T^2(p, n-1). \]

这里\(A/(n-1)\)就是样本自协方差矩阵\(S\)。

联合协方差阵：设\(X_{(\alpha)}\sim N_p(\mu_1,\Sigma)\)，\(Y_{(\alpha)}\sim N_p(\mu_2, \Sigma)\)，即两总体有相同的协方差阵，则样本联合协方差阵为

\[S_\text{pooled}=\frac{A_1+A_2}{m+n-2}. \]

这里

\[A_1=\sum_{\alpha=1}^{m}(X_{(\alpha)}-\bar X)(X_{(\alpha)}-\bar X)',\\ A_2=\sum_{\alpha=2}^{n}(Y_{(\alpha)}-\bar Y)(Y_{(\alpha)}-\bar Y)',\\ A_1+A_2\sim W_p(m+n-2, \Sigma). \]

方差未知但相等时，两样本均值检验：

\[H_0:\mu_1=\mu_2\Leftrightarrow H_1:\mu_1\ne \mu_2. \]

利用

\[\bar X\sim N_p\left(\mu_1,\frac{\Sigma}{m}\right),\bar Y\sim N_p\left(\mu_2,\frac{\Sigma}{n}\right),\\ \bar X-\bar Y\sim N_p\left(\mu_1-\mu_2,\left(\frac{1}{m}+\frac{1}{n} \right)\Sigma\right), \]

得到

\[\sqrt{\frac{mn}{m+n}}[\bar X-\bar Y-(\mu_1-\mu_2)]\sim N_p(0,\Sigma),\\ A_1+A_2\sim W_p(m+n-2,\Sigma), \]

所以

\[\frac{mn}{m+n}(\bar X-\bar Y)'\left(\frac{A_1+A_2}{m+n-2} \right)^{-1}(\bar X-\bar Y)\stackrel{H_0}{\sim }T^2(p,m+n-2),\\ \frac{mn}{m+n}(\bar X-\bar Y)S_{\text{pooled}}^{-1}(\bar X-\bar Y)\stackrel{H_0}{\sim }T^2(p,m+n-2). \]

方差不等但样本数量相等时，构造成对数据即

\[Z_{(\alpha)}=X_{(\alpha)}-Y_{(\alpha)}\sim N_p(\mu_1-\mu_2,\Sigma_1+\Sigma_2)\xlongequal{def}N_p(\mu_Z,\Sigma_Z). \]

将其转化为一元问题。

多元均值检验：对于\(k\)个总体\(N_p(\mu_i,\Sigma)\)，\(i=1,2,\cdots,k\)，它们的协方差阵相等，检验问题是

\[H_0:\mu_1=\cdots=\mu_k\Leftrightarrow H_1:\text{otherwise}. \]

离差阵分解：定义\(\bar X\)为所有样本的均值，\(\bar X^{(i)}\)为第\(i\)组的样本均值，则

\[T=A+B,\\ n=\sum_{i=1}^k n_i, \\ \bar X=\frac{1}{n}\sum_{i=1}^{k}\sum_{\alpha=1}^{n_j}X_{(\alpha)}^{(i)}, \\ \bar X^{(i)}=\frac{1}{n_i}\sum_{\alpha=1}^{n_i}X_{(\alpha)}^{(i)}. \\ T=\sum_{i=1}^{k}\sum_{\alpha=1}^{n_i}(X_{(\alpha)}^{(i)}-\bar X)(X_{(\alpha)}^{(i)}-\bar X)',\\ A = \sum_{i=1}^{k}\sum_{\alpha=1}^{n_i}(X_{(\alpha)}^{(i)}-\bar X^{(i)})(X_{(\alpha)}^{(i)}-\bar X^{(i)})',\\ B=\sum_{i=1}^k n_i(\bar X^{(i)}-\bar X)(\bar X^{(i)}-\bar X)'. \]

有相关结论：

1. \(A\)为组内离差阵，如果\(A_i\)为第\(i\)组的离差阵，则

   \[A_i\sim W_p(n_i-1, \Sigma),\\ A=\sum_{i=1}^{k}A_i\sim N(n-k,\Sigma). \]

2. \(B\)为组间离差阵，如果\(H_0\)成立，有\(B\sim W_p(k-1,\Sigma)\)，且 \(A\)与\(B\)相互独立。

3. \(T\)为总离差阵，如果\(H_0\)成立，有\(T=A+B\sim W_p(n-1,\Sigma)\)。

因此，构造威尔克斯\(\Lambda\)统计量，有

\[\Lambda=\frac{|A|}{|A+B|}=\frac{|A|}{|T|}\stackrel{H_0}{\sim}\Lambda(p, n-k, k-1). \]

1-5 其他假设检验

似然比检验法：要检验的假设是\(H_0\)，似然函数为\(L(\Theta)\)，则似然比统计量为

\[\lambda=\frac{\max_{H_0}L(\Theta_0)}{\max L(\Theta)}\in[0, 1]. \]

这里\(\Theta\)为参数空间，\(\Theta_0\)为\(H_0\)下参数空间，则

\[-2\ln\lambda \stackrel{\text{approx}}{\sim }\chi^2(f),\\ f=\dim(\Theta) - \dim(\Theta_0). \]

协方差阵为单位阵的检验：对于单总体\(N_p(\mu,\Sigma)\)，要检验的假设是

\[H_0:\Sigma = I_p\Leftrightarrow H_1:\Sigma \ne I_p. \]

构造似然比检验统计量为

\[\lambda = \frac{\max_{\mu,\Sigma_0}L(\mu,\Sigma_0)}{\max_{\mu,\Sigma}L(\mu,\Sigma)}=\frac{L(\bar X,I_p)}{L(\bar X, \frac{A}{n})}\\ -2\ln \lambda \stackrel{\text{approx}}{\sim} \chi^2\left(\frac{p(p+1)}{2} \right). \]

一般协方差阵的检验：要检验的假设是

\[H_0:\Sigma=\Sigma_0\Leftrightarrow H_1:\Sigma\ne \Sigma_0, \]

可以使用某个线性变换：\(\Sigma_0^{-1/2}X\sim N(\Sigma_0^{-1/2}\mu,I_0)\)，转变成单位阵的检验。

检验如下的假设：

\[H_0:\Sigma=\sigma^2\Sigma_0\Leftrightarrow H_0:\Sigma\ne\sigma^2\Sigma_0. \]

这里\(\sigma^2\)未知但\(\Sigma_0\)已知，似然比检验统计量为

\[\lambda=\frac{\max_{\mu,\sigma}L(\mu,\sigma^2\Sigma_0)}{\max L(\mu,\Sigma)}=\frac{L(\bar X, \frac{1}{np}{\rm tr}(\Sigma_0^{-1}A)\cdot \Sigma_0)}{L(\bar X,\frac{A}{n})},\\ -2\ln\lambda \stackrel{\text{approx}}{\sim} \chi^2\left(\frac{p(p+1)}{2}-1 \right). \]

多总体协方差阵的检验：有\(k\)个总体\(X_k\sim N(\mu_i,\Sigma_i)\)，每一个总体中抽取样本容量为\(n_i\)，一共\(n\)个样本，检验假设为

\[H_0:\Sigma_1=\Sigma_2=\cdots=\Sigma_k\Leftrightarrow H_1:\text{otherwise}. \]

联合似然函数为

\[L(\mu_1,\Sigma_1,\cdots,\mu_k,\Sigma_k)=\prod_{i=1}^kL_i(\mu_i,\Sigma_i). \]

似然比检验统计量为（注意这里的\(A\)是组内离差阵而非总离差阵\(T\)）

\[\lambda=\frac{\prod_{i=1}^kL_i(\bar X^{(i)},\frac{A}{n})}{\prod_{i=1}^k L_i(\bar X_i^{(i)},\frac{A_i}{n_i})}=\frac{|\frac{A}{n}|^{-n/2}}{\prod_{i=1}^k|\frac{A_i}{n_i}|^{-n_i/2}},\\ -2\ln\lambda \stackrel{\text{approx}}{\sim }\chi^2\left(\frac{p(p+1)(k-1)}{2} \right). \]

一般说来，\(\mu\)的极大似然估计为\(\bar X\)，每一组组内协方差阵\(A_i\)的极大似然估计为\(\frac{A_i}{n_i}\)，联合协方差阵的极大似然估计为\(\frac{A}{n}\)。也可以进行无偏修正，令

\[\hat{\Sigma}_i=\frac{A_i}{n_i-1},\quad \hat{\Sigma}_\text{pooled}=\frac{A}{n-k}. \]

独立性检验：检验两个总体\(X^{(1)}\)、\(X^{(2)}\)是否相互独立，即检验\(\Sigma_{12}=O\)是否成立。用似然比检验法，

\[\lambda=\frac{L(\bar X^{(1)},\frac{A_{11}}{n},\bar X^{(2)},\frac{A_{22}}{n})}{L(\bar X^{(1)},\bar X^{(2)},\frac{A}{n})}=\frac{|\frac{A_{11}}{n}|^{-n/2}|\frac{A_{22}}{n}|^{-n/2}}{|\frac{A}{n}|^{-n/2}}=\left(\frac{|A|}{|A_{11}A_{22}|} \right)^{n/2},\\ -2\ln \lambda \stackrel{\text{approx}}{\sim}\chi^2(r(p-r)). \]

如果要检验\(k\)个总体的是否相互独立，则

\[\lambda=\frac{L(\bar X^{(1)},\frac{A_{11}}{n},\cdots,\bar X^{(k)},\frac{A_{kk}}{n})}{L(\bar X^{(1)},\cdots,\bar X^{(k)},\frac{A}{n})}=\frac{\prod_{i=1}^n|\frac{A_{ii}}{n}|^{-n/2}}{|\frac{A}{n}|^{-n/2}}=\left(\frac{|A|}{\prod_{i=1}^n |A_{ii}|} \right)^{n/2} \\ -2\ln\lambda\stackrel{\text{approx}}{\sim}\chi^2\left(\frac{p(p+1)}{2}-\sum_{i=1}^k\frac{p_i(p_i+1)}{2} \right) \]

要注意的是，此时各个分块的极大似然估计为

\[\hat\Sigma_{ii}=\frac{A_{ii}}{n}, \]

与之前在多总体检验中的不同（分母、阶数都不同）。

Part 2：归类

2-1 距离判别

基本思想：样本距离哪个总体最近，就判断它属于哪一个总体。

度量方式：马氏距离，设\(G\)为\(p\)元总体，均值为\(\mu\)，自协方差矩阵\(\Sigma\)，则样品\(X\)与\(G\)的马氏距离为

\[d^2(X,G)=(X-\mu)'\Sigma^{-1}(X-\mu) \]

实际应用中如果\(\mu\)和\(\Sigma\)是未知量，则从样本均值\(\bar X\)和样本自协方差阵\(S\)出发计算马氏距离。

联合协方差阵：在判别分析中，我们常常假定总体有相同的自协方差矩阵\(\Sigma\)，从而判别函数是线性的。\(k\)个总体的联合协方差阵是

\[S_\text{pooled}=\frac{A_1+A_2+\cdots+A_k}{n-k}. \]

此时的线性判别函数为

\[W(X)=\left(X-\frac{1}{2}(\bar X^{(1)}+\bar X^{(2)})\right)'S^{-1}(\bar X^{(1)}-\bar X^{(2)}). \]

其最后一部分是\(\bar X^{(1)}-\bar X^{(2)}\)，这意味着，如果\(W(X)>0\)，则判\(X\in G_1\)；如果\(W(X)<0\)，则判\(X\in G_2\)。

错判概率：假设\(X^{(i)}\sim N_p(\mu_i,\Sigma)\)并假定其已知，则

\[D\xlongequal{def}(\mu_1-\mu_2)'\Sigma^{-1}(\mu_1-\mu_2),\\ \begin{aligned} \mathbb{P}(2|1)=&\mathbb{P}\left((\mu_1-\mu_2)'\Sigma^{-1}X<\frac{1}{2}(\mu_1-\mu_2)'\Sigma^{-1}(\mu_1+\mu_2)\bigg|X\sim N_p(\mu_1,\Sigma) \right)\\ =&\mathbb{P}\left((\mu_1-\mu_2)'\Sigma^{-1}(X-\mu_1)<-\frac{1}{2}D\bigg|X\sim N_p(\mu_1,\Sigma) \right)\\ =&\mathbb{P}\left(N(0,1)<-\frac{1}{2}\sqrt{D}\bigg|X\sim N_p(\mu_1,\Sigma) \right)\\ =&\Phi(-\frac{1}{2}\sqrt{D})<\frac{1}{2}. \end{aligned} \]

这表明线性判别的错判概率不会高于50%.

2-2 贝叶斯判别

贝叶斯判别的目标是使得错判损失的期望值最小。

• 先验概率：\(q_i\)，表明样本出自某个类\(i\)的概率。
• 错判损失：\(L(j|i)\)，表明样本出自\(i\)类却被错判成\(j\)类的损失大小，常常\(L(j|i)=1-\delta_{j-i}\)。
• 判别法：一种划分，\(\mathbb{R}^n=\{D_1,D_2,\cdots,D_k\}\)，将样本空间划分成\(k\)块。

贝叶斯判别：设有\(k\)个总体\(G_1,\cdots,G_k\)，每一类\(G_i\)的联合密度为\(f_i(X)\)，先验概率为\(q_i\)，错判损失为\(L(j|i)\)，贝叶斯判别要使得

\[g(D)=\sum_{i=1}^kq_i\sum_{j=1}^k\mathbb{P}(j|i)L(j|i)=\sum_{i=1}^k q_i\sum_{j=1}^kL(j|i)\int_{D_j^n}f_i(X){\rm d}X. \]

将样品判别到第\(i\)类的平均损失：

\[h_i(X)=\sum_{j=1}^k q_jf_j(X)L(i|j) \]

贝叶斯判别的依据是，在\(h_1(X),\cdots,h_k(X)\)中，选择判别平均损失最小的一个。如果\(L(j|i)=1-\delta_{j-i}\)，则

\[h_i(X)<h_j(X)\Leftrightarrow q_if_i(X)>q_jf_j(X), \]

所以等价于选择\(q_if_i(X)\)最大的那个\(i\)作为贝叶斯判别的归类，看起来就像是来自哪个总体的概率最大，就归为哪一类。

同损失正态总体：

\[q_if_i(X)=\frac{q_i}{(2\pi)^{p/2}|\Sigma_i|^{1/2}}\exp\left(-\frac{1}{2}(X-\mu_i)'\Sigma^{-1}(X-\mu_i) \right),\\ \ln(q_if_i(X))=\ln|q_i|-\frac{1}{2}\ln|\Sigma_i|-d^2(X,G_i). \]

所以基于贝叶斯判别的广义平方距离为

\[\bar d^2(X,G_i)=d^2(X,G_i)+\ln|\Sigma_i|-2\ln q_i . \]

它没有考虑损失函数的设置。

2-3 费希尔判别

费希尔判别的方法是对样本做某个线性变换，使不同类之间的差距尽可能大，表现为组间离差阵（组间平方和）占总离差阵（总离差平方和）的比例大。

考虑双总体\(X^{(1)}\sim N_p(\mu_1,\Sigma_1)\)和\(X^{(2)}\sim N_p(\mu_2,\Sigma_2)\)，在没有给定样本的情况下，视为各组使用相同的样本容量，则组内离差阵和组间离差阵分别为

\[A=\Sigma_1+\Sigma_2,\\ B=(\mu_1-\bar \mu)'(\mu_1-\bar \mu)+(\mu_2-\bar\mu)'(\mu_2-\bar\mu),\\ \bar\mu=\frac{1}{2}(\mu_1+\mu_2). \]

如果给定了各组样本，则可以计算组内离差阵\(A\)和组间离差阵\(B\)。目标是使用线性变换\(a\)，让\(a'Ba/a'Aa\)最大，并规定\(a'Aa=1\)，变为如下的规划问题：

\[\max {a'Ba}=\Delta(a),\\ \text{s.t. } a'Aa=1. \]

解得\(a'Ba\)为\(A^{-1}B\)的特征值，\(a\)为相应的特征向量（使得\(a'Aa=1\)，但实际应用中是否归一化对问题没有影响）。

结论：设\(A^{-1}B\)的非零特征值为\(\lambda_1\ge\lambda_2\ge \cdots\lambda_r>0\)，相应的满足约束的特征向量为\(l_1,l_2,\cdots,l_r\)，则前\(l\)个线性判别函数为

\[l_1'X,\cdots,l_l'X, \]

其判别能力为

\[P_l=\frac{\lambda_1+\cdots+\lambda_l}{\lambda_1+\cdots+\lambda_r}. \]

实际应用时，常取使得\(P_l\ge 0.7\)的\(l\)值。

2-4 系统聚类

点间度量：闵可夫斯基距离（绝对值距离、欧氏距离、切比雪夫距离），兰氏距离，马氏距离，斜交空间距离等。

系统聚类步骤：

1. 初始每个样品视为一类，计算样品间距离
2. 每一次选择两个距离最小的类，合并，总类数减少1
3. 重新计算新类与其他类的距离
4. 直到总类数为1时，停止合并过程
5. 绘制谱系聚类图，决定分类个数

聚类方法：以下定义\(D_{pq}\)为类\(G_p,G_q\)间的距离，\(\bar X_i\)为\(G_i\)的均值，\(n_i\)为\(G_i\)内的样本数量，\(D_{rk}\)为由\(G_p,G_q\)合并后的新类\(G_r\)与其他类\(G_k\)间的距离。

1. 最短距离法(single)：定义类间距离为，两类中相距最近的样品之间的距离。

   \[D_{pq}=\min_{i\in G_p,j\in G_q} d^2_{ij},\\ D_{rk}=\min\{D_{pk},D_{qk} \}. \]

2. 最长距离法(complete)：定义类间距离为，两类中相距最远的样品之间的距离。

   \[D_{pq}=\max_{i\in G_p,j\in G_q} d^2_{ij},\\ D_{rk}=\max\{D_{pk},D_{qk} \}. \]

3. 中间距离法(median)：基于递推的距离定义法，根据参数的不同而不同。

   \[D^2_{rk}=\frac{1}{2}(D_{pk}^2+D_{qk}^2)+\beta D_{pq}^2,\quad \beta\in[-\frac14, 0]. \]

   常取\(\beta=-1/4\)，此时\(D_{rk}\)就是以\(D_{pk},D_{qk},D_{pq}\)为边的三角形中，\(D_{pq}\)边上的中线长度。

4. 重心法(centroid)：每一类的重心即均值，定义类间距离类重心间的距离。

   \[D_{pq}=d^2(\bar X_p, \bar X_q),\\ D_{rk}^2=\frac{n_p}{n_r}D_{pk}^2+\frac{n_q}{n_r}D_{qk}^2-\frac{n_pn_q}{n_r^2}D_{pq}^2. \]

5. 类平均法(average)：用两类中所有样品之间距离的平均作为类间距离。

   \[D_{pq}=\frac{1}{n_pn_q}\sum_{i\in G_p,j\in G_q} d^2_{ij},\\ D_{rk}^2=\frac{n_p}{n_r}D^2_{pk} + \frac{n_q}{n_r}D^2_{qk}. \]

6. MCQ相似分析法(mcquitty)：使用如下的递推公式。

   \[D^2_{rk}=\frac{D_{pk}^2+D_{qk}^2}{2}. \]

7. 离差平方和法(ward)：类间平方距离视为合并类后会增加的离差平方和，它是实际应用中效果最好的方法。

   \[D_{pq}^2=W_r-(W_p+W_q), \\ D_{rk}^2=\frac{n_p+n_k}{n_r+n_k}D^2_{pk}+\frac{n_q+n_k}{n_r+n_k}D^2_{qk}-\frac{n_k}{n_r+n_k}D^2_{pq}. \]

统一表达式：

\[D_{rk}^2=\alpha_pD_{pk}^2+\alpha_qD_{qk}^2+\beta D_{pq}^2+\gamma|D_{pk}^2-D_{qk}^2|. \]

应用此公式时，重心法和Ward法要求使用欧氏距离。

Part 3：降维

3-1 主成分分析

基本思想：用总体的线性组合形成新的变量来代替原有变量，要求尽可能保留多的信息，数量尽可能少。

主成分：设\(X=(X_1,\cdots,X_p)'\)为\(p\)维随机向量，\(Z_i=a_i'X\)为\(X\)的第\(i\)主成分，如果

\[a_i'a_i=1,\\ \forall j<i,\quad {\rm Cov}(Z_i, Z_j) = a_i'\Sigma a_j=0,\\ \mathbb{D}(Z_i)=\max_{a'a=1, a'\Sigma a=0}\mathbb{D}(a'X). \]

求解第一主成分：即求解以下的线性规划问题：

\[\max a_1'\Sigma a_1',\\ \text{s.t. } a_1'a_1 = 1. \]

其拉格朗日函数是

\[L(a_1,\lambda)=a_1'\Sigma a_1-\lambda(a_1'a_1-1), \]

解得\(\lambda_1=a_1'\Sigma a_1'\)，且\(\lambda_1\)是\(\Sigma\)的最大特征值，\(a_1\)是对应的单位正交向量。

主成分的求法：设\(X=(X_1,\cdots,X_p)'\)且\(\mathbb{D}(X)=\Sigma\)，\(\Sigma\)的特征值为\(\lambda_1\ge \lambda_2 \ge \cdots \ge\lambda_p\ge 0\)，相应的单位正交向量为\(a_1,\cdots,a_p\)，则\(X\)的第\(i\)主成分为\(Z_i=a_i'X\)。

• 令\(Z=(Z_1,\cdots,Z_p)'\)，则\(Z=A'X\)，这里\(A\)是正交阵，即

  \[A = \begin{bmatrix} a_1 & a_2 & \cdots & a_p \end{bmatrix}=\begin{bmatrix} a_{11} & a_{12} & \cdots & a_{1p} \\ a_{21} & a_{22} & \cdots & a_{2p} \\ \vdots & \vdots & & \vdots \\ a_{p1} & a_{p2} & \cdots & a_{pp}. \end{bmatrix} \]

  矩阵\(A\)中，每一列代表一个主成分系数，即一个特征向量。

• \(\mathbb{D}(Z)={\rm diag}(\lambda_1,\cdots,\lambda_p)\)，且\(\lambda_1\ge \lambda_2\ge \cdots \lambda_p\)。这也是\(Z\)的各个分量依次是\(X\)的第\(i\)主成分的充要条件。

总方差（总惯量）：令\(\mathbb{D}(X)=\Sigma=(\sigma_{ij})_{p\times p}\)，则\(\sum_{i=1}^p \sigma_{ii}\)称为原总体的总方差（总惯量），即各分量方差之和。有

\[\sum_{i=1}^p \sigma_{ii} = \sum_{i=1}^p \lambda_i. \]

因子载荷：把主成分\(Z_k\)与原始变量\(X_i\)之间的相关系数，称为因子载荷量，即

\[\rho(Z_k, X_i)=\frac{{\rm Cov}(a_k'X, e_i'X)}{\sqrt{\mathbb{D}(Z_k)\mathbb{D}(X_i)}}=\frac{\lambda_ke_i'a_k}{\sqrt{\lambda_k\sigma_{ii}}}=\frac{\sqrt{\lambda_k}a_{ik}}{\sqrt{\sigma_{ii}}}.\\ A'\Sigma A=\Lambda\Rightarrow \Sigma=A\Lambda A',\sigma_{ii}=\sum_{k=1}^p\lambda_k a_{ik}^2. \\ \sum_{k=1}^p\rho^2(Z_k, X_i)=\sum_{k=1}^p\frac{\lambda_k a_{ik}^2}{\sigma_{ii}}=1,\\ \sum_{i=1}^p \sigma_{ii}\rho^2(Z_k, X_i)=\lambda_k. \]

贡献率：贡献率分为两种，一种是主成分整体的贡献率，一种是主成分对某个变量的贡献率。

1. 主成分\(Z_k\)的贡献率：

   \[\frac{\lambda_k}{\lambda_1+\cdots +\lambda_k}. \]

2. 前\(m\)个主成分\(Z_1,\cdots,Z_m\)的累计贡献率：

   \[\frac{\sum_{k=1}^m\lambda_k}{\sum_{i=1}^p\lambda_i}. \]

3. 主成分\(Z_k\)对变量\(X_i\)的贡献率：

   \[\rho^2(Z_k,X_i)=\frac{\lambda_k a_{ik}^2}{\sigma_{ii}}. \]

4. 前\(m\)个主成分\(Z_1,\cdots,Z_m\)对变量\(X_i\)的贡献率：

   \[\nu_i^{(m)}=\sum_{k=1}^m \rho^2(Z_k, X_i)=\sum_{k=1}^m\frac{\lambda_ka_{ik}^2}{\sigma_{ii}}. \]

   由前面的性质，\(p\)个主成分对任何变量的贡献率都是1。

3-2 主成分分析实操

样本主成分的处理，一般会使用标准化以后的数据，此时得到的自协方差矩阵就是相关阵\(R\)，

\[R=\frac{1}{n-1}X'X. \]

得到\(R\)的按特征值大小排列的单位正交向量\(a_1,\cdots,a_p\)，\(A=(a_1,\cdots,a_p)\)。

主成分得分：第\(t\)个样品\(X^{(t)}=(x_{t1},x_{t2},\cdots,x_{tp})'\)，主成分得分向量为

\[Z_{(t)}=A'X^{(t)}=\begin{bmatrix} z_{t1} \\ z_{t2} \\ \vdots \\ z_{tp} \end{bmatrix}, \]

所有样品的主成分得分矩阵为

\[ Z=\begin{bmatrix} Z_{(1)}' \\ Z_{(2)}' \\ \vdots \\ Z_{(n)}' \end{bmatrix}=\begin{bmatrix} z_{11} & z_{12} & \cdots & z_{1p} \\ z_{21} & z_{22} & \cdots & z_{2p} \\ \vdots & \vdots & & \vdots \\ z_{n1} & z_{n2} & \cdots & z_{np} \end{bmatrix}. \]

3-3 因子模型

正交因子模型：设\(X=(X_1,\cdots,X_p)'\)是可观测的随机向量，\(\mathbb{E}(X)=\mu\)，\(\mathbb{D}(X)=\Sigma\)。设\(F=(F_1,\cdots,F_m)\)，\(m<p\)是不可观测的随机变量，称为公共因子；\(\varepsilon=(\varepsilon_1,\cdots,\varepsilon_p)'\)称为特殊因子。正交因子模型是这样的模型：

\[X=\mu+AF+\varepsilon, \]

这里

1. 公共因子满足\(\mathbb{E}(F)=0\)，\(\mathbb{D}(F)=I_m\)，即因子是标准化互不相关的。
2. 特殊因子满足\(\mathbb{E}(\varepsilon)=0\)，\(\mathbb{D}(\varepsilon)={\rm diag}(\sigma_1^2,\cdots,\sigma_p^2)\xlongequal{def}D\)，即方差是对角阵；
3. 公共因子与特殊因子互不相关，即\({\rm COV}(F, \varepsilon)=O\)。
4. \(A=(a_{ij})_{m\times p}\)是待估的系数矩阵，称为因子载荷阵。\(a_{ij}\)是变量\(X_i\)在因子\(F_j\)上的载荷，简称为因子载荷。

方差分解：可以把原始变量的方差分解为公共方差与特殊方差。

\[\mathbb{D}(X)=\Sigma=AA'+D\Rightarrow \Sigma-D=AA'. \]

载荷分析：

\[{\rm COV}(X, F)=A\Rightarrow{\rm Cov}(X_i, F_j)=a_{ij}. \]

因子载荷\(a_{ij}\)反映了变量\(X_i\)在因子\(F_j\)上的相对重要性。

1. 变量共同度：\(A\)中第\(i\)行元素的平方和称为变量\(X_i\)的共同度，即

   \[h_i^2=\sum_{j=1}^m a_{ij}^2. \]

   有\(\mathbb{D}(X_i)=h_i^2+\varepsilon_i^2\)，即分量\(X_i\)的方差可以分为公因子方差与特殊方差两部分。\(h_i^2\)描述了变量\(X_i\)对所有公因子的依赖程度，它越大，变量\(X_i\)就越能由因子描述。

2. 公因子贡献：\(A\)中第\(j\)列元素的平方和称为公因子\(F_j\)对\(X\)所有分量的总影响，即

   \[q_j^2=\sum_{i=1}^p a_{ij}^2. \]

   如果\(q_j^2\)越大，则\(F_j\)对\(X\)的贡献越大。如果把载荷矩阵每一列的\(q_j^2\)算出来，就能够找到最有贡献的因子。

载荷矩阵\(A\)的参数估计：由\(p\)个相关变量的观测数据得到协方差阵\(S\)和相关阵\(R\)。

1. 主成分解：设\(S\)的特征值为\(\lambda_1\ge\cdots\ge\lambda_p\ge0\)，相应的单位正交向量为\(l_1,\cdots,l_p\)，则\(S\)有谱分解式

   \[S=\sum_{i=1}^n \lambda_il_il_i'. \]

   如果最后几个特征值比较小，则近似有

   \[\begin{aligned} S&\approx \sum_{i=1}^m\lambda_il_il_i'+D\\ &=(\sqrt{\lambda_1}l_1,\cdots,\sqrt{\lambda_m}l_m)\begin{pmatrix} \sqrt{\lambda_1}l_1 \\ \vdots \\ \sqrt{\lambda_m}l_m \end{pmatrix}+\begin{pmatrix} \sigma_1^2 \\ & \sigma_2^2 \\ && \ddots \\ &&& \sigma_p^2 \end{pmatrix}\\ &\xlongequal{def}AA'+D,\\ \sigma_{i}^2&=s_{ii}-\sum_{t=1}^m a_{it}^2. \end{aligned} \]

   主成分解的每一列与主成分分析中使用的主成分，仅相差\(\sqrt{\lambda_i}\)倍数。事实上也经常用相关阵\(R\)代替\(S\)计算主成分解。

2. 主因子解：从相关阵\(R\)出发，需要特殊方差的初始估计量\((\hat{\sigma}_i^*)^2\)，得到初始共同度的估计量\((h_i^*)^2=1-(\hat{\sigma}_i^*)^2\)，定义约相关阵为

   \[R^*=R-D=\begin{bmatrix} (\hat{\sigma}_1^*)^2 & r_{12} & \cdots & r_{1p} \\ r_{21} & (\hat{\sigma}_2^*)^2 & \cdots & r_{2p} \\ \vdots & \vdots & & \vdots \\ r_{p1} & r_{p2} & \cdots & (\hat{\Sigma}_p^*)^2 \end{bmatrix}. \]

   从这个约相关阵出发，如同主成分解一样构造近似分解\(R^*=AA'\)，然后重新计算特殊方差：

   \[\hat{\sigma}_i^2=1-\sum_{t=1}^m a_{it}^2. \]

   再将这个方差作为特殊方差的初始估计量，反复迭代计算，直到得出稳定界。

3-4 因子模型实操

求主成分解的步骤：

1. 由样本数据阵，计算相关数字特征：\(\bar X, A, R\)。

2. 求\(R\)的特征值和标准化特征向量，特征值为\(\lambda_1\ge \cdots \ge\lambda_p\ge 0\)，对应的单位特征向量为\(l_1,\cdots,l_p\)。

3. 确定因子个数\(m\)，一般使用特征值之和不小于0.7的\(m\)作为因子个数。

4. 令\(a_i=\sqrt{\lambda_i}l_i\)，得到载荷矩阵\(A=(a_1,\cdots,a_m)\)。

5. 求特殊因子方差，得到对角阵\(D\)：

   \[\sigma_{ii}^2=1-\sum_{t=1}^m a_{it}^2. \]

6. 解释因子。

正交旋转：找到一个正交矩阵\(\Gamma_{m\times m}\)，使得旋转因子为\(Z=\Gamma' F\)，构建出的因子模型为

\[X = A\Gamma Z+\varepsilon. \]

• 目标：使得因子载荷矩阵\(A\)的方差尽可能大，即希望每一列的数值是分散的。

• 方差度量（用不上吧）：

  \[V=\sum_{j=1}^m V_m=\frac{1}{p^2}\left\{\sum_{j=1}^m\left[p\sum_{i=1}^p\frac{a_{ij}^4}{h_i^4}-\left(\sum_{t=1}^p\frac{a_{tj}^2}{h_t^2} \right)^2 \right] \right\}. \]

因子得分：把公共因子表示为变量的线性组合，或反过来对每一个样品估计公共因子的估计值，是对不可观测的随机向量的估计。现假设\(X=AF+\varepsilon\)，\(\Sigma=AA'+D\)，对于每一个观测样本，\(X\)、\(A\)和\(D\)是已知的。

1. 加权最小二乘法：因子得分为

   \[\hat{F}=(A'D^{-1}A)^{-1}(A'D^{-1})X. \]

   这样算出来的因子得分称为Bartlett因子得分。如果使用主成分法估计，因子得分常使用不加权的最小二乘法，即

   \[\hat{F}=(A'A)^{-1}A'X. \]

2. 回归法：因子得分为

   \[\hat{F}=A'R^{-1}X=A'(AA'+D)^{-1}X. \]

   这样算出来的因子得分称为Thompson因子得分。

3-5 典型相关变量求解

典型相关分析：研究两组随机变量\(X=(X_1,\cdots,X_p)'\)、\(Y=(Y_1,\cdots,Y_q)'\)之间的相关性。

1. \(p=q=1\)时，相关性为相关系数：

   \[\rho(X, Y)=\frac{{\rm Cov}(X, Y)}{\sqrt{\mathbb{D}(X)\mathbb{D}(Y)}}. \]

2. \(p>1, q=1\)时，相关性为全相关系数：

   \[R=\sqrt{\frac{\Sigma_{XY}\Sigma_{XX}^{-1}\Sigma_{YX}}{\sigma_{YY}}}. \]

3. \(p\ne 1,q\ne 1\)时，研究两个线性组合之间的相关系数：

   \[V=\boldsymbol{a}_p'X,\quad W=\boldsymbol{b}_q'Y,\\ \rho(V, W). \]

典型相关系数：设\(X=(X_1,\cdots,X_p)'\)，\(Y=(Y_1,\cdots,Y_q)'\)，将其连接为\(Z=(X', Y')'\)，\(p\le q\)，满足\(\mathbb{E}(Z)=0\)，\(\mathbb{D}(Z)=\Sigma>0\)。如果存在\(a_1=(a_{11},\cdots,a_{1p})'\)和\(b_1=(b_{11},\cdots,b_{1q})'\)，使得

\[\rho(a_1'X, b_1'Y)=\max_{\mathbb{D}(a'X)=\mathbb{D}(b'Y)=1}\rho(a'X, b'Y), \]

则称：\(a_1X, b_1Y\)是第一对典型相关变量，相关系数\(\rho(a_1X, b_1Y)\)是第一个典型相关系数。

如果存在\(a_k=(a_{k1},\cdots,a_{kp})'\)和\(b_k=(b_{k1},\cdots,b_{kq})'\)使得

1. \(a_k'X\)与\(b_k'Y\)与前面的\(k-1\)对典型相关变量都不相关；
2. \(\mathbb{D}(a_k'X)=\mathbb{D}(b_k'Y)=1\)；
3. \(a_k'X\)与\(b_k'Y\)的相关系数最大；
   则称\(a_k'X, b_k'Y\)是第\(k\)对典型相关变量，\(\rho(a_k'X, b_k'Y)\)为第\(k\)个典型相关系数。这里\(k\le p\)。

协方差阵分块：\(\mathbb{D}(Z)=\Sigma\)将分解为

\[\Sigma=\begin{bmatrix} \Sigma_{11} & \Sigma_{12} \\ \Sigma_{21} & \Sigma_{22} \end{bmatrix}. \]

求解典型相关变量：求解第一对典型相关变量是求解以下的规划问题：

\[\max a_1'\Sigma_{12} b_1,\\ \text{s.t. }a_1'\Sigma_{11} a_1=b_1'\Sigma_{22} b_2=1. \]

拉格朗日函数为

\[L(a_1,b_1,\lambda_1,\lambda_2)=a_1'\Sigma_{12}b_1-\frac{\lambda_1}{2}(a_1'\Sigma_{11}a_1-1)-\frac{\lambda_2}{2}(b_1'\Sigma_{22}b_1-1),\\ \frac{\partial L}{\partial a_1}=\Sigma_{12}b_1-\lambda_1\Sigma_{11}a_1=0,\\ \frac{\partial L}{\partial b_1}=\Sigma_{21} a_1-\lambda_2\Sigma_{22}b_1=0. \]

用\(a_1'\)，\(b_1'\)左乘以上方程，得到

\[a_1'\Sigma_{12}b_1=\lambda_1=\lambda_2\xlongequal{def}\lambda. \]

代回上式得到

\[\left\{ \begin{array}l -\lambda \Sigma_{11} a_1 + \Sigma_{12} b_1 = 0, \\ \Sigma_{21} a_1 - \lambda \Sigma_{22} b_1 = 0, \end{array} \right. \]

方程组具有非零解的条件是：

\[\begin{vmatrix} -\lambda \Sigma_{11} & \Sigma_{12} \\ \Sigma_{21} & -\lambda\Sigma_{22} \end{vmatrix}= 0. \]

为了求解此\(\lambda\)，常作以下变换：将拉格朗日偏导中的二式左乘\(\Sigma_{12}\Sigma_{22}^{-1}\)后，代入一式，得到（另一个式子同理）

\[(\Sigma_{11}^{-1}\Sigma_{12}\Sigma_{22}^{-1}\Sigma_{21}-\lambda^2 I_p)\alpha = 0.\\ (\Sigma_{22}^{-1}\Sigma_{21}\Sigma_{11}^{-1}\Sigma_{12}-\lambda^2 I_p)\beta = 0. \]

即\(\lambda^2\)是\(\Sigma_{11}^{-1}\Sigma_{12}\Sigma_{22}^{-1}\Sigma_{21}\)和\(\Sigma_{22}^{-1}\Sigma_{21}\Sigma_{11}^{-1}\Sigma_{12}\)的公共特征值，且非零特征值至多为\(p\)个。

求解定理：设\(T=\Sigma_{11}^{-1/2}\Sigma_{12}\Sigma_{22}^{-1/2}\)，\(p\)阶方阵\(TT'\)的特征值依次是\(\lambda_1^2\ge \lambda_2^2\ge\cdots\ge \lambda_p^2>0\)，对应的单位特征向量为\(l_1,\cdots,l_p\)，则第\(k\)对典型相关变量为

\[a_k=\Sigma_{11}^{-1/2}l_k,\quad b_k=\frac{1}{\lambda_k}\Sigma_{22}^{-1}\Sigma_{21} a_k. \]

如果以上任何矩阵不是可逆的，则使用其广义逆矩阵。

矩阵\(A\)的广义逆矩阵\(D\)是唯一的，满足

\[ADA=A,\quad DAD=D, \quad (AD)'=AD, \quad (DA)'=DA. \]

典型相关变量的相关性：设\(V_k, W_k\)分别是\(A, B\)的第\(k\)对典型相关变量，定义\(V=(V_1,\cdots, V_p)'\)，\(W=(W_1,\cdots,W_p)'\)，则

\[\mathbb{D}\begin{bmatrix} V \\ W \end{bmatrix}=\begin{bmatrix} I_p & \Lambda \\ \Lambda & I_p \end{bmatrix},\\ \Lambda = {\rm diag}(\lambda_1, \lambda_2, \cdots, \lambda_p). \]

典型结构：记\(A_{p\times p}=(a_1,\cdots, a_p)\)，\(B_{q\times p}=(b_1,\cdots, b_q)\)，则\(V=A'X\)，\(W=B'Y\)，有

\[{\rm COV}(X, V)={\rm COV}(X, A'X)=\Sigma_{11}A,\\ {\rm COV}(X, W)={\rm COV}(X, B'Y)=\Sigma_{12}B,\\ {\rm COV}(Y, V)={\rm COV}(Y, A'Y)=\Sigma_{21}A,\\ {\rm COV}(Y, W)={\rm COV}(Y, B'Y)=\Sigma_{22}B. \]

3-6 典型相关实操

\(\Sigma\)未知时，会使用其无偏估计\(S\)作为替代，也常常从相关阵\(R\)出发导出典型相关变量：

\[a_k=R_{11}^{-1/2} l_k, b_k=\frac{1}{\lambda_k}\Sigma_{22}^{-1}\Sigma_{21}a_k,\\ T = R_{11}^{-1/2}R_{12}R_{22}^{-1/2},\\ {\rm eigenvalue}(TT')=\lambda_k^2,\quad {\rm eigenvector}(TT')=l_k,\\ V_k=a_k'X_\text{scale}, \quad W_k=b_k'Y_\text{scale}. \]

相关假设检验：

1. \(\Sigma_{12}=O\)：如果两组样本不相关，则样本的典型相关分析没有意义。（回顾1-5节）
2. \(\lambda_k=0\)：检验是否需要使用第\(k\)组典型相关变量。

典型相关得分：将样本\(Z_{(t)}=(X_{(t)}',Y_{(t)}')'\)的值代入典型相关变量中，则第\(t\)个样本的第\(i\)对典型变量的得分为\((v_{ti}, w_{ti})\)

\[v_{ti}=a_i'(X_{(t)}-\bar X), \\ w_{ti}=b_i'(Y_{(t)}-\bar Y). \]

典型冗余：记

\[R(X, V)={\rm COV}(X, V)=R_{11}A=(r(X_i, V_j))_{p\times p}, \\ R(Y, W)={\rm COV}(Y, W)=R_{22}B=(r(Y_i, W_j))_{q\times p}, \\ R(X, W)={\rm COV}(X, W)=R_{12}B=(r(X_i, W_j))_{p\times p}, \\ R(Y, V)={\rm COV}(Y, V)=R_{21}A=(r(Y_i, V_j))_{q\times p}. \]

则称第\(k\)个典型变量\(V_k, W_k\)解释本组变量\(X, Y\)的总变差的百分比为

\[R_d(X;V_k)=\frac{1}{p}\sum_{i=1}^p r^2(X_i, V_k); \\ R_d(Y;V_k)=\frac{1}{q}\sum_{i=1}^q r^2(Y_i, W_k). \]

称前\(m\)个典型变量\(V_k, W_k\)解释本组变量\(X, Y\)的总变差累计百分比为

\[R_d(X;V_1,\cdots,V_m)=\frac{1}{p}\sum_{k=1}^m\sum_{i=1}^p r^2(X_i, V_k); \\ R_d(Y;W_1, \cdots, W_m)=\frac{1}{q}\sum_{k=1}^m\sum_{i=1}^p r^2(Y_i, W_k). \]

第\(k\)个典型变量\(V_k, W_k\)解释另一组变量\(Y, X\)的总变差的百分比为

\[R_d(X;W_k)=\frac{1}{p}\sum_{i=1}^p r^2(X_i, W_k)=\lambda_k^2R_d(X; V_k), \\ R_d(Y; V_k)=\frac{1}{q}\sum_{i=1}^q r^2(Y_i, V_k)=\lambda_k^2R_d(Y; W_k). \]

其中\(R_d(X;W_k)\)称为第一组典型变量的冗余测度，\(R_d(Y; V_k)\)称为第二组典型变量的冗余测度。冗余测度的大小表现这对典型变量能够对另一组变差相互解释的程度大小，体现了两组变量之间的相关程度。由另一个典型相关变量解释的总变差，总与自身相关变量解释的总变差距离\(\lambda_k^2\)倍。