步骤 07:坏导检测(第一次)

import mne
import numpy as np

# ====== 方法 A：LOF（Local Outlier Factor） ======
# mne-bids-pipeline 默认的自动检测算法
# 真实函数: mne.preprocessing.find_bad_channels_lof（New in MNE 1.7）
from mne.preprocessing import find_bad_channels_lof

bad_channels = find_bad_channels_lof(
    raw,
    n_neighbors=20,
    threshold=1.5,         # 阈值，越小越激进
    picks=['eeg'],
    verbose='INFO',
)
print(f"LOF 检测到的坏导: {bad_channels}")

# ====== 方法 B：基于功率谱的自动检测（无需额外依赖） ======
# 适合超平坦或超噪通道
psd = raw.compute_psd(picks='eeg')
psd_data = psd.get_data(power=True)
log_power = 10 * np.log10(psd_data)

# 自定义阈值（dB）
log_power_mean = log_power.mean(axis=1)
median_power = np.median(log_power_mean)
mad = np.median(np.abs(log_power_mean - median_power))
threshold_high = median_power + 5 * mad
threshold_low = median_power - 5 * mad

flat_channels = [psd.ch_names[i] for i, p in enumerate(log_power_mean) if p < threshold_low]
noisy_channels = [psd.ch_names[i] for i, p in enumerate(log_power_mean) if p > threshold_high]
print(f"平导: {flat_channels}, 噪导: {noisy_channels}")

# ====== 方法 C：交互式人工标记（推荐！） ======
# 弹出 GUI 窗口，点击通道名可标记坏导
raw.plot(block=True, bad_color='red')

# ====== 应用：把检测结果写入 raw.info['bads'] ======
raw.info['bads'] = bad_channels  # 加入自动检测结果
# raw.info['bads'].extend(flat_channels)  # 手动补充
print(f"当前坏导: {raw.info['bads']}")

# 验证 1：打印坏导列表
print(f"坏导数: {len(raw.info['bads'])}/{len(raw.ch_names)}")
print(f"坏导列表: {raw.info['bads']}")

# 验证 2：坏导占比应在 5%~15%
bad_ratio = len(raw.info['bads']) / sum(1 for ch in raw.info['chs'] if ch['kind'] == 2)
assert 0 <= bad_ratio <= 0.20, "坏导比例异常，请人工复核"

# 验证 3：可视化（坏导在 plot 中会被标红或标灰）
raw.plot(duration=5, n_channels=30)

# 验证 4：检查 raw.info['projs'] 是否仍为 []
# 坏导检测本身不投影，仅标记
assert len(raw.info['projs']) == 0, "本步骤不应创建 SSP 投影"