3. CCA（典型相关分析）原理：数学基础、目标函数、迭代过程
好，咱们开始聊CCA。很多人一听“典型相关分析”就觉得头大，觉得这名字太唬人。其实说白了，它要干的事很单纯：在两个数据集里，分别找到一组“投影方向”，让投影之后的两组数据相关性最大。
举个直观的例子。你手头有同一个组织的RNA-seq和ATAC-seq数据，一个测基因表达，一个测染色质开放性。虽然测的东西不一样，但它们描述的是同一个生物学过程。CCA要做的，就是在这两组数据里各找一根“翻译尺子”，把两个模态的信息拉到同一个坐标系下，然后告诉你：“看，这两个数据在后面前几维是高度协同的。”
我早期做CITE-seq分析时就遇到过这问题——蛋白数据和mRNA数据尺度差太多，硬拼在一起完全不对。后来用了CCA，才找到它们共享的那个“生物学信号”空间。嗯，这就是它最有魅力的地方。
3.1 数学基础：从协方差到典型相关
要理解CCA，咱们得从最基本的统计量开始。
假设我们有两个矩阵X（n×p）和Y（n×q）。n是细胞数，p和q分别是各自的特征数（比如基因数）。
CCA的核心是三个协方差矩阵：
- Σ_xx：X内部的协方差矩阵
- Σ_yy：Y内部的协方差矩阵
- Σ_xy：X和Y之间的交叉协方差矩阵
你想想看，如果Σ_xy全是0，那X和Y就彻底没关系，CCA也没必要做了。CCA要找的，就是让这个“跨模态相关性”最大化的投影。
3.2 目标函数：最大化相关系数
好，目标函数来了。CCA想找到一对向量 a（p维）和 b（q维），使得投影后的变量 u = Xa 和 v = Yb 的相关系数最大。
数学上，这个相关系数可以写成：
ρ = a' Σ_xy b / sqrt( (a' Σ_xx a) * (b' Σ_yy b) )注意，分子是a和b通过交叉协方差矩阵的“耦合”，分母是两个投影变量的方差的乘积。这个形式其实跟Pearson相关系数一模一样，只是把变量换成了线性组合。
CCA最妙的约束是：要求a' Σ_xx a = 1 且 b' Σ_yy b = 1。为啥？因为如果不约束，你可以把a和b同时乘一个大数，让相关系数看起来很大，但这没有意义。约束之后，问题就变成了一个有条件的优化问题。
最终，CCA可以转化为求解以下两个耦合的特征方程：
Σ_xx^{-1} Σ_xy Σ_yy^{-1} Σ_yx a = λ a
Σ_yy^{-1} Σ_yx Σ_xx^{-1} Σ_xy b = λ b我建议你记住这个形式：CCA = 广义特征值分解。很多工具包（包括Seurat里的RunCCA）底层就是这么算的。
3.3 迭代过程：不止一对
那为什么说CCA是“迭代”的呢？因为一组a和b不够用，我们需要多对投影方向。
第一对投影方向（a¹, b¹）让u¹和v¹的相关系数最大。那么第二对呢？
要求在与第一对正交的子空间里，再找一对（a², b²），使得u²和v²的相关系数次大。以此类推。
这个过程在数学上等价于：
- 求解特征矩阵：计算M = Σ_xx^{-1} Σ_xy Σ_yy^{-1} Σ_yx
- 特征分解：对M做特征分解，得到特征值λ₁ ≥ λ₂ ≥ ... ≥ λ_k和对应的特征向量a¹, a², ..., aᵏ
- 导出b：对每个aⁱ，bⁱ = (1/√λᵢ) Σ_yy^{-1} Σ_yx aⁱ
- 截断：通常只保留前K对（比如K=15或30），因为后面的典型相关系数σ已经很小，没有生物学意义了
我在项目中遇到过一件事：有时K取太大，引入了噪声，批次校正反而变差了。我后来一般会画一个“典型相关系数碎石图”，选拐点之前的维度。说白了，跟PCA选主成分一个道理。
3.4 CCA在单细胞批次校正中的角色
你可能会问：CCA原理我们懂了，它怎么用到批次校正里？
其实很简单：
- 视批次为两个“模态”：批次1的数据看做X，批次2的数据看做Y
- 共享的生物学信号：用CCA找到的前K对典型变量(u¹, v¹)、(u², v²)等，就是跨批次共享的“低维嵌入”
- 对齐后校正：在CCA空间里，批次1和批次2对应的细胞（同一类型）应该落在相近的位置。然后可以用这个低维空间做匹配，或者直接投影回原空间做校正
嗯，这就是Seurat里经典的“CCA+MNN”策略的核心——先用CCA找共享空间，再用MNN做精细的对齐匹配。
3.5 一张图总结CCA框架
3.6 代码层面看CCA
理论讲完了，我们看看代码里怎么玩。其实核心计算就几行：
# 假设 X 和 Y 已标准化（n×p, n×q）
# 计算协方差
Sxx <- cov(X)
Syy <- cov(Y)
Sxy <- cov(X, Y)
# 组装广义特征值问题的核心矩阵
M <- solve(Sxx) %*% Sxy %*% solve(Syy) %*% t(Sxy)
# 特征分解
eig <- eigen(M)
a_vectors <- eig$vectors[, 1:K] # 前K个投影向量（对应于X）
eigenvalues <- eig$values[1:K]
# 导出Y的投影向量
b_vectors <- solve(Syy) %*% t(Sxy) %*% a_vectors %*% diag(1/sqrt(eigenvalues))
# 投影得到典型变量
U <- X %*% a_vectors # 批次1的典型变量
V <- Y %*% b_vectors # 批次2的典型变量嗯，你可能会觉得这代码看着有点“重”。别担心，实际项目里没人手写这些——Seurat的RunCCA()、scipy的sklearn.cross_decomposition.CCA都封装好了。但理解底层逻辑能帮你调参，比如知道为什么lambda参数能防止过拟合，或者为什么K不能取太大。