# Python实战：用MNE和NetworkX构建脑电复杂网络（附完整代码） 脑科学研究正在经历一场从局部到整体的范式转变。过去我们关注单个脑区的激活模式，现在则更重视不同脑区之间的动态交互。这种转变催生了一个新兴领域——脑网络分析，它借助复杂网络理论来量化大脑各区域之间的连接模式。对于Python开发者而言，MNE和NetworkX这两个强大的工具库为我们提供了从原始脑电数据到复杂网络分析的全套解决方案。 本文将带你完整走通这个技术流程：从EEG数据预处理开始，到计算连接矩阵，再到构建复杂网络并进行图论分析。我们不仅会讲解核心概念，更重要的是提供可直接运行的代码示例，让你能够快速将理论转化为实践。无论你是神经科学研究者还是对脑机接口感兴趣的开发者，这些技术都将为你打开一扇新的大门。 ## 1. 数据准备与预处理 任何脑电分析的第一步都是确保数据质量。原始EEG信号中混杂着各种噪声，包括眼动、肌电、电源干扰等。MNE-Python提供了完整的预处理流水线，让我们能够高效地清洗数据。 首先安装必要的库： ```bash pip install mne numpy scipy networkx matplotlib ``` 加载示例EEG数据并进行基本预处理： ```python import mne from mne.datasets import sample # 加载示例数据 data_path = sample.data_path() raw_fname = data_path / 'MEG/sample/sample_audvis_filt-0-40_raw.fif' raw = mne.io.read_raw_fif(raw_fname, preload=True) # 选取EEG通道 raw.pick_types(eeg=True, meg=False) # 基本预处理流程 raw.filter(1, 30) # 带通滤波1-30Hz raw.notch_filter(50) # 陷波滤波去除电源干扰 # 分段为epochs events = mne.find_events(raw) epochs = mne.Epochs(raw, events, tmin=-0.2, tmax=0.5, baseline=(-0.2, 0)) ``` 预处理的关键步骤包括： - **滤波处理**：保留与认知活动相关的频段（通常1-30Hz） - **伪迹去除**：使用ICA或回归方法消除眼动和肌电干扰 - **分段对齐**：根据实验事件标记切分数据 - **基线校正**：消除直流偏移和慢波漂移 > 提示：对于实际研究，建议使用至少64导的EEG系统，并确保电极阻抗低于5kΩ。采样率应不低于250Hz以捕捉高频振荡。 ## 2. 计算功能连接矩阵 功能连接矩阵是构建脑网络的基础，它量化了不同脑区之间的交互强度。常用的连接指标包括： | 指标类型 | 计算方法 | 方向性 | 适用场景 | |---------|---------|--------|---------| | 相干性 | 频域相关性 | 无向 | 研究特定频段的同步活动 | | 相位锁定值 | 相位差一致性 | 无向 | 神经振荡耦合分析 | | Granger因果 | 时间序列预测 | 有向 | 信息流向推断 | 下面以Granger因果为例展示计算过程： ```python from mne_connectivity import spectral_connectivity_epochs # 定义要分析的频段 freq_bands = {'theta': (4, 7), 'alpha': (8, 12), 'beta': (13, 30)} # 计算Granger因果 conn = spectral_connectivity_epochs( epochs, method='gc', fmin=4, fmax=30, gc_n_lags=20, # 自回归模型阶数 n_jobs=1 ) # 获取连接矩阵 (n_nodes × n_nodes × n_freqs) conn_matrix = conn.get_data() ``` 连接矩阵的可视化对于理解脑区交互模式很有帮助： ```python import matplotlib.pyplot as plt # 平均所有频率的连接强度 mean_conn = conn_matrix.mean(axis=-1) plt.figure(figsize=(10, 8)) plt.imshow(mean_conn, cmap='hot', interpolation='nearest') plt.colorbar(label='Granger Causality') plt.title('Average Functional Connectivity Matrix') plt.xlabel('EEG Channels') plt.ylabel('EEG Channels') plt.show() ``` ## 3. 构建复杂脑网络 获得连接矩阵后，我们可以将其转化为图结构进行分析。NetworkX库提供了丰富的图论算法实现。 首先需要确定网络构建的阈值策略： - **绝对阈值**：保留强度超过固定值的连接 - **相对阈值**：保留前x%最强的连接 - **最小生成树**：保证网络连通性的最简结构 ```python import networkx as nx import numpy as np def threshold_matrix(matrix, percentile=90): """保留前percentile%的连接""" threshold = np.percentile(matrix, percentile) binary = matrix > threshold return matrix * binary # 应用阈值 thresholded = threshold_matrix(mean_conn, 85) # 构建有向加权图 G = nx.from_numpy_array(thresholded, create_using=nx.DiGraph) # 添加节点标签 labels = {i: ch_name for i, ch_name in enumerate(epochs.ch_names)} nx.set_node_attributes(G, labels, 'label') ``` 网络可视化可以帮助我们直观理解连接模式： ```python pos = nx.spring_layout(G, seed=42) # 定义节点位置 plt.figure(figsize=(12, 10)) nx.draw_networkx_nodes(G, pos, node_size=500, node_color='skyblue') nx.draw_networkx_edges(G, pos, width=1.5, edge_color=[d['weight'] for u,v,d in G.edges(data=True)], edge_cmap=plt.cm.Reds) nx.draw_networkx_labels(G, pos, labels, font_size=10) plt.title('Brain Network Topology') plt.axis('off') plt.show() ``` ## 4. 图论指标计算与分析 图论提供了丰富的指标来量化网络特性。以下是几个核心指标及其神经科学意义： ### 4.1 节点中心性分析 识别网络中的关键枢纽节点： ```python # 计算加权度中心性 degree_centrality = nx.degree_centrality(G) # 计算特征向量中心性 eigen_centrality = nx.eigenvector_centrality(G, max_iter=1000) # 结果排序 top_nodes = sorted(degree_centrality.items(), key=lambda x: x[1], reverse=True)[:5] print("Top hub nodes:") for node, score in top_nodes: print(f"{labels[node]}: {score:.3f}") ``` ### 4.2 全局网络指标 评估整体网络特性： ```python # 转换为无向图计算全局指标 G_undir = G.to_undirected() metrics = { '平均聚类系数': nx.average_clustering(G_undir), '全局效率': nx.global_efficiency(G_undir), '小世界性': nx.sigma(G_undir, niter=10, seed=42), '模块化': nx.algorithms.community.modularity( G_undir, nx.algorithms.community.greedy_modularity_communities(G_undir)) } print("\n全局网络指标:") for name, value in metrics.items(): print(f"{name}: {value:.3f}") ``` ### 4.3 社区检测 识别功能模块： ```python from networkx.algorithms import community # 使用Louvain算法检测社区 partition = community.louvain_communities(G_undir, seed=42) print("\n检测到的社区:") for i, comm in enumerate(partition): print(f"Community {i+1}: {[labels[n] for n in comm]}") ``` ## 5. 实际应用与注意事项 将这套方法应用于实际研究时，有几个关键点需要考虑： - **数据质量**：EEG信号的信噪比直接影响连接估计的准确性 - **参数选择**：Granger因果的模型阶数、网络构建的阈值都需要合理设置 - **统计验证**：需要通过置换检验等方法评估结果的显著性 - **生理解释**：图论指标的神经生物学意义需要结合具体研究问题 一个完整的分析流程可能还需要： ```python # 组间比较示例 def compare_groups(group1_conn, group2_conn, n_perm=1000): """使用置换检验比较两组网络指标""" # 实现细节省略 pass # 动态网络分析 def sliding_window_analysis(raw, window_size=1.0, step=0.2): """计算时变网络特性""" # 实现细节省略 pass ``` 在实际项目中，我发现最常遇到的挑战是参数选择。例如Granger因果的模型阶数，太小会导致模型欠拟合，太大则会引入噪声。经过多次实验，我发现对于采样率250Hz的EEG数据，20-30阶通常能取得较好平衡。另一个实用技巧是在网络可视化时，使用脑地形图作为背景，这样能更直观地理解连接模式的解剖基础。