说明
转到末尾下载完整示例代码。
SourceEstimate 数据结构#
源估计通常称为 STC(Source Time Courses,源时间过程),由源定位方法获得。源定位方法用于求解所谓的“逆问题”。MNE 提供了多种求解方法:dSPM、sLORETA、LCMV、MxNE 等。
源定位就是把 EEG/MEG 传感器数据投影到位于个体受试者大脑解剖结构中的三维“源空间”中。因此,数据经过变换,使每个传感器位置上记录的时间序列映射到“源空间”中每个空间位置上的时间序列,而源估计正是在该源空间中定义的。
STC 对象包含源随时间变化的幅值。它只存储激活的幅值,而不存储源的位置。要访问位置信息,你需要使用计算正向算子(常缩写为 fwd)所用的源空间(常缩写为 src)。
关于正向建模的更多细节,请参阅 头模型与正向计算;关于使用 dSPM、sLORETA 或 eLORETA 进行源定位的示例,请参阅 使用 MNE、dSPM、sLORETA 和 eLORETA 进行源定位。
源估计有多种不同的形式:
mne.SourceEstimate:用于皮层约束的源空间。
mne.VolSourceEstimate:用于体积源空间。
mne.VectorSourceEstimate:用于皮层约束且源激活为矢量(强度与方向)的源空间。
mne.MixedSourceEstimate:用于由皮层约束源与体积源组合而成的源空间。
说明
(矢量)SourceEstimate 属于表面表示,通常与 FreeSurfer 表面表示配合使用。
我们先来了解 mne.SourceEstimate 究竟是什么。首先搭建环境并加载一些数据:
# Authors: The MNE-Python contributors.
# License: BSD-3-Clause
# Copyright the MNE-Python contributors.
from mne import read_source_estimate
from mne.datasets import sample
print(__doc__)
# Paths to example data
sample_dir_raw = sample.data_path()
sample_dir = sample_dir_raw / "MEG" / "sample"
subjects_dir = sample_dir_raw / "subjects"
fname_stc = sample_dir / "sample_audvis-meg"
加载并检查示例数据#
该数据集包含一项视听任务中的源估计数据。数据已使用 dSPM 方法映射到由 FreeSurfer 获得的膨胀皮层表面表示上,并凸显出听觉皮层中一个明显的峰值。
我们来看看它是什么样子。
stc = read_source_estimate(fname_stc, subject="sample")
# Define plotting parameters
surfer_kwargs = dict(
hemi="lh",
subjects_dir=subjects_dir,
clim=dict(kind="value", lims=[8, 12, 15]),
views="lateral",
initial_time=0.09,
time_unit="s",
size=(800, 800),
smoothing_steps=5,
)
# Plot surface
brain = stc.plot(**surfer_kwargs)
# Add title
brain.add_text(0.1, 0.9, "SourceEstimate", "title", font_size=16)

SourceEstimate(stc)#
源估计包含源空间所定义的空间位置上激活的时间序列。在 FreeSurfer 表面(由三维三角剖分构成)的语境下,我们可以把膨胀大脑表示上的每个数据点称为一个顶点(vertex)。如果每个顶点代表一条时间序列的空间位置,那么时间序列和空间位置就可以写成一个矩阵:每个顶点(行)在多个时间点(列)上都被赋予一个值,这个值就是信号在给定空间点和时间点上的强度。存储在 stc.data 中的正是这个矩阵。
我们来看一下它的形状
shape = stc.data.shape
print(f"The data has {shape} vertex locations with {shape} sample points each.")
The data has (7498, 25) vertex locations with (7498, 25) sample points each.
可以看到,stc 携带了 7498 条长度为 25 个采样的时间序列。这些时间序列属于 7498 个顶点,而这些顶点又代表皮层表面上的位置。那么这些顶点值来自哪里呢?
FreeSurfer 把左右半球分开,并分别为左半球和右半球创建表面表示。表面位置的索引存储在 stc.vertices 中。这是一个包含两个整数数组的列表,每个数组索引 FreeSurfer 网格中的某个特定顶点。因此,数值 42 会映射到网格中索引为 42 的顶点的 x、y、z 坐标。关于如何访问 mne.SourceSpaces 对象中的位置信息,请参阅下一节。
由于两个半球总是分别表示,上面介绍的两个属性也可以通过选择相应半球来获得。方法是加上正确的前缀(lh 或 rh)。
shape_lh = stc.lh_data.shape
print(
f"The left hemisphere has {shape_lh} vertex locations with {shape_lh} sample points"
" each."
)
The left hemisphere has (3732, 25) vertex locations with (3732, 25) sample points each.
由于我们没有改变时间表示,改变的只是所选顶点的子集,因此只有矩阵的行数发生了变化。我们可以检查 stc.lh_data 与 stc.rh_data 的行数之和是否等于之前的数值。
is_equal = stc.lh_data.shape[0] + stc.rh_data.shape[0] == stc.data.shape[0]
print(
"The number of vertices in stc.lh_data and stc.rh_data do "
+ ("not " if not is_equal else "")
+ "sum up to the number of rows in stc.data"
)
The number of vertices in stc.lh_data and stc.rh_data do sum up to the number of rows in stc.data
确实如此,而且正如细心的读者已经猜到的那样,顶点也是如此。stc.lh_vertno 对应 FreeSurfer 的左侧膨胀表面表示,stc.rh_vertno 则对应右侧膨胀表面表示。
与 SourceSpaces(src)的关系#
如上所述,src 携带了从 stc 到表面的映射。每个半球的表面都由三角网格构成,构成一个面的每个三角形由 3 个顶点组成。由于 src 是两个源空间(左半球和右半球)的列表,我们可以先选择源空间,再访问相应的数据。构成左半球的面可通过 src[0]['tris'] 访问,其中索引 \(0\) 代表左半球,\(1\) 代表右半球。
src[0]['tris'] 中的值指向 src[0]['rr'] 中的行索引,这里存放的是表面网格的实际坐标。因此,每个顶点索引值都会从这里选取一个坐标。此外,src[0]['vertno'] 存储的数据与 stc.lh_vertno 相同,除非使用 mne.inverse_sparse.mixed_norm() 这类稀疏求解器,因为此时只有一小部分顶点真正具有非零激活。
换句话说,stc.lh_vertno 等于 src[0]['vertno'],而 stc.rh_vertno 等于 src[1]['vertno']。因此,stc.lh_data 中的第 N 条时间序列分别对应 stc.lh_vertno 和 src[0]['vertno'] 中的第 N 个值,这些值又反过来把该时间序列映射到表面上的特定位置,该位置由 src[0]['rr'] 中的笛卡尔坐标集合 stc.lh_vertno[N] 表示。
我们获取数据的峰值幅值,以及对应的顶点和时间点索引:
peak_vertex, peak_time = stc.get_peak(hemi="lh", vert_as_index=True, time_as_index=True)
其中第一个值表示使用的是 stc.lh_vertno 或 stc.lh_data 中的哪个顶点,第二个值表示哪个时间点索引。我们可以利用这些信息得到对应峰值的表面顶点索引及其数值。
我们也可以使用与开头相同的 surfer_kwargs 对此进行可视化。
brain = stc.plot(**surfer_kwargs)
# We add the new peak coordinate (as vertex index) as an annotation dot
brain.add_foci(peak_vertex_surf, coords_as_verts=True, hemi="lh", color="blue")
# We add a title as well, stating the amplitude at this time and location
brain.add_text(0.1, 0.9, "Peak coordinate", "title", font_size=14)

总结#
stc 是 MNE-Python 中的一个类,表示从源估计获得的变换后时间序列。对于两个半球,数据分别以 \(m \times n\) 矩阵的形式存储在 stc.lh_data 和 stc.rh_data 中,其中 \(m\) 是该半球对应的空间位置数,\(n\) 是时间点数。
stc.lh_vertno 和 stc.rh_vertno 分别对应 src[0]['vertno'] 和 src[1]['vertno']。这些是表面表示上位置的索引。
表面的网格坐标分别存储在左半球的 src[0]['rr'] 和右半球的 src[1]['rr'] 中。可按上述逻辑访问三维坐标:
>>> lh_coordinates = src[0]['rr'][stc.lh_vertno]
>>> lh_data = stc.lh_data
或者:
>>> rh_coordinates = src[1]['rr'][src[1]['vertno']]
>>> rh_data = stc.rh_data
脚本总运行时间:(0 分钟 3.834 秒)