一、浅谈功能脑网络

首先，甲学员从他人那里获取了每个脑区的信号序列。其次，计算任意两个信号序列间的相关（皮尔逊相关）。这样，把脑区视为节点，相关值视为边，功能网络就构建好啦！

红色脑区的信号序列（均值）和蓝色脑区的信号序列（均值）做一下相关，就是红区和蓝区间的FC值，也就是红蓝脑区间的连边值。

这样，包含两个节点（红色脑区和蓝色脑区）和一个连边的最小“脑网络”已经构建好了。当我们把任意两个脑区信号序列的相关计算出来后，一个完整的脑网络就计算好了！

二、不同模态脑网络的构建

不同模态脑网络有五种构建方式：功能脑网络、结构脑网络、白质纤维束网络、加权脑网络、二值网络

功能脑网络

（1）功能脑网络：首先，获取了每个脑区的信号序列。其次，计算任意两个信号序列间的相关（皮尔逊相关）。把脑区视为节点，相关值视为边，连接边和节点就构造好了功能网络

结构脑网络

（2）结构脑网络的2种构建方式：

方式1：假定我们有100名被试，计算出红色脑区和蓝色脑区的平均灰质密度值（GMV），得到了100个红区的GMV和100个蓝区的GMV。把红区的100个GMV值和蓝区的100个GMV值做相关，就得到了红区和蓝区的相关值，也就是红区和蓝区的连边。当我们把所有脑区间的GMV相关计算出来后，就得到了所谓的结构共变脑网络。

方式2：假设每个被试都每隔三个月扫描一次T1像，至今扫描了12次。那么对每个被试的红蓝脑区，可分别计算出12个GMV值。对红蓝脑区的12个GMV值做相关即得到红蓝脑区间的连接值。

白质纤维束脑网络

（3）白质纤维束脑网络：首先，基于DTI成像，使用纤维追踪技术（以确定性纤维追踪为例）可以追踪出两个脑区间的纤维束。“如果两个脑区间存在纤维束，就认为它们之间存在连接，连接强度即为纤维束数量。也就是说，在确定性纤维追踪中，以脑区为节点，纤维束数量值为边构建白质脑网络。同样地，每个人都可以得到一个纤维束脑网络。

加权网络

（4）加权网络：每个脑网络都可以表示为邻接矩阵的形式，横轴和纵轴都是脑区编号，横纵交叉的地方就是相应脑区间的连接值，这种网络就是加权网络。”

二值网络

（5）二值网络：设定一个阈值，当连接值大于这个阈值时，就视为1（有连接）；小于这个阈值时，就视为0（无连接），这样就得到一个只区分有连接或者无连接的二值网络。其实，我们很多脑网络指标（出入度、局部效率等）都是在二值网络上计算的

（只不过把time points换成了subject）

“你看，功能脑网络既然你已经会了，现在看一下结构脑网络。“

• 假定我们有100名被试，并计算出了红色脑区和蓝色脑区的平均灰质密度值（GMV），得到了100个红区的GMV和100个蓝区的GMV。

• 同样地，把红区的100个GMV值和蓝区的100个GMV值做相关，就得到了红区和蓝区的相关值，也就是红区和蓝区的连边。

• 当我们把所有脑区间的GMV相关计算出来后，就得到了所谓的结构共变脑网络。只不过，这个脑网络一组被试才能构建一个
  

（同一个被试，每隔三个月扫描一次T1像）

“当然啦，结构脑网络也可以仿照功能脑网络的方式构建。看我画的第二张图，

• 假设每个被试都每隔三个月扫描一次T1像，至今扫描了12次。那么对每个被试的红蓝脑区，可分别计算出12个GMV值，
• 对红蓝脑区的12个GMV值做相关即得到红蓝脑区间的连接值。
• 通过这种方式，每个被试都可以得到对应的结构脑网络，是不是很类似于功能脑网络？”

（红蓝脑区间追踪出了5条纤维束，依然是荡漾的画风）

现在看一下白质纤维束脑网络。首先，基于DTI成像，使用纤维追踪技术（以确定性纤维追踪为例）可以追踪出两个脑区间的纤维束。”

（红蓝脑区间的连接值为5）

“如果两个脑区间存在纤维束，就认为它们之间存在连接，连接强度即为纤维束数量。也就是说，在确定性纤维追踪中，我们可以以脑区为节点，纤维束数量值为边构建白质脑网络。同样地，每个人都可以得到一个纤维束脑网络，是不是也很简单？”

（加权脑网络的邻接矩阵表示）

“至此，你已经知道了4种脑网络构建方式。这里我再补充一个概念：加权网络和二值网络。每个脑网络都可以表示为邻接矩阵的形式，参考上图，横轴和纵轴都是脑区编号，横纵交叉的地方就是相应脑区间的连接值（颜色随连接值大小变化），这种网络就是加权网络。”

（二值网络的邻接矩阵表示）

“我们可以设定一个阈值，当连接值大于这个阈值时，就视为1（有连接）；小于这个阈值时，就视为0（无连接），这样就得到一个只区分有连接或者无连接的二值网络。其实，我们很多脑网络指标（出入度、局部效率等）都是在二值网络上计算的…”

三、趣谈散点图与相关系数

说起相关系数，从字面上的含义就可看出，就是两个信号之间的相关性。但是你真正理解内在的机理吗？

结论放在最前面：相关系数，其实就是通过散点图来的。

所有的一切，由这个图说起：

图1：Ref: JamesAH, BMJ, 1995, 311: 1668.

有一个人，他测量了一组人的“量表”。其中这个“量表”包含着年龄和耳朵长度。这样子他就得到一个二维小表格如下图示：

图2：萌萌哒的二维量表小表格

然后他在坐标纸上面进行打点，X轴坐标设置为年龄，y轴坐标设置为耳朵长度。然后每一行就是一个点，也就是说：每一个点对应着一个被试信息。

紧接着，他就拿手来比划，画出一根能最好拟合这个散点趋势的线（拟合或最小二乘法）。这样他就发现：年龄越大，耳朵越长。Ps：怪不得如来佛耳朵如此长，连起来可以绕地球一圈。

图3：散点图与拟合线（橙色），左：正相关；中：不相关；右：负相关

其实我告诉你，现在这根橙色拟合线的趋势就是相关性。如果这根线是朝着右上角走，就是正相关；如果这根线是朝着右下角走，就是负相关；如果这根线水平，就代表着不相关。

但是理想很美好，现实很残酷。真正拿到数据进行计算相关系数，多多少少会存在一定的相关性，真正不相关的例子太少太少。前一阵子有一篇文章说：中国三峡大坝是影响日本地震的原因。该文说这个相关性还是非常非常显著的。

那么问题来了：相关系数的计算怎么会有显著性呢？

多图警示！

我做了一堆图，上面的这个例子取得不够恰当，我应该让这些图拟合的斜率是一致的（相关系数一致），但是不要在意这些细节。小伙伴们有没有看到图上的P值，那个p值就是相关的显著性。我们可以很轻松的发现：只要散点的点，越靠近拟合曲线，那么显著性越强。散点越分散，显著性越差。

以上几点细节部分特此做一个说明：

1、画散点图的时候，有白点，有黑点。图中的白点是剔除的，黑点是选取的。

2、剔除点是根据三倍方差以上的点进行剔除（图中色彩斑斓的圈，就是n倍方差边界线）

3、相关系数和斜率存在一定关系，可以说，斜率越靠近1，相关系数越大。（相关系数取值-1~1），其中0代表不相关，1代表正相关，-1代表负相关。

4、最佳拟合线的画法：在这里三种方法是等价关系，得到的数值是一样的。（最小二乘法 = 一次拟合 = 一次回归）

5、本文所指相关，指代皮尔逊相关

其中：

皮尔逊相关ρ，协方差Cov(x,y)，标准差σx σy

注：在公式内并无显著性水平计算，显著性解释是作者领悟的。在matlab中，计算相关系数是有显著性输出的。此显著性并未通过多重比较较正。有关校正，敬请看后面推送。

现在说了这么多，让我来告诉你，一些在脑科学领域用散点图来解释的本质：

1、功能连接：功能连接最早的定义就是皮尔逊相关，而功能连接就是两个脑区时间点的散点图

2、结构上的协变连接：协变连接是用得最早的，在磁共振出现之前，前人研究PET就是采用协变连接。简单说，就是A、B两个脑区之间的散点图。

3、回归：有没有发现回归也是这样子的？

四、脑电信号频域变换

脑电信号频域变换的原理

在脑电数据处理中，我们经常要做信号的时、频域变换。许多初学者已经做过一些时、频域分析，但是仍然对时、频域信号没有一个直观的理解。那么，究竟什么是时域信号？什么又是频域信号呢？

首先，理解时域信号很简单。时域信号是什么？就是以时间为横轴，数据值为纵轴的信号呀。比如这样的：

对于我们的脑电信号，我们看到的每个通道的脑电波形就是时域信号。我们经常听到的事件相关电位（ERP），也是时域信号，只不过是在某个事件（比如实验刺激）发生后，脑电信号会出现一个波形罢了。比如这样的：

在事件（0秒处）发生的300ms后，有一个正的波形，这个可能就是P300成分（Positive， 300ms）。我们把300ms称为潜伏期，波形的高度，即350uV，称为幅值。

其实，在现实世界中，我们收集到的各种各样的数据，大都是时域信号。直到一位大神出现，他就是傅里叶！

傅大神发现了什么呢？他发现周期函数，像下图这样的，都可以用一系列的正弦（和余弦）波叠加近似地表示出来。想象一下，一个方波，它可以由一个大波浪，加一个小波浪，加一个小小波浪，加一个小小小波浪。这样表达出来！如下图所示。

而且，不仅仅是方波，任何周期函数，像下图这样的，都可以用一系列（无限个）的正弦和余弦波来叠加表示。

也就是说，

我们在中学数学中已经学习到，一个正弦（余弦）波，可以由频率、幅值、相位三者来决定

所谓频率，就是1秒内，这个波形走了几个周期；所谓幅值，就是波形的高度；所谓相位，就是在0时刻时，波形走到了一个周期的哪个进度。所以，现在，我们知道原信号可以由一系列不同频率的正弦（余弦）波来构成，如下图：

那么，我们以频率为横轴，幅值为纵轴，就画出了原信号的频谱图，如下图所示。同理，也可以画出相位谱。

然而呢，大部分信号都是非周期信号，非周期信号的傅里叶变换如何呢？不要着急，傅大爷也给出了答案。

非周期信号由于没有周期，所以把原信号做傅里叶变换后，在频域上，需要用到所有频率（而不是特定的一些频率点）才能把原信号表示出来。

更进一步，现实中我们能采集到的信号都是有限非周期离散信号，于是又有人基于离散傅里叶变换提出了FFT（快速傅里叶变换）算法，可以把我们采集到的离散信号近似地表示成一系列频率点（不超过采样率的一半）的波形的集合。所以，现在我们明白了，所谓的频域信号，就是把原信号以频率为横轴表示出来。

在对每个通道的信号做了频域变换后，我们就可以考察特定频率或者频率带上的信息了。正是由于我们发现不同频率带上脑电信号的幅值在特定环境下会产生相应变化，我们才定义了Alpha波（8_{13Hz）、Beta波（13}30Hz）等脑电波段，到这里，整个过程是不是很熟悉了？

有关快速傅里叶的讲解：https://blog.csdn.net/enjoy_pascal/article/details/81478582

五、fMRI中的FDR校正

当我们招完被试，收完数据，做完预处理，统计，过五关斩六将，认为自己马上将要发SCI,走上人生巅峰的时候，不好意思，你还需要面对最残酷无情的对手：多重比较校正，比如FDR校正。

（玻璃脑表示，你随意算，有激活区算你赢！）

那什么是多重比较，什么是FDR校正呢？多重比较是统计学中的术语。当我们进行多次统计检验后，假阳性的次数就会增多，所以要对假阳性进行校正。

先解释一些概念：以任务态数据为例，先看一个体素。那么原假设（H0）：该体素没有激活。从原假设基础上进行统计推断，得到P=0.02 < 0.05。所以我们拒绝原假设，然后推断出该体素激活。请注意：我们推断的这个结果有可能是错的，也就是说有可能错误地拒绝了原假设。这种犯错误的概率称之为假阳性率。这个例子里面假阳性率=2%，也就是说该体素激活的这种推断有2%的概率是错的。我们一般显著性水平设置为P=0.05。那么在1000次假设检验中，假阳性的结果就有1000*0.05=50次。我们希望的是假阳性是越少越好，所以要对假阳性进行校正。这个过程就称之为多重比较校正。

那什么是FDR校正呢？

先解释三个指标： Vaa,Via, Da=Vaa+Via。这个三个指标表示的是在V（比如这里V指体素个数）次测试中，它们各自出现了多少次。Vaa表示该体素本来就应该激活（Truly active），我们进行统计推断发现它确实激活了（Declared active）。Via表示该体素本来不应该激活（Truly inactive），我们进行统计推断发现它居然激活了（Declared active）。Da表示我们发现的激活的体素的个数。从上面可以看出，Vaa是我们想要的结果，而Via不是我们想要的结果（这里Via就是假阳性次数），因此我们要控制Via出现的次数。FDR公式如下：

或者，我们可以用中文表达该公式：

FDR校正的是在检测出的激活的体素中（Da），伪激活（假阳性）的体素的个数（Via）。这里要注意和FWE校正的区别。FWE校正的是在所有被检测的体素中，假阳性体素的个数。

最后举个例子：假设大脑总共有50000个体素（V=50000），通过假设检验发现有20000个体素的P<0.05，也就是说Da=20000。FDR和FWE的校正水平都设为0.05。那么FDR说的是在20000个激活的体素中，假阳性的体素不超过20000*0.05=1000个。FWE说的是在总共50000个体素中（包括检测到的激活的体素和不激活的体素），假阳性的体素不超过50000*0.05=2500个。FWE比FDR严格，是因为FWE控制的是全脑中假阳性次数，FDR控制的是我们发现的激活（Declared active）的体素中假阳性次数。

总结起来就三句话：
（1）当同一个数据集有n次（n>=2）假设检验时，要做多重假设检验校正
（2）对于Bonferroni校正，是将p-value的cutoff除以n做校正，这样差异基因筛选的p-value cutoff就更小了，从而使得结果更加严谨
（3）FDR校正是对每个p-value做校正，转换为q-value。q=p*n/rank，其中rank是指p-value从小到大排序后的次序。

举一个具体的实例：

我们测量了M个基因在A,B,C,D,E一共5个时间点的表达量，求其中的差异基因，具体做法：
（1）首先做ANOVA，确定这M个基因中有哪些基因至少出现过差异
（2）5个时间点之间两两比较，一共比较5*4/2=10次，则多重假设检验的n=10
（3）每个基因做完10次假设检验后都有10个p-value，做多重假设检验校正（n=10），得到q-value
（4）根据q-value判断在哪两组之间存在差异

通过T检验等统计学方法对每个蛋白进行P值的计算。T检验是差异蛋白表达检测中常用的统计学方法，通过合并样本间可变的数据，来评价某一个蛋白在两个样本中是否有差异表达。
但是由于通常样本量较少，从而对总体方差的估计不很准确，所以T检验的检验效能会降低，并且如果多次使用T检验会显著增加假阳性的次数。
例如，当某个蛋白的p值小于0.05（5%）时，我们通常认为这个蛋白在两个样本中的表达是有差异的。但是仍旧有5%的概率，这个蛋白并不是差异蛋白。那么我们就错误地否认了原假设（在两个样本中没有差异表达），导致了假阳性的产生（犯错的概率为5%）。
如果检验一次，犯错的概率是5%；检测10000次，犯错的次数就是500次，即额外多出了500次差异的结论（即使实际没有差异）。为了控制假阳性的次数，于是我们需要对p值进行多重检验校正，提高阈值。

方法一.Bonferroni

“最简单严厉的方法”

例如，如果检验1000次，我们就将阈值设定为5%/ 1000 = 0.00005；即使检验1000次，犯错误的概率还是保持在N×1000 = 5%。最终使得预期犯错误的次数不到1次，抹杀了一切假阳性的概率。
该方法虽然简单，但是检验过于严格，导致最后找不到显著表达的蛋白（假阴性）。

方法二.FalseDiscovery Rate

“比较温和的方法校正P值”

FDR（假阳性率）错误控制法是Benjamini于1995年提出的一种方法，基本原理是通过控制FDR值来决定P值的值域。相对Bonferroni来说，FDR用比较温和的方法对p值进行了校正。其试图在假阳性和假阴性间达到平衡，将假/真阳性比例控制到一定范围之内。例如，如果检验1000次，我们设定的阈值为0.05（5%），那么无论我们得到多少个差异蛋白，这些差异蛋白出现假阳性的概率保持在5%之内，这就叫FDR＜5%。
那么我们怎么从p value 来估算FDR呢，人们设计了几种不同的估算模型。其中使用最多的是Benjamini and Hochberg方法，简称BH法。虽然这个估算公式并不够完美，但是也能解决大部分的问题，主要还是简单好用！

FDR的计算方法

除了可以使用excel的BH计算方法外，对于较大的数据，我们推荐使用R命令p.adjust。

1.我们将一系列p值、校正方法（BH）以及所有p值的个数（length§）输入到p.adjust函数中。
2.将一系列的p值按照从大到小排序，然后利用下述公式计算每个p值所对应的FDR值。
公式：p * (n/i)， p是这一次检验的pvalue，n是检验的次数，i是排序后的位置ID（如最大的P值的i值肯定为n，第二大则是n-1，依次至最小为1）。
3.将计算出来的FDR值赋予给排序后的p值，如果某一个p值所对应的FDR值大于前一位p值（排序的前一位）所对应的FDR值，则放弃公式计算出来的FDR值，选用与它前一位相同的值。因此会产生连续相同FDR值的现象；反之则保留计算的FDR值。
4. 将FDR值按照最初始的p值的顺序进行重新排序，返回结果。
最后我们就可以使用校正后的P值进行后续的分析了。

https://blog.csdn.net/zhu_si_tao/article/details/71077703

六、模板（mask）

在这里笔者讲述下mask，与mask使用的技巧与原则。

mask是谁，mask要做什么？

1、模板(mask )往往是与ROI联系在一起的

(region of interest感兴趣区，如果你一定要问我感兴趣区是什么，我觉得我们不在一个频道上，放手吧，我们俩是不可能的)

在小时候填写机答题卡的时候，老师改卷就是在正确答案上面挖洞，然后数个数。如图1示。

图1：挖过洞的正确答案机读卡（其实我觉得全部涂D似乎得分更高，当然了，我绝对没有这么做过，因为会给零分，你至少要涂个C上去）

老师改卷的正确答案就是这张顶层挖了洞的答题卡，就是一个mask。在mask的作用下，你的正确答案（ROI）就一目了然的显现了出来。

总结：在脑影像研究中，mask把我们不关心的区域屏蔽了。（比如你在外面各种浪，想发个朋友圈显摆下你出众，不凡，有品位的生活的时候，是不是要选择不给你老板看呢？这时候你就给你的朋友圈加了个MASK，如果你一定要说你敢不屏蔽你老板，我只能说：大哥在上，请受小弟一拜，你很社会，你是老板）

2、mask作用的原理

用数学的表达式来描述：就是"乘"。

用逻辑学的表达式描述：就是"与"。

图2：mask在计算中处理的机理

这样就把感兴趣的区域信号“520” 提取出来了。

3、常见的mask

根据我们使用的需要，把mask做成：二值模板、多值模板、概率模板。

1）二值模板：

二值模板就是简单的卡一个阈值，白色的区域置为数字1，黑色区域设置为数字0。得到的区域就是我们想要得到的我们所感兴趣的区域。如果需要研究大脑灰质，则把白质、脑脊液与大脑外的区域滤除即可

​ 图3：成年人大脑灰质模板

谁也不能阻挡我作为皇家（ROYAL）高端职业玩家研究脑岛的热情！如果想单独提取脑岛区域作为一个mask呢？想把脑岛作为一个种子点那该怎么办呢？那就要利用下面的多值模板了。

2）多值模板：

多值模板的制作过程可以是解剖形态学分割，也可以是功能相近的区域，相较于二值模板，多值模板利用多个数值来定义不同的脑区，其实简单粗暴的理解就是用一个模板来表示多个二值模板，在这里笔者都摆上来给大伙儿看看。

​ 图4：AAL模板

如图4展现的AAL模板，每一个颜色的区块代表一个不同的区域。不同颜色代表一个不同的数字。然后特定数字对应特定脑区。

​ 图5：yeo的7网络模板

与其说上面的结构模板是按照大脑皮层沟沟回回来分割，功能模板（如图5）就是按照大脑功能区来划分的。我们如果要提取特定脑区的ROI，可以查表（如表1），然后经过**一系列非常高端非常娴熟非常恐怖的操作（这就是我们高端职业玩家与普通玩家的区别之所在，普通玩家看看就行了，千万不要模仿，妄想达到我这个层次）**可以得到的单独一个或几个脑区的二值模板。

表1：多值模板数值与脑区对应表

注：功能模板的制作与表格的来源可以参考这篇文献 The Human Brainnetome Atlas ANew Brain Atlas Based on Connectional Architecture

ref：YeoBT, Krienen FM, Sepulcre J, Sabuncu MR, Lashkari D, Hollinshead M, Roffman JL,Smoller JW, Zollei L., Polimeni JR, Fischl B, Liu H, Buckner RL. Theorganization of the human cerebral cortex estimated by intrinsic functionalconnectivity. J Neurophysiol 106(3):1125-65, 2011.

3）概率模板：

在定义什么是灰质，什么是白质的时候我们应该额外的小心。有时候灰质和白质的边界在MRI上面表现得并不是非灰即白。（要不然为什么freesurfer估计的表面与cat12估计的表面存在一定出入呢）在灰质与白质交界的地方存在一些模糊区域。具体想把大脑灰白质中间这部分模糊区域分割成灰质还是白质需要自己定义。这就诞生了概率模板。

图6：灰质概率模板

图7：白质概率模板

看着图6灰质概率模板，然后看看最初我讲的灰质二值模板。是不是有几分类似？灰质的二值模板就是选取灰质概率大于某一个特定值制作成的。概率模板在spm预处理的空间标准化这一步中起到了关键性的作用。如果您想用儿童的数据就请在网上找找儿童的模板来进行空间标注化。用成年人的模板虽然能够计算出差不多的结果。但是请注意：如果计算儿童数据用成人模板，很有可能被审稿人提问噢。

七、假设检验和效果量

在MRI脑影像领域，统计是几乎必不可少的一环。很多软件如SPM，FSL都可以进行统计分析。我们习惯了点点点。但这背后的机理是什么？

首先要有原假设H0（一般假设某种效应不存在）和备择假设H1（假设效应存在）。然后统计推断得到P值，如果P值很小（比如P<0.05），说明这件事情发生是小概率事件，原假设不成立，从而判定效应存在（是不是与以前学的反证法过程类似）。如果P值过大（如P>0.05），这时候接受原假设。如果P=0.05怎么办？只能说这是一个尴尬的P值，有人称在P在0.05 附近时为边缘显著（当然也有人不认可这种边缘显著）。这里还有个有趣的问题，为什么要以0.05作为显著性的临界值？

下面进入严肃的问题，来了解统计学中一些重要的概念。

举一个实际的例子。环保局要检查某个工厂污染排放是否有问题。假定污染排放量的上限是3。原假设是该工厂污染排放没有问题，环保局派人进行抽样调查，发现该工厂的污染排放量是4，那么我们是否就可以下结论说该工厂有问题。不是！我们还需假设检验，得到P值，如果P>0.05，我们就认为污染排放量4是由于随机抽样误差引起的（刚好抽到了污染多的地方）。如果P<0.05，说明该工厂污染严重。这里注意下。我们说该工厂污染严重有一定几率是错的，即该工厂没有污染，而环保局认为它污染严重（冤枉别人），这种错误称之为I型错误（也叫假阳性率）。还有一种情况，是该工厂污染严重，而环保局认为它没有污染（包庇工厂），这种错误称之为II型错误。具体看下图：

上面图表有两种正确的结果。一是工厂没污染，环保局鉴定过后确实没污染；二是工厂有污染，环保局鉴定后确实污染。这两种结果无需过多关注。我们更感兴趣的是I型错误和II型错误。

那么哪种错误更严重？对于环保局来说，肯定不能冤枉别人，所以应考虑控制I型错误。II型错误的后果是：工厂继续污染，没有得到惩罚，周围百姓继续忍受污染。对于周围百姓来讲，要控制II型错误。那么一个理想的方案是把I和II型错误都控制很小，然而现实是不可能的！！！！！比如要把P控制在P<0.0000000000001，这样我们才拒绝H0（非常小心求证）。那么要找1000条污染证据才能让P达到这样小。但事实上，结果我们只找到20条证据，这时候自己都会对自己说：证据这么少，这个工厂应该没有污染吧！看，II型错误显著上升了。那么有没有办法在其他条件一定的情况下，降低II型错误呢? 唯一的办法就是增加样本量（样本量增多，就有可能找到更多的证据）！！

下面介绍Power。Power=1 - II型错误。II型错误是工厂确实污染，环保局认为没污染。那么Power就是工厂确实污染，环保局认为工厂也污染（正确打击了这种危害性工厂）。所以Power指的是对真实存在的差异正确检测出来的能力。Power越大说明检测差异的能力越大。一种统计方法，即使差异再小，它都能把该差异检测出来，就说该统计方法的Power很大。比如比较两组人的ALFF，如果该统计方法的power=0.8，就是说10个脑区有真实差异，我就能检测出来8个。

下面介绍效果量。

当我们辛辛苦苦收集完数据，统计结果也显著（P值那是相当小），觉得非常perfect的时候，突然审稿人来了一句：请报一下研究的效果量！。你不觉会问：这是什么东东？

效果量，英文名为effectsize。假设对两组数据的均数差异进行统计推断，会得到统计值T值和P值，如果P<0.05，那么就说该差异显著。问题是这样的显著性差异在实际中有没有用？统计推断会受样本影响。比如调查男女身高的差异，在重庆收集了一批样本，发现男性身高显著高于女性。那么这种结论能否推广到其它城市？显然不能。统计推断还会受样本大小的影响。比如研究某治疗方法对治疗抑郁症是否有效，实际结果是实验组比控制组平均高4分，两组人数都是12人，标准差都是8。可以计算P>0.05，不显著。但当两组的人数增加到100（均数差异和标准差不变），差异极其显著。而下结论说该治疗方法有显著效果是不令人信服的。也就是说通过增大样本量达到的统计显著可能并没有实际效果。如果P值很小，但是效果量也很小，就说明即使该治疗方法效果显著，但并不能在实际当中使用。只有那种P值小，效果量也大的治疗方法才能推广使用。

所以效果量反应的是该差异在实际上是否“显著”（不受样本容量大小的影响），而**P值只反应该差异在统计上是否显著。**比如对于男女人数的显著差异（假设男人数>女人数），如果效果量大，表明随便往哪条大街上一站，就能看到男人多于女人。如果效果量很小，那么男人多于女人这种现象可能只限于某局部区域（如某某理工类高校！！！）。正因为效果量重要，所以美国心理学会1994年就发出通知，要求公开发表的研究报告需包含效果量的测定结果。

图2.Cohen’s d图示例

下面介绍几种效果量的计算方法：

八、组水平标准化

在处理脑影像数据时，常有一些刚入门的选手提问：高端玩家啊，一批数据该如何标准化？既然你提问了，我什么也不说也不好，这一期《大话脑成像》，我们专门来讨论这个问题！另外强调：本文所有“标准化”均指单组数据的Z-score变换，并非预处理时图像匹配到标准模板那个步骤！

在MRI领域，标准化常常伴随着数据处理与统计。评价脑影像仅仅采用单一被试来描述往往显得不大合理，所以在分析数据时常常采用组分析方法（batch analysis）。把评价单独被试变换到一组被试上需要引入标准化。在不同的脑影像处理方法中，指标的量纲往往不尽相同，在各个指标之间水平相差巨大的时候如何进行统一分析？避免因为量纲不同带来的假结果？在重测时如何保证数据的可重复性？这就需要对数据进行标准化。

目前标准化的方法非常多，不同的标准化方法带来的评价结果会产生不同的影响，但是在数据标准化方法的选择上并没有标准。这就需要我们了解各种各样标准化方法的机理与可能产生的问题，方便我们有需要的时候进行合理选择。组水平标准化（standardlization）最主要的目的是对一组数据进行比例、放缩变换，把有量纲的数据变成无量纲的数据。**无量纲的数据处理的好处在于不用考虑数据的物理含义，可以在不同单位或者量级的数据之间进行加减或者比较。**所以在这里我们引入今天的重点：Z-score——Z分数化。(敲黑板，划重点)

为了减少各位读者的阅读负担，我们对数学公式进行形象解释（你看高端玩家对你们多好）。**Z分数化：一组数据减去均值除以标准差。**经过Z-score标准化后，组内的数据服从标准正态分布（如图1：均值为0，组标准差为1）。经过Z-score标准化的洗礼，所有的数据变成了没了单位的纯数量值。数据做过Z-score标准化后，把标准差变成一个单位的距离，非0的Z-score值表示距离平均（也就是0）的距离。查看图1，基本99.9%的数据都落在3倍标准差以内。

图1：标准正态分布图（from Wikipedia）

注：

Z-score的数学定义：

其中： Z 为 Z分数，x 为一组数据值，μ为均值（或样本均值），δ为标准差（或样本标准差）

在使用Z-score标准化时，有两个问题需要注意：

（1）此"组"非彼"组"

在做Z-score标准化时，心里一定要清楚我们所选择的“组”是什么（定义）。这里的**“组”就是指一系列值的意思**。举个例子：如果对一组人的某个ROI进行Z-score标准化时，“组”就是这组被试·一个ROI信号值序列。如果对一个人的3维或4维脑图进行Z标准化时，“组”就是该被试·大脑区域·某个时间点的所有数据值。注意这里所指大脑区域，表示非大脑区域的信号一定不能混入其中，需要选取这个被试的大脑模板。

有读者可能会产生疑问：为什么非大脑区域对Z-score标准化影响巨大呢？我们把问题想得极端一些就明白了。如下图4×6的矩阵代表整个大脑区域，其中红色框部分为大脑区域，两个非零元素代表大脑区域数值。根据定义：减去均值除以标准差，加上非大脑区域，均值产生了明显降低，故导致结果错误。

图2：极端化的大脑区域与外界区域

（2）Z-score 与 Fisher-Z 的区别

很多人在Z-score 与 Fisher’s z 变换上面傻傻分不清楚，包括某些文章上面所得到的结果就是错误的。因为某些文章把Z变换用成了fisher z。Z-score 与Fisher-Z其实没有任何联系，没有任何联系，没有任何联系。只不过名称相似，又常常同时在磁共振数据处理中使用，难免混淆。Z-score，又称Z分数化，“大Z变换”，Fisher-z，又称Fisher z-transformation，“小z变换”。

Fisher’s z 变换，主要用于皮尔逊相关系数的非线性修正上面。因为普通皮尔逊相关系数在0-1上并不服从正态分布，相关系数的绝对值越趋近1时，概率变得非常非常小。相关系数的分布非常像断了两头的正态分布。所以需要通过Fisherz-transformation对皮尔逊相关系数进行修正，使得满足正态分布。

图3：fisher‘s ztransformation（from Wikipedia）

注：

相关系数定义：

fisher‘s z transformation：

其中：r 为相关系数，ln为自然对数。

最后，我们用一句话进行总结：Z-score 标准化，用于一组数据去量纲，变换后得到数据均值为0，标准差为1.

九、由 ALFF 说开去

我有两个中心采集的数据，只不过一个中心采集的是正常人的数据，另外一个中心采集的是病人的数据。我能不能算ALFF，并直接拿两个不同机器的fMRI的ALFF值相减呢？

要回答这个问题，我们可以把上面这个问题拆分成以下几个小问题，如果把小问题（小目标）解决了，大问题也就迎刃而解了：

1、什么是频率域？

2、什么是ALFF？

3、为什么要（为什么不要）算ALFF？

4、ALFF能不能直接相减？

小目标1：什么是频率域？

频率域，简单来说：只是另外一种描述信号的手段，与时域类似，只不过我们更容易理解时域罢了。

更简单来说：就是一个公式

图1：信号，时域，频域公式

频域方面的计算永远不是脑电的代名词，相反，我们往往可以师夷长技以制夷。更换一种看问题的角度往往能够更好探究人脑的内部机制。例如把脑电的相关知识运用在磁共振脑成像领域中。当然在运用的时候需要严格把关底层的原理，否则就容易做成了四不像。本文最后提供部分频域方面的参考文献。

理解上面这个图，频域方面就可以PASS了。

小目标2：什么是ALFF？

ALFF的全称： Amplitude of Low Frequency Fluctuations（低频振幅）

笔者找到推广ALFF的第一篇文章，截了一个图，并稍微做了些修改，为了方便理解把它分成了好几个大块。（记住这个图，后面的大话脑成像系列笔者还会再拿这张图说其他事 ^_）

图2：ALFF计算过程

ALFF具体是怎样来的，死磕文章的定义：Amplitude of Low Frequency Fluctuations（低频振幅）

a）首先需要看“ Low Frequency （低频）”。实现磁共振信号低频保留当然使用低通滤波器(注1)。立即送入低通滤波器low-pass filter 实现A->B的步骤。

b）时域上面看腻了以后，想看看频域上有没有什么出彩的地方。所以需要把时域上的信号转换到频域上看看。简单来说就是换一个视角，把时域信号转换到频谱上面看看(注2)。在这里需要注意一点。B图和C图的横坐标和纵坐标分别代表什么？ B图和C图的纵坐标都代表着信号的强度。而横坐标不相同，一个代表着时间，另外一个代表着频率。这样就实现了B->C。

c）对C中每一个点的纵坐标数值进行开平方取平方根，得到D。

d）对这根线取平均（各个离散的点计算出均值）得到E。其中 E 就是计算出来的ALFF值（求出振幅）

e）剧本不是这样写的，怎么能把“F”丢掉了呢！其实呢，不要在意这些细节，F图是在E的基础上除以全脑每一个体素ALFF的均值，进行了除均值的标准化。

一句话，死磕定义，就明白了ALFF是什么了。

注1：a）部分为了方便理解，以低通滤波器来描述，这样描述得不够严谨。ALFF内部计算是使用带通滤波器，因为在ultra-slow频段上的生理意义存在争议。部分文章认为极低频信号属于低频漂移，随磁共振机器温升而变化。而反对者认为极低频信号存在生理意义，它可能反应生物内部自发振荡节律（生物钟）。

注2：b）部分高端玩家会问，按照图1的总结，时域转换为频率域时除了幅频谱以外还应有相位谱才对，而在图2中却没有相位谱。绝大数文章做FFT的时候都只考虑了其中的幅频谱，而丢掉了相位谱。换句话说，计算ALFF将时域转换为频域时，丢失了部分信息。部分文献为了捡回这部分信息，采用了动态ALFF的方法。但是笔者认为这种方法没有直接用相位谱来得好。

小目标3：为什么要（为什么不要）计算ALFF？

a）为什么要：

为什么要用ALFF来刻画大脑呢？还是死磕定义：**低频振幅。振幅可以代表着大脑活动强度。**而BOLD机制告诉我们，**大脑激活和不激活可以通过BOLD信号强度来反映。**根据小目标2的计算过程，大脑中的每一个体素都能够计算出一个ALFF值，而每一个ALFF的值代表着大脑中该体素的活动强度。

图3：BOLD机制

比较ALFF的差异无非是比较两个脑区之间信号的活动强度，换做在癫痫发作时，可以是某个特定脑区的异常激活。

b）为什么不要：

既然想看看低频振幅。为什么不直接滤波后取平均，得到一个振幅数值作为笔者定义的ALFF呢？为什么要转换到频域再平均呢？

小目标4：ALFF能不能直接相减？

问ALFF能不能直接相减前，想想看，BOLD信号能不能直接相减。理解不了？看看下图就了然了。

图4：不同机器（不同扫描参数）下BOLD信号差异

什么？曲线是一样的？不，你看看纵坐标。不同的扫描参数会影响磁共振的对比度和强度。不同参数的BOLD信号不能比较的，这样比较是没有意义的。

你可能会问：那么我算ALFF岂不是不能比较了？

我的回答：仔细看看小目标2中，为什么作者需要多加一个E -> F 的步骤，作者并不会平白无故再多加一个除均值标准化的！他的用处体现在这里啦。经过这个步骤，组间就可以进行比较了

笔者有一个问题：既然想看看低频振幅。为什么不直接滤波后取平均，得到一个振幅数值作为笔者定义的ALFF呢？为什么要转换到频域再平均呢？

频率，频域相关参考：

Wang Y F, Liu F, Long Z L,et al. Steady-state BOLD response modulates low frequency neural oscillations[J].Scientific Reports, 2014, 4: 7376.

Wang Y, Zhu L, Zou Q, etal. Frequency dependent hub role of the dorsal and ventral right anteriorinsula[J]. NeuroImage, 2018, 165: 112-117.

ALFF的相关参考：

Yu-Feng Z, Yong H,Chao-Zhe Z, et al. Altered baseline brain activity in children with ADHDrevealed by resting-state functional MRI[J]. Brain and Development, 2007,29(2): 83-91.

Zou Q H, Zhu C Z, Yang Y,et al. An improved approach to detection of amplitude of low-frequencyfluctuation (ALFF) for resting-state fMRI: fractional ALFF[J]. Journal ofneuroscience methods, 2008, 172(1): 137-141.

十、计算机存取MRI影像的那些事

昨天，小芳（隔壁村的）问笔者：为什么我输出不了超过256个大脑区域？

宏观来讲，普通玩家对脑影像分析处理的步骤无非： 读取 -> 分析处理 -> 输出（写入保存）。这三个步骤会变的也就是中间这个步骤：分析处理。结合自己的问题来找到自己特定的分析处理方法。但是今天的重点主要是前、后两个步骤——数据读取与结果保存。通过MATLAB底层函数读取一个功能磁共振影像或者结构磁共振影像。

这些功能在基础的软件包里面已经集成了，但是对于喜欢脑 （no ）洞 （zuo ）大 （no ）开（ die ）的我们来说，普通的软件包点点点已经不够了，这时候就需要一些深入底层的操作。从内部剖析磁共振数据了。

从磁共振机器里面拿出来的数据分成两个部分：头+脑影像。如下图：

图1：磁共振脑影像数据结构图

可以用下图来理解：一个火车头，火车头里面装载着这个图像的信息，这些信息包含着层厚，层数，体素大小等等描述后面数据的各个信息。而每一层的脑影像就存放在后面的多个车厢里面。换句话说，那个火车头就像目录，而后面的车厢存放的是每一层大脑切片的信息。如图2所示：在计算机里面，每一层（大脑切片）的信息就是一张二维的图片，而数字的图片在计算机里面用离散的数值进行储存。每一个数值代表着当前切片颜色的深浅。

图2：每一个车厢里面装载的信息

上面的前期铺垫工作已经做好，读到这里，想必各位读者了解了磁共振脑影像在计算机中存储的数据结构形式，接下来进行非常高端非常娴熟的操作。

首先我们在MATLAB下面安装好SPM与 rest软件，基本的rest软件操作想必各位都会。如果你想从中挖出一些更加炫酷狂拽的东西，可以操作rest软件包里面底层的函数对脑影像进行分析。

1、读取：

比如对于今天我们先读取一个脑影像：运用rest_ReadNiftiImage函数。该函数基本用法如下：

​ [Data, Head] = rest_ReadNiftiImage(fileName,volumeIndex)或

​ [Data, Head] = rest_ReadNiftiImage(fileName)

别慌，我来讲解一下。

你可以理解为：

​ [火车车身，火车头 ] <= 函数名 （文件存放路径）

其中：

​ 文件存放路径：你的C盘下存放一个AAL模板（名字叫aal.nii），如果你用windows系统那么文件路径就是 C:\aal.nii ，如果你是用linux或者mac系统的话，那么就是 /home/[your name]/aal.nii

​ 火车车身：你读取出来的数据（aal模板是3维的数据，你可能看到他的维度是91×109×91的）

​ 火车头：读取出来数据的头文件。（不用管它是什么啦）

2、写入（保存）

我们对这个火车车身的数据进行了一些非常恐怖的操作后（比如高阶相关、换成统计参数值、计算结构协变网络等等，恐怖吗？害怕吗？）要对结果进行保存，该怎么办呢？其实这里很简单，如果是全脑体素水平的操作，仍然把那个冒着黑烟的火车头拿过来，这次拉的货物就不是新鲜的猪肉了，而是已经炖好的猪肉了。换句话说：旧瓶装新酒。

保存方式如下：

rest_WriteNiftiImage(Data,Head,filename)

看见了没有，data已经换成了新的，而Head还是咱们当时读入的火车头，filename 就是我们需要保存的位置与文件名。这个文件名可以自己取，运用这个函数会自动生成。

总结：对于写入和保存，读取文件的时候拉来了一火车的数据，这一火车的数据包含着：火车头，和很多车厢的数据。我们对各个车厢的数据进行处理，处理完毕后继续装回车厢，挂上火车头继续开走。

在解答以上问题之前，我问各位一个问题：我们的磁共振脑影像数据大小为什么有的一样，有的不一样？换句话问：数据大小由什么决定的？

图3：不同影像数据，不同大小。

看过网络电视的朋友们肯定知道一点：高清4K和流畅画质是两个截然不同的概念。一般看着超清的视频消耗的流量更多，而看普通标清的画质的视频消耗的流量更少。于此进行类比，画面的分辨率决定着数据的大小。在磁共振这里也是类似——为什么一个3维脑图，T1加权的结构像更加的清晰，而功能像更加的模糊？这个是因为T1像分辨率更加高。

结论1： 分辨率更高的影像数据更大

但是我们的问题仍然没有解决，我明明存放的是有小数的，为什么全部变成了整数？

在这里补充一点计算机的背景知识：

​ 在计算机中用 char 型的数据用1个字节保存，一个字节8位二进制。表示范围2^8=256。换句话说，char型数据只能保存 -128~127 之间的数值，而unsigned char就抹去了负数部分，描述了 0~255范围之间的数值。（绿皮有窗火车）

​ 在计算机中 short 型数据 用2个字节保存，也就是16位二进制 2^16 = 65,536个数值，但是注意，他与char型一样，都只能够描述整数。（普通红皮蓝皮车）

​

​ 在计算机中 double 型数据 用8个字节保存，自己算算可以表示多少数据。double可以描述小数，但是注意：这个小数是有16位精度限制的。也就是说double是用科学计数法来表