3. PCoA原理详解：距离矩阵、特征分解与坐标变换

这一章，我们来啃硬骨头——PCoA的数学原理。很多教程上来就甩公式，看得人头晕。我换个讲法，先告诉你PCoA到底在干什么。

说白了，PCoA就是给数据找个“能看得见”的投影。跟PCA不一样的是，它不关心变量，只关心样本之间的距离。距离算好了，剩下的就是特征分解的活儿。

3.1 距离矩阵——PCoA的起点

PCoA的第一步，是先把原始数据变成距离矩阵。这一步很关键，你选什么距离，结果就偏向什么。

怎么理解距离矩阵？

假设你有n个样本。距离矩阵D是一个n×n的方阵，第i行第j列的元素d_ij表示样本i和样本j之间的“不相似度”。对角线上全是0——自己跟自己当然没距离。

常用的距离度量

3.2 核心推导——从距离到内积

拿到距离矩阵D之后，我们不能直接用，因为D不是内积矩阵。PCoA的核心思想是：在欧氏空间中，距离和内积是可以互相转换的。

这一转换需要两步：

1. 平方距离矩阵：对每个元素取平方，得到D²。
2. 双中心化（double centering）：减去行均值、列均值，再加上总均值。

公式长这样：

B = -½ J D² J

其中 J = I - (1/n) · 1·1ᵀ
     I 是n维单位矩阵
     1 是n维全1列向量

这个B矩阵，就是我们要的内积矩阵。它的第i行第j列元素，等于样本i和样本j在PCoA空间里的内积。

3.3 特征分解——找主方向

有了B矩阵，接下来就熟悉了：特征分解。

对B做特征分解：

B = V Λ Vᵀ

其中 Λ = diag(λ₁, λ₂, ..., λₙ)
     V = [v₁, v₂, ..., vₙ] 是特征向量矩阵

特征值λ₁ ≥ λ₂ ≥ ... ≥ λₙ，对应的特征向量就是PCoA的主轴方向。

特征值有什么用？

• 特征值大小代表该轴对总方差的贡献
• 特征值占比 = 该轴解释的方差比例
• 前两个轴加起来能解释多少？这个值越大，降维效果越好

3.4 坐标变换——从高维到低维

特征分解之后，我怎么拿到样本在低维空间的坐标？

很简单。样本在第k个主轴的坐标就是：

坐标矩阵 = V · diag(√λ₁, √λ₂, ..., √λₖ)

具体来说，样本i在第k轴上的坐标 = λₖ · vₖ[i]

这里vₖ[i]是第k个特征向量的第i个分量。为什么开根号？因为我们要从内积空间回到欧氏空间，开根号保证坐标的尺度是对的。

通常我们只取前2-3个轴来画图。为什么不取更多？因为超过3维就看不懂了，可视化意义不大。

用R做这一套，其实很简单：

# R代码示例
library(vegan)

# 计算Bray-Curtis距离
dist_mat <- vegdist(data, method = "bray")

# 执行PCoA（也叫MDS）
pcoa_res <- cmdscale(dist_mat, k = 2, eig = TRUE)

# 提取坐标
scores <- as.data.frame(pcoa_res$points)
colnames(scores) <- c("PCoA1", "PCoA2")

Python里也不复杂：

# Python代码示例
from skbio.stats.ordination import pcoa
from skbio.diversity import beta_diversity

# 计算Bray-Curtis距离矩阵
dm = beta_diversity("bray_curtis", data, ids)

# 执行PCoA
pcoa_results = pcoa(dm)

# 获取前两轴坐标
scores = pcoa_results.samples.iloc[:, :2]

3.5 实际应用中的几个坑

实战下来，有几点我想提醒大家注意：

1. 负特征值：如果距离不是欧氏的，特征分解可能出现负特征值。解决办法：加一个足够大的常数让矩阵变成正定。
2. 解释率别太当真：PCoA前两轴的解释率通常不如PCA高，特别是用非欧氏距离的时候。我见过很多人盯着20%的解释率不放，其实没必要。
3. 距离选择影响太大：同一个数据，用欧氏距离和Bray–Curtis距离，出来的图可能完全不一样。我个人习惯先试两三种距离，看看哪种能体现已知的生物学分组。