Hello，大家好！
这里是行上行下，我是喵君姐姐~
脑电图(EEG)是测量人脑神经元活动最原始的技术之一，特别是在临床诊断中具有无可争议的价值。它能够识别癫痫和睡眠障碍，以及对感觉传导通路功能障碍进行评估。然而，EEG的缺点主要在于其空间分辨率较差，这使得在推断头皮上产生的神经元活动的大脑区域位置时较为困难。
因此，研究者们旨在找到一个更精确、更可靠的解决方案以定位EEG活动产生位置，并能够结合精确的头部解剖信息和复杂的源定位算法，将其转换为一个真正的神经成像模式的高密度脑电图系统。
本期，我们就来介绍将采集到的EEG信号转换为三维神经活动成像时的一系列步骤。图1列举了EEG源定位中常用的一些软件和安装包。这里我们使用的是一款独立免费的软件：Cartool。
PS：以下一些翻译可能还有待商榷，欢迎大家留言批评指正！
原始EEG数据受到许多非生理源(电源线、不良电极接触等)和生理源(心跳、肌肉活动、出汗、运动等)的伪迹污染。由于数据集较大，我们通常会采取自动对伪迹进行识别和消除的方式。但由于造成伪迹的因素较多且不确定，因此在自动识别之后还需要人工检查。下面我们就来主要地介绍EEG预处理的一系列步骤。
首先，我们需要使用时间滤波器（Temporal Filtering），以去除数据中那些被认为是非生理源或与当前研究无关的频率。在Cartool中，我们采用的是2阶无限脉冲响应（IIR）巴特沃思滤波器。
它的优点在于具有最佳的平通带响应，能够限制局部最大值位置上出现伪迹。另外，巴特沃斯低通和高通滤波器都能够消除相位移动，这确保了局部最大值保持在合理的范围，而与它们的频率无关。
滤波后，我们便可以对数据进行下采样，即降低采样率。它可以通过减少样本量来降低后续处理所需的内存，而不会丢失任何信息。奈奎斯特定理要求向下采样不低于最高剩余频率（highest remaining frequency）的两倍。
然而，在实际中，由于滤波器的截止频率不可能完全清晰，为了保持一些额外的时间分辨率，最终的采样频率应该选择在最高剩余频率的4倍左右。在Cartool中，下采样主要是级联积分梳状(CIC)滤波器来完成的，它在离线状态下也能用于计算。
独立成分分析（ICA）一般用于检测和校正伪迹，特别是眼球运动、眨眼和心跳伪迹等。它能够计算出反映特定事件的独立成分以及伪迹，从而对伪迹进行排除。
在Cartool中，数据检测是半自动的。我们可以检测出明显偏离平均振幅的电极。这些电极可能是由于接触不良而导致的。因此，在计算ICA时要对这些电极进行标记并忽略。而在ICA校正后，检测到的不良电极需要采用三维或球面样条算法（spherical spline algorithm）来进行插值。
由于Cartool计算ICA还并不全面，最好的方法是使用EEGLAB来计算ICA。
在获得EEG源成像之前，我们需要继续对头皮电位图进行观察和处理。即使在人工电极插补和移除无关的ICA成分后，瞬时事件仍可能在短时间内损坏少数电极，并对源定位产生影响。
在Cartool中，我们介绍了一个空间滤波器，它是一种瞬时滤波器，通过空间平滑地图来去除局部异常值，而不丢失其地形特征。
方法如下:
所以我们有一个间隔加权平均值。对于每个电极e:
在经过上面几个步骤的处理之后，可能还仍然存在残留的瞬时伪迹，如肌肉运动、出汗、眨眼等。因此我们需要对某个“不良时段”进行人工标记，并删除。
在Cartool中，我们可以对轨迹进行简单统计，估计每条轨迹与它自己的基线偏离的程度，从而识别这些“不良时段”。它主要是计算瞬时值(在给定时间点下电极的绝对值、方差、偏度和峰度)和时间的交叉卷积，这是一种估计信号中的噪声比较方便的方法。
所有这些偏离估计量被合并到一个单一的复合估计量，明显偏离平均值的时间段被突出显示。通过视觉检查，我们可以决定是否应该将某段标记为“不良时段”。尽管Cartool的许多工具允许跳过它们，但将它们进行移除会便于后期处理。
头部模型是计算脑电正解（forward solution）的模型。正解决定了大脑中的某个特定源对头皮上每个电极的影响程度。我们可以使用参与者的单个MRI来构建头部模型。
由于从源到头皮表面的电场会被这些隔室(尤其是颅骨)减弱，因此在头部模型中适当地纳入头部形状和电导率参数对于EEG源重建至关重要。我们需要对MRI进行预处理，以得到不同脑室(如皮肤，颅骨，脑脊液，大脑)之间的电导率参数。
在预处理时，电极必须和记录过程中电极放置的位置一致，以保证源定位的正确性。
EEG源成像的头部模型是基于MRI的，为了正确提取与头部形状、颅骨厚度和灰质体积等有关的信息，需要进行以下几个步骤：重新采样和重新定位、颅骨剥离、偏场校正，以及灰质和白质分离。下面介绍如何使用Cartool实现这些步骤。
第一步需要对MRI重新采样并重新定位。如果采集是各向异性的，通常情况下，体素的大小在所有三维中都是不相等的（图4A）。这对于任何的三维处理，都是非常有害的。
在Cartool中，可以通过一些线性调整来简单地对最低分辨率的轴进行上采样，以达到3个轴的最高分辨率。Cartool默认将MRI转换为X、Y、Z三个轴的右-前-上(“RAS”)方向(右手系统)，类似于MNI模板大脑。这是通过适当旋转90度来实现的（图4B），由此使得三维中的体素大小一致。
图5展示了一个最佳切割面，并确定了原点
在分析过程中，我们需要得到一个标准化的个体头部模型。因此，我们的下一步是剥离颅骨以分离颅骨、脑脊液和大脑（图6）。
MRI通常在空间中存在不均匀性，这被称为偏差场。在这种情况下，脑组织（如灰质）会在扫描仪中的不同物理位置上存在不同的值（图7左），这会阻碍脑组织的分离。因此在Cartool中，我们需要对其进行校正（图7右）。
MRI处理的最后一步便是白质与灰质分离。这是因为EEG源定位通常将源空间限制在包含突触的灰质中。Cartool能够使用混合高斯函数来估计灰质、白质和脑脊液的整体强度分布，从而提取gray mask（图8）。
通过灰质分离而获得体积称为solution space，上面的电极活动能够被定位。solution space通常包含3000-6000个单独的solution point。
solution point提取是源定位中十分重要的一步。我们需要注意以下几个方面：
1）左右分布
不对称的solution point分布会对源定位产生影响，因为它会给大脑的一侧赋予更多的权重，并将源定位到错误的一侧。因此我们需要保证它们在大脑左右半球是平均分布的。
2）临近最小值
由于需要在solution space中计算一个离散的拉普拉斯矩阵，因此在Cartool中，每个solution point的附近都应该有足够的相邻solution point。
3）连续性
由于源活动无法在不存在的solution point上重建，导致部分大脑区域缺乏精度。因此我们需要保证solution point的连续性（图9）。
4）solution point的数量
通常最佳solution point的数量在3000-6000，数量太少或太多都存在不同的利弊（图10）。
一些研究表明，电极数量较少时，不仅会导致分辨率模糊，而且还会产生不正确的定位。一定范围内的高密度的电极分布更能帮助我们精准定位。
除了电极的数量，电极的位置在头部的覆盖也起着重要的作用。标准的10-20系统不包括头部下方的电极，这不利于正确记录颞叶基部下方和前部的活动，而起源于或传播自颞中结构的活动是最大的（图11）。因此，建议临床常规标准10-20系统中每侧至少增加3个下电极。
对于MRI，在被试头部表面正确安放电极非常重要。最精确的方法就是让被试戴上脑电帽，进入扫描仪中，然后把电极引起的伪迹在MRI图像上标记出来（图12）。
在Cartool中，我们需要通过显示和手动调整电极阵列的整体三维形状来实现三维电极阵列的协同配准。一旦调整完成，Cartool还提供了最后一个有用的功能:将电极“粘合”在头部，这是一种虚拟功能。通过激活这种虚拟的粘合，所有电极都将完美地垂直投射到头皮上最近的位置（图13）。
为了计算导程电场（lead field），必须继续创建一个头部模型，该模型必须尽可能真实地结合头部的形状和不同组织的电导率参数，包括大脑中的电流源和头皮上的电位。最常用的模型是BEM和FEM。
在Cartool中，还有一种LSMAC模型。尽管三者的复杂度不一致，但它们都能够提供相似的源定位。下面我们主要地介绍一下LSMAC模型。
对于LSMAC头部模型，在每个电极下，颅骨的内外边界能够被自动地确定，并对整体电阻率值进行局部校正，减少了脑电信号逆建模的误差来源。头骨的边界是通过分析从大脑中心到头皮电极的辐射线的灰色水平来确定的。图14为在3个电极上估算颅骨半径的示例。
由于这种测量方法还存在一定的不确定性，Cartool要求提供被试年龄，并绘制除了不同年龄下的头骨厚度变化（图15A）。
另外，颅骨电阻率也会随着年龄的变化而变化。因此在Cartool中，给出了颅骨相对于邻近组织(大脑、脑脊液、头皮)的相对导电性与年龄之间的函数（图15B）。
在Cartool中，计算了三种线性源分布模型：LORETA、LAURA、NM。我们在一些实验和临床研究中验证了这些方法，并将其与颅内记录、电皮层刺激、fMRI和术后临床结果进行了比较。分析如下：
正则化方法（Tikhonov regularization）通常被用于确定方程组的情况下，如在转换lead field时。正则化程度越高，结果越平滑，对噪声的敏感性越低。然而，过多的正则化，过度平滑的结果，将降低定位的精度。
为了处理所有情况，Cartool计算了13个矩阵，这些矩阵从0(无)增加到12(对于非常嘈杂的数据)，并乘以一个常数∝（取决于所选择的逆模型），可以得到各个矩阵的正则化因子。
我们可以根据实际的EEG信号选择一个最佳矩阵。∝R代表正则化因子，R代表矩阵，K代表Lead Field，W代表一些特定的逆权重因子（图16）
对于给定的EEG数据，最优正则化的定义是：位于solution point与正则化因子R的函数拐点处。
在EEG进程中，我们可以观察到每一个solution point的功率会发生显著变化。为了可靠地估计被试大脑活动随时间的波动，并在不同被试之间进行比较，我们需要找到一种方法来校正这种功率变化。
Cartool提供了一种标准化方法来计算校正因子。在标准化过程之后，电流密度的功率在所有solution point上都可以进行比较，且噪声分量是正态分布的。
首先，图17A展示了标准化之前的3个solution point信号，它们的平均幅值存在差异；图17B为标准化之后的信号，可以看到它们的平均幅值已经在相同的范围内。
同样，标准化之前的背景活动主要集中在左侧，且不一致（图18A）；标准化之后的背景活动都集中在1上（图18B）。
图19为所有solution point的源幅值分布，在标准化之后，它们的最大幅值范围（红色区域）更为集中。
逆矩阵与脑电信号相乘的输出结果在于每个solution point上都有等效的偶极子。由于每个偶极子都是一个三维矢量，它也可以代表x、y和z方向上的振幅，即标量。大多数时候，我们所导出的都是振幅，而不是偶极子。
图18展示了分布源的实际矢量结果(左上)和相应的振幅值(右上，底部为横向切片)（图20）。
EEG源定位越来越多地用于推断在头皮上观察到的活动的大脑区域。在临床应用方面，EEG源定位广泛地应用于癫痫，目的是定位癫痫区以进行药物治疗。
除了临床意义外，癫痫EEG源定位也为评估不同头部源模型的性能和精度提供了独特的可能性，因为颅内记录或术后结果可作为“黄金标准”。
在实验研究中，EEG源成像已经成为涉及感觉、运动和认知功能的不同脑区域定位的标准，最常用于事件相关电位。EEG源成像通过连接性测量，也成为了定义大规模网络动态的一种方式，它是功能连接分析的先决条件。
本文从细节上描述了从多通道头皮EEG记录到评估潜在神经源分布的一系列步骤，解释了每个步骤背后的逻辑，执行这些步骤需要满足的条件，以及如何在Cartool中实现。这对我们具有极大的启发意义。
PS：在本文首发的微信公众号 行上行下 的后台回复关键词“脑电溯源”即可获得所述原文及参考文献啦！
参考文献
