### 如何使用Python进行EEG微状态分析 #### 创建和预处理EEG数据集 为了执行EEG微状态分析，首先需要加载并预处理原始EEG数据。这通常涉及读取文件、设置参考电极以及过滤信号噪声。 ```python import mne from mne.datasets import sample data_path = sample.data_path() raw_fname = data_path + '/MEG/sample/sample_audvis_raw.fif' raw = mne.io.read_raw_fif(raw_fname, preload=True) # 设置新的参考电极 raw.set_eeg_reference('average') # 使用平均重参考[^2] # 应用带通滤波器去除不需要的频率成分 raw.filter(1., 40.) ``` #### 提取全局场功率（Global Field Power, GFP） 微状态分类依赖于识别由GFP定义的时间点上的局部最大值。这些时间点对应着不同大脑活动模式下的瞬态稳定期。 ```python evoked = raw.average() # 计算GFP作为通道间方差的标准偏差 gfp = evoked.std(axis=0).flatten() times = evoked.times * 1e3 # 时间轴转换成毫秒表示法 plt.figure(figsize=(10, 6)) plt.plot(times, gfp, label='GFP') plt.xlabel('Time (ms)') plt.ylabel('Amplitude ($\mu V$)') plt.title('Global Field Power over Time') plt.legend() plt.show() ``` #### 定义微状态类别并通过聚类算法确定最优数量 通过K-means或其他无监督学习技术对基于GFP峰值时刻的地图进行分组，从而找到最佳的微状态数目。 ```python from sklearn.cluster import KMeans import numpy as np # 获取所有GFP峰位置处的空间分布矩阵 peak_indices = np.where(gfp >= np.percentile(gfp, 95))[0] maps_at_peaks = evoked.data[:, peak_indices].T kmeans = KMeans(n_clusters=4, random_state=0).fit(maps_at_peaks) labels = kmeans.labels_ centers = kmeans.cluster_centers_ print(f'Found {len(set(labels))} distinct microstates.') ``` #### 可视化微状态地图 一旦确定了各个微状态，则可以绘制它们各自代表性的头皮拓扑图以便直观理解每种状态下活跃区域的位置。 ```python for i, center in enumerate(centers): fig, ax = plt.subplots() im = ax.imshow(center.reshape((10, 11)), cmap='RdBu_r', interpolation='nearest', extent=[-np.pi, np.pi, -np.pi / 2., np.pi / 2.]) ax.axis('off') cbar = fig.colorbar(im, shrink=.8) cbar.ax.set_ylabel('Voltage (\u00B5V)', rotation=-90, va="bottom") plt.title(f'Microstate Map #{i+1}') plt.show() ``` #### 分析微状态序列特性 最后一步是对检测到的一系列微状态随时间变化的趋势做进一步探索，比如持续时间和发生概率统计等。 ```python transitions = [] current_label = labels[0] start_time = times[peak_indices][0] for idx, lbl in zip(peak_indices[:-1], labels[:-1]): next_idx = peak_indices[idx+1] if lbl != current_label or abs(next_idx-start_time)>10: # 超过10 ms视为新事件 transitions.append([current_label, start_time, times[next_idx]]) current_label = lbl start_time = times[next_idx] transition_df = pd.DataFrame(transitions, columns=['State', 'Start(ms)', 'End(ms)']) print("Transition DataFrame:") print(transition_df.head()) ```