梯度计算方法

传统脑科学研究关注离散的脑区和模块但近年来越来越多研究者开始将脑网络映射到一个连续的低维流形表示也就是梯度 (gradient)。人脑功能网络中有两个梯度轴一个是从初级皮层到跨模态皮层一个是从视觉皮层到感觉运动皮层 (Margulies et al., 2016)。梯度的计算 (Wael et al., 2020)1.数据一个n*p的矩阵其中n是脑区数p是每个脑区的特征2.构建一个相似性矩阵相关性可以通过高斯核余弦相似性皮尔逊相关和斯皮尔曼相关等。2.降维所谓的梯度其实就是从高维(p)中提取代表性的低维嵌入所以需要进行降维。有三种方式线性降维使用PCA非线性降维使用拉普拉斯特征映射 (Laplacian Eigenmap, LM) 以及扩散映射 (Diffusion Maps, DM)。PCA见帖子脑网络分析基础3A——PCA_脑网络基础学习-CSDN博客拉普拉斯特征映射首先构建图拉普拉斯矩阵L D-AA是相似性矩阵D是相似性矩阵的对角矩阵。当然L也可以是归一化拉普拉斯矩阵Lnorm D^(1/2)LD(1/2)。然后求解广义特征值问题得到其中g是特征向量λ是特征值 (Ghojogh et al., 2023)。这个公式怎么理解首先基于拉普拉斯矩阵或者归一化拉普拉斯矩阵可以得到下面的评价图的平滑度的公式如果让这个公式最小化则高相关的两个脑区a和b (Aa,b较大) 的特征x(a)和x(b)越接近图越平滑。我们要得到合适的低维嵌入就需要得到是x^TLx最小化的x。同时不能让x的尺度任意变化于是加上一个约束x^TDx 1用拉普拉斯乘子法求极值点 (Ghojogh et al., 2023)得到进一步又gDg1所以gLg λ也就是λ越小图越平滑。这里要注意第一个特征向量特征值不能要因为此时g 0λ也等于0。除了令g0的最小的特征值对应令图最平滑的特征向量以此类推。扩散映射 (diffusion mapping, DM)首先构建Wα矩阵即Wα D^(-α) A D^(-α)其中α属于[0,1]一般取0.5类似于构建广播矩阵。然后构建扩散算子Pα Dα^(-1) Wα其中Dα^(-1)是Wα的度或说强度矩阵这里也就是给Wα做行归一化。再之后对Pα进行特征值分解。去掉最大的特征值对应特征向量因为特征向量所有元素全为1最大特征值也为1。得到特征向量g。g的每一个元素g_i的意义是从节点 i 随机扩散一步后所到达节点的 g 值期望为λg_iλ为特征值。如果λ越接近最大值1自然g向量所对应的图是更平滑的因为每个元素i扩散一步到的地方和自身相差无几。最终使用λ^t g作为梯度向量即为什么要使用λ^t g作为梯度向量首先P gk​λk​ gk然后有P^tgk (λk)^t gk。P^tgk就是gk元素扩散t步所到达的地方。t也就为控制扩散时间的参数。梯度对齐对不同的被试计算得到的梯度可能不能用于直接比较因为可能不同被试的梯度向量的符号是相反的以及可能存在多个重复的特征值需要进行对齐。第一种方法Procrustes analysis给定一个源梯度Gs和目标梯度GtProcrustes 分析旨在寻找一个正交线性变换ψ将源梯度空间 Gs​ 对齐到目标梯度空间 Gt​。在估计对齐变换之前还可以先对数据进行中心化centering和归一化normalizing从而进一步处理两个梯度空间之间的位置差异和数值大小差异。这个方法是一个迭代的方法第一次迭代把所有被试的梯度对齐到参考梯度之后的迭代把所有被试对齐之后的梯度平均作为参考梯度再进行对齐。具体的对齐方法首先构建T^T * X矩阵其中T是参考梯度矩阵X是要对齐到参考梯度的梯度矩阵。然后对这个矩阵进行奇异值分解得到V和U然后要对齐的梯度矩阵右乘V * U实现对齐。第二种方法联合表示首先构建一个大相似性矩阵对角线上的A矩阵是每个被试内部的相似性矩阵Aij就是i被试和j被试之间的相似性。然后对这个大相似性矩阵求梯度这样不同被试的梯度在同一个空间中。检验梯度的显著性研究者经常会把梯度和其它数据例如皮层厚度基因表达强度等做相关或去回归。但因为许多脑影像模态都具有明显的空间自相关spatial autocorrelation直接使用线性回归或类似统计方法可能会得到存在偏差的检验统计量。也就是哪怕不存在什么真正的相关由于空间自相关性也会产生虚假相关。所以使用spin test如果数据能映射到球面上就将皮层表面映射到球面然后随机旋转这个球面从而重新分配脑区位置上的数值或者Moran spectral randomization (MSR)这样可以在不破坏空间相关性的前提下构建零模型然后检验真实的梯度和其它模态属性相关是否真的显著这样得到的显著性就不能被空间自相关所解释。梯度的计算对齐和显著性检验可以使用brain space这个工具箱。