偏相關

在機率論和統計學中，偏相關（Partial correlation）度量了兩個隨機變量在移除了其他控制變量（Confounding variables）的影響後的相關程度。

在確定兩個感興趣變量之間的數值關係時，如果存在與這兩個變量都相關的混雜變量，直接計算它們的相關係數可能會產生誤導性結果（偽相關）。通過計算偏相關係數控制混雜變量，可以避免這種誤導性資訊。這正是多元線性回歸中包含其他右側變量的動機；雖然多元回歸給出了效應值的無偏結果，但它並不直接給出一個數值來度量兩個感興趣變量之間關係的強度。

例如，給定關於個人消費、收入和財富的經濟數據，如果希望分析消費與收入之間的關係，而在計算相關係數時未能控制「財富」這一變量，可能會得到誤導性的結果。因為收入可能與財富相關，而財富又與消費相關；直接測量的消費與收入之間的相關性可能會被這些其他相關性所「污染」。使用偏相關可以避免這個問題。

與普通的相關係數一樣，偏相關係數的取值範圍也是從 -1 到 1。

- 值 -1 表示在控制某些變量後，存在完全的負相關（即一個變量值越高，另一個變量值越低）。
- 值 1 表示完全的正線性關係。
- 值 0 表示沒有線性關係。

如果隨機變量的聯合分布是多元正態分布、其他橢圓分布、多元超幾何分布、多元負超幾何分布、多項分布或狄利克雷分布，則偏相關與條件相關一致，但一般情況下兩者並不相同。

## 正式定義

形式上，給定一組 $n$ 個控制變量 $\mathbf{Z} = \{Z_1, Z_2, \dots, Z_n\}$，變量 $X$ 和 $Y$ 之間的偏相關係數 $\rho_{XY\cdot\mathbf{Z}}$ 定義為：$X$ 關於 $\mathbf{Z}$ 的線性回歸的殘差 $e_X$ 與 $Y$ 關於 $\mathbf{Z}$ 的線性回歸的殘差 $e_Y$ 之間的相關係數。

一階偏相關（即 $n=1$ 時）是相關係數與可移除相關係數之積的差，除以可移除相關係數的相異係數（coefficient of alienation）之積。

## 計算方法

### 使用線性回歸

計算樣本偏相關的一種簡單方法是求解兩個相關的線性回歸問題，並計算殘差之間的相關性。

設 $X$ 和 $Y$ 為取實數值的隨機變量，$\mathbf{Z}$ 為 $n$ 維向量值隨機變量。設 $x_i, y_i, \mathbf{z}_i$ 表示從實隨機變量 $X, Y, \mathbf{Z}$ 的聯合分布中抽取的 $N$ 個獨立同分布（i.i.d.）觀測值中的第 $i$ 個（其中 $\mathbf{z}_i$ 增加了一個 1 以允許回歸中的常數項）。

求解線性回歸問題相當於找到 $(n+1)$ 維回歸係數向量 $\mathbf{w}_X^*$ 和 $\mathbf{w}_Y^*$，使得：

$$\mathbf{w}_X^* = \arg\min_{\mathbf{w}} \left\{ \sum_{i=1}^N (x_i - \langle\mathbf{w}, \mathbf{z}_i \rangle)^2 \right\}$$

$$\mathbf{w}_Y^* = \arg\min_{\mathbf{w}} \left\{ \sum_{i=1}^N (y_i - \langle\mathbf{w}, \mathbf{z}_i \rangle)^2 \right\}$$

其中 $N$ 是觀測數量，$\langle\mathbf{w}, \mathbf{z}_i \rangle$ 是向量 $\mathbf{w}$ 和 $\mathbf{z}_i$ 之間的標量積。

殘差為：

$$e_{X,i} = x_i - \langle\mathbf{w}_X^*,\mathbf{z}_i \rangle$$

$$e_{Y,i} = y_i - \langle\mathbf{w}_Y^*,\mathbf{z}_i \rangle$$

樣本偏相關係數即為這些新生成的殘差值之間的樣本相關係數：

$$\begin{align} \hat{\rho}_{XY\cdot\mathbf{Z}}&=\frac{N\sum_{i=1}^N e_{X,i}e_{Y,i}-\sum_{i=1}^N e_{X,i}\sum_{i=1}^N e_{Y,i}} {\sqrt{N\sum_{i=1}^N e_{X,i}^2-\left(\sum_{i=1}^N e_{X,i}\right)^2}~\sqrt{N\sum_{i=1}^N e_{Y,i}^2-\left(\sum_{i=1}^N e_{Y,i}\right)^2}}\\ &=\frac{N\sum_{i=1}^N e_{X,i}e_{Y,i}} {\sqrt{N\sum_{i=1}^N e_{X,i}^2}~\sqrt{N\sum_{i=1}^N e_{Y,i}^2}}. \end{align}$$

在第一個表達式中，減號後的三項都等於 0，因為它們包含普通最小二乘法回歸的殘差和（必定為 0）。

#### 示例

考慮三個變量 $X, Y, Z$ 的以下數據：

X｜Y｜Z

2｜1｜0

4｜2｜0

15｜3｜1

20｜4｜1

計算變量 $X$ 和 $Y$ 之間的皮爾遜相關係數約為 0.970；而使用上述公式計算 $X$ 和 $Y$ 之間的偏相關係數（控制 $Z$）則為 0.919。以下是使用 R語言 進行的計算：

```R
> x <- c(2, 4, 15, 20)
> y <- c(1, 2,  3,  4)
> z <- c(0, 0,  1,  1)

# 将 x 对 z 进行回归并计算残差
> res_x <- lm(x ~ z)$residuals

# 将 y 对 z 进行回归并计算残差
> res_y <- lm(y ~ z)$residuals

# 计算残差的相关性
> cor(res_x, res_y)
# [1] 0.919145

# 显示这与 x 和 y 之间的直接相关性不同
> cor(x, y)
# [1] 0.9695016
```

### 使用遞推公式

求解線性回歸問題在計算上可能非常昂貴。實際上，$n$ 階偏相關（即 $|\mathbf{Z}| = n$）可以很容易地從三個 $(n-1)$ 階偏相關計算得出。零階偏相關 $\rho_{XY\cdot\varnothing}$ 定義為普通的相關係數 $\rho_{XY}$。

對於任意 $Z_0 \in \mathbf{Z}$，有：

$$\rho_{XY\cdot \mathbf{Z} } = <pre> \frac{\rho_{XY\cdot\mathbf{Z}\setminus\{Z_0\}} - \rho_{XZ_0\cdot\mathbf{Z}\setminus\{Z_0\}}\rho_{Z_0Y\cdot\mathbf{Z}\setminus\{Z_0\}}} {\sqrt{1-\rho_{XZ_0\cdot\mathbf{Z}\setminus\{Z_0\}}^2} \sqrt{1-\rho_{Z_0Y\cdot\mathbf{Z}\setminus\{Z_0\}}^2}}$$

如果單純地將此計算實現為遞歸演算法，會產生指數級的時間複雜度。然而，由於該計算具有重疊子問題性質，使用動態規劃或簡單地緩存遞歸調用的結果可以產生 $\mathcal{O}(n^3)$ 的複雜度。

特別地，當 $Z$ 是單個變量時，公式簡化為：

$$\rho_{XY\cdot Z}= <pre> \frac{\rho_{XY}-\rho_{XZ}\rho_{ZY}} {\sqrt{1-\rho_{XZ}^2}\sqrt{1-\rho_{ZY}^2}}$$

### 使用矩陣求逆

偏相關也可以用聯合精度矩陣（Precision matrix）來表示。考慮一組隨機變量 $\mathbf{V} = \{X_1,\dots, X_n\}$。我們想要計算在給定其他所有變量 $\mathbf{V} \setminus \{X_i,X_j\}$ 的情況下，變量 $X_i$ 和 $X_j$ 之間的偏相關。

假設（聯合/完整）協方差矩陣 $\Sigma = (\sigma_{ij})$ 是正定矩陣，因此是可逆矩陣。定義精度矩陣為 $\Omega = (p_{ij}) = \Sigma^{-1}$，則：

$$\rho_{X_i X_j\cdot \mathbf{V} \setminus \{X_i,X_j\}} = - \frac{p_{ij}}{\sqrt{p_{ii}p_{jj}}}$$

計算此式需要求協方差矩陣 $\Sigma$ 的逆 $\Sigma^{-1}$，時間複雜度為 $\mathcal{O}(n^3)$。值得注意的是，只需一次矩陣求逆即可獲得 $\mathbf{V}$ 中所有變量對之間的偏相關係數。

## 解釋

### 幾何解釋

<image>

設三個變量 $X, Y, Z$（其中 $Z$ 是「控制」或「額外變量」）是從 $n$ 個變量 $\mathbf{V}$ 的聯合機率分布中選出的。進一步設 $\mathbf{v}_i, 1 \le i \le N$ 是從 $\mathbf{V}$ 的聯合機率分布中抽取的 $N$ 個 $n$ 維獨立同分布觀測值。

考慮由觀測值形成的 $N$ 維向量 $\mathbf{x}$、$\mathbf{y}$ 和 $\mathbf{z}$。

可以證明，來自 $X$ 對 $\mathbf{Z}$ 的線性回歸的殘差 $e_{X,i}$，如果視為一個 $N$ 維向量 $\mathbf{e}_X$（在附圖中記為 $\mathbf{r}_X$），與 $\mathbf{Z}$ 生成的向量 $\mathbf{z}$ 的標量積為零。這意味著殘差向量位於垂直於 $\mathbf{z}$ 的 $(N-1)$ 維超平面 $S_{\mathbf{z}}$ 上。

同理也適用於 $Y$ 的殘差向量 $\mathbf{e}_Y$。所需的偏相關就是向量 $\mathbf{x}$ 和 $\mathbf{y}$ 在垂直於 $\mathbf{z}$ 的超平面上的投影 $\mathbf{e}_X$ 和 $\mathbf{e}_Y$ 之間夾角 $\phi$ 的餘弦。

### 作為條件獨立性測試

假設所有涉及的變量都是多元高斯分布的，偏相關 $\rho_{XY\cdot\mathbf{Z}}$ 為零若且唯若 $X$ 在給定 $\mathbf{Z}$ 的條件下與 $Y$ 條件獨立。這一性質在一般情況下不成立。

為了檢驗樣本偏相關 $\hat{\rho}_{XY\cdot\mathbf{Z}}$ 是否意味著總體偏相關不為 0，可以使用偏相關的費希爾變換（Fisher's z-transform）：

$$z(\hat{\rho}_{XY\cdot\mathbf{Z}}) = \frac{1}{2} \ln\left(\frac{1+\hat{\rho}_{XY\cdot\mathbf{Z}}}{1-\hat{\rho}_{XY\cdot\mathbf{Z}}}\right)$$

零假設是 $H_0: \rho_{XY\cdot\mathbf{Z}} = 0$，針對雙尾備擇假設 $H_A: \rho_{XY\cdot\mathbf{Z}} \neq 0$。如果下式成立，則拒絕 $H_0$：

$$\sqrt{N - |\mathbf{Z}| - 3}\cdot |z(\hat{\rho}_{XY\cdot\mathbf{Z}})| > \Phi^{-1}(1-\alpha/2)$$

其中 $\Phi$ 是具有零均值和單位標準差的高斯分布的累積分布函數，$\alpha$ 是顯著性水平，$N$ 是樣本大小。

## 半偏相關

半偏相關（Semipartial correlation，或 part correlation）統計量與偏相關統計量類似；兩者都比較兩個變量在控制某些因素後的變異。然而，計算半偏相關時，第三個變量僅針對 $X$ 或 $Y$ 中的一個保持常數，而不是同時針對兩者。偏相關比較的是一個變量的唯一變異與另一個變量的唯一變異，而半偏相關比較的是一個變量的唯一變異（移除了與 $Z$ 相關的變異）與另一個變量的未過濾變異。

半偏相關可以被視為更具有實際相關性，因為它相對於因變量（響應變量）的總變異進行了縮放。相反，它在理論上的用處較小，因為它對自變量的唯一貢獻的作用不那麼精確。

$X$ 與 $Y$ 的半偏相關的絕對值總是小於或等於 $X$ 與 $Y$ 的偏相關。原因如下：假設從 $X$ 中移除了 $X$ 與 $Z$ 的相關性，得到殘差向量 $e_x$。在計算半偏相關時，$Y$ 仍然包含唯一方差和由於其與 $Z$ 關聯而產生的方差。但是 $e_x$ 由於與 $Z$ 不相關，只能解釋 $Y$ 方差中唯一的一部分，而不能解釋與 $Z$ 相關的那部分。相比之下，對於偏相關，只有 $e_y$（$Y$ 方差中與 $Z$ 無關的部分）需要被解釋，因此原本 $e_x$ 無法解釋的那類方差在 $e_y$ 中已經不存在了。

## 時間序列分析中的應用

在時間序列分析中，時間序列的偏自相關函數（Partial autocorrelation function, PACF，有時也稱為「偏相關函數」）定義為滯後 $h$ 時的偏相關：

$$\varphi(h)= \rho_{X_0X_h\,\cdot\, \{X_1,\,\dots\,,X_{h-1} \}}$$

該函數用於確定自回歸（AR）過程的適當滯後長度。

## 帶有收縮估計的偏相關

當樣本量小於變量數量（即高維數據設置）時，估計偏相關具有挑戰性。在這種情況下，樣本協方差 $\hat{\Sigma}$ 是病態的，求其逆矩陣 $\hat{\Omega}$ 會出現問題。

收縮估計（Shrinkage estimation）方法可以改善 $\hat{\Sigma}$ 或 $\hat{\Omega}$ 的估計，從而產生更可靠的偏相關估計。一個例子是 Ledoit-Wolf 收縮估計量：

$$\hat{\Sigma}^{[\lambda]} = \lambda T + (1 - \lambda) \Sigma$$

其中 $\hat{\Sigma}$ 是樣本協方差矩陣，$T$ 是目標矩陣（例如對角矩陣），收縮強度 $\lambda\in (0,1)$。

Ledoit-Wolf 收縮下的偏相關估計為：

$$\hat{P}_{ij}^{[\lambda]} = \frac{\hat{\Omega}_{ij}^{[\lambda]}}{\sqrt{\hat{\Omega}_{ii}^{[\lambda]} \hat{\Omega}_{jj}^{[\lambda]} }} <pre>$$

其中 $\hat{\Omega}_{ij}^{[\lambda]}$ 是 $\hat{\Sigma}_{ij}^{[\lambda]}$ 的逆矩陣元素。這種方法應用於金融和基因組學等多個領域。

## 參見

- 線性回歸
- 條件獨立
- 復相關 (Multiple correlation)
- 偏自相關函數 (Partial autocorrelation function)