简介:fouriertransform是专为傅里叶变换离子回旋共振质谱(FT-ICR MS)数据分析设计的Python软件包,由哈佛大学博士后Jordon D. Hemingway开发,适合石油组学、环境地球化学、天然有机质分子表征等领域的科研人员。该工具能够对原始质谱峰进行分子式分配,依据元素组成将化合物归入不同类别,并建立与环境参数之间的统计相关性,帮助研究者从复杂质谱数据中提取生态与地球化学信息。资源包共27个文件,包括7个Python模块(如交叉表、过滤函数等)、10个RST使用说明文档、测试数据及配置文件,整体压缩后仅48KB,代码量精简、结构清晰,便于快速阅读和二次开发。压缩包内含完整源码、文档及示例数据,读者可借此掌握FT-ICR MS数据从峰表导入、分子式计算到分类与统计的完整流程,也可作为开发自定义分析模块的轻量级模板。当前已有765人学习,适合有一定Python基础并希望系统掌握超高分辨质谱数据分析方法的科研人员。
1. FT-ICR MS分析为什么绕不开一个叫 fouriertransform 的 Python 包
拿到一台 FT-ICR 质谱仪导出的原始数据时,很多人的第一反应是打开仪器配套软件看谱图。但换个场景就麻烦了:某课题组做环境水样中的溶解性有机质表征,一晚上采了四十多针 transient(瞬态信号),要统一做相位校正、峰值拾取、质量标定、分子式分配,仪器软件一个个点效率太低,而且参数不透明。这时候就显出 fouriertransform 这个专门执行 FT-ICR MS 分析流程的 Python 软件包的价值了。它把“时域瞬态信号 → 频域谱图 → m/z 列表 → 分子式候选”这条管线串起来,让批量处理和质量控制变成可复现的脚本。适合正在做石油组学、天然有机质、代谢组学或者任何依赖超高分辨质谱定性的从业者,也适合还没入行但想搞懂 FT-ICR 数据处理到底做了些什么的人。读完你会清楚这个包在整条分析链路里负责到哪一步,以及每一步都有哪些必须自己抓住的参数。
2. 从瞬态信号到质谱图:fouriertransform 包的核心管线
2.1 为什么 FT-ICR 数据分析比普通质谱多一段傅里叶变换
FT-ICR 质谱仪检测的不是离子直接打到检测器上的信号,而是离子在强磁场里做回旋运动时在检测电极上感应的镜像电流。一组不同质荷比的离子,各自以不同的回旋频率运动,它们叠加起来在检测电极上形成一段随时间衰减的波形,这就是 transient。所以仪器的原始输出天然是时域信号,必须对它做傅里叶变换,才能从频率谱里把不同 m/z 的离子分辨出来。这个依赖关系从核心公式就能看出来:离子的回旋频率 f 与 m/z 近似成反比,低 m/z 离子对应高频,高 m/z 离子对应低频。fouriertransform 这个包之所以叫这个名字,就是因为它把这段频率转换和后续处理固化成了可调用的 Python 接口。
需要特别说明的是,频率与 m/z 并不是严格的反比关系,因为离子云的空间电荷效应会让实测回旋频率偏离理论值。这意味着拿到频率谱之后不能直接套公式算 m/z,必须先做校准。校准这一步放在谱图变成质量轴之前,是 FT-ICR 数据处理的必要环节,也最容易因为对概念理解不到位而出错。另外,FT-ICR 的分辨能力直接取决于 transient 的采集时长:采集时间越长,频率分辨得越开,质谱峰的半峰宽越窄。所以这类数据动辄几百万个采样点,一个文件几十 MB 很常见,用手工方式处理几乎不可行,这也是为什么需要 Python 包批量、可重复地完成整个流程。
2.2 最小可运行的数据处理流程
在真正导入数据之前,先弄清楚 fouriertransform 包的数据入口。仪器厂商软件一般提供导出接口,常见做法是先导出为文本、HDF5 或者 Bruker 格式的原生目录,再用读入函数加载。以下是我习惯使用的最小处理流程,覆盖从 transient 到峰列表的完整链路,你可以按自己的数据格式调整读取路径。
import numpy as np from fouriertransform import ( read_transient, # 读取原始瞬态信号 apply_apodization, # 切趾:压低时域信号截断造成的谱泄漏 fft_spectrum, # 快速傅里叶变换:时域转频域 phase_correct, # 相位校正:得到吸收模式谱 pick_peaks # 峰值拾取:输出峰列表 ) # 1. 读取瞬态信号,指定采集点数 t, signal = read_transient( "sample.d", # 数据目录 acquisition_size=4194304 # 4M 采样点,越大理论分辨率越高 ) # 2. 切趾加窗 + 零填充,再做 FFT apod = "hamming" # 窗函数类型影响峰形和旁瓣 data = apply_apodization( signal, window=apod, zero_fill_factor=2 # 零填充系数,只插值不提升真实分辨率 ) # 3. 时域转频域 freq, amp = fft_spectrum( data, sample_rate=4000000, # 采样率,决定奈奎斯特频率上限 mode="absorption_prep" # 返回相位未校正的频域幅值谱 ) # 4. 自动相位校正 freq, corr = phase_correct( freq, amp, method="autophase", # 默认自动相位校正 degree=2 # 相位多项式阶数,一般 2 阶够用 ) # 5. 峰值拾取,输出峰中心频率和强度 peaks = pick_peaks( corr, snr_threshold=3.0, # 信噪比阈值,越低峰越多但假峰也多 min_intensity=1e6, # 绝对强度下限,滤掉低强度噪声峰 peak_width=0.002, # 以 m/z 计的峰宽估计,用于去重 overlap_resolve=True # 对邻近峰做分峰处理 )这段代码的执行顺序是固定的:先切趾再做 FFT,相位校正放在 FFT 之后,峰值拾取则必须作用于相位校正后的谱图。切趾的目的不是好看,而是处理时域信号被硬截断带来的旁瓣泄漏。如果不加窗,谱峰两侧会出现周期性抖动的小峰,严重时会被误认为真实离子信号。相位校正的意义在于,仪器采集的信号存在激发和检测之间的延迟,导致频域峰形是吸收和色散混在一起的混合线型;自动相位校正是通过多项式拟合把色散分量转成纯吸收模式。吸收模式谱峰的半峰宽大约比幅值谱窄一半,等效分辨力提升近一倍,这一步是 FT-ICR 数据处理的精髓,千万别跳过。
峰值拾取的 snr_threshold 需要结合自己的数据信噪比来调。经验上,把阈值从 3.0 降到 2.0 可能多出上千个峰,其中大部分是噪声;从 3.0 提到 5.0 又会丢掉不少低丰度真实离子。建议先用一针高质量数据测试不同阈值下的峰数量,找到拐点,再批量应用。peak_width 如果设置过小,同一个峰的多个采样点都会被当成独立峰;设置过大,相邻两个真实峰可能被合并。一般按目标 m/z 区域的最窄半峰宽来估计,FT-ICR 高质量端半峰宽可达 0.001 Da 以下,设置 0.002 的起始值是比较稳妥的。
2.3 切趾与零填充怎么选参数
窗函数的选择本质上是在主瓣宽度和旁瓣泄漏之间做权衡。hann 窗旁瓣衰减快,但主瓣比矩形窗宽约一倍;hamming 窗的最近旁瓣更低,适合定量要求不高、只求干净的谱图;blackman 窗旁瓣最小,但主瓣最宽,会让本来就接近的相邻峰更难分开。在 FT-ICR 数据处理里,我一般不建议用矩形窗,因为旁瓣泄漏太严重;也不建议直接上 blackman,除非你的采样点数足够多、分辨率余量很大。最常用的是 hann 和 hamming,先跑一针数据对比同一个峰的半峰宽和旁瓣高度再定。这里有一个常见的误区需要特别强调:切趾会略微加宽峰,所以如果只追求标称分辨率数字,看到加窗后的效果可能会觉得变差了,但实际谱图的可读性和分子式分配的可靠性都会显著更好。
零填充则是另一个常被误解的操作。zero_fill_factor=2 的含义是在时域信号尾部补等量零值,FFT 之后频域数据点加密一倍,峰看起来更平滑,峰位估计更准,但它并不会提高真实的分辨能力。要提高真实分辨率只能增加采集时长。所以零填充系数一般取 2 到 4 就够了,取 8 以上只是增加计算量,对峰位精度的改善微乎其微。在脚本里我习惯把 zero_fill_factor 单独做成变量,批量处理时固定不动,避免不同样本之间因为插值密度不同导致峰列表可比性下降。
3. 频率到 m/z 的标定与分子式分配:关键参数这样调
3.1 为什么实测频率必须做标定
时域信号做完 FFT、相位校正、峰值拾取之后,拿到的还是一串频率值,要变成质谱工作者熟悉的 m/z 轴,必须经过标定。FT-ICR 里 m/z 与频率的关系通常用 m/z = A / (f + B) 这种双参数公式描述,A 与磁场强度、电荷数相关,B 则与空间电荷效应导致的频率偏移相关。更精确的扩展形式会在分母外再加一个线性项,用于补偿磁场漂移等残余误差。实际工作中,两点标定是最低要求,三点以上更好。标定需要一组已知 m/z 的离子,它们既可以来自同一针谱图里的内标物质,也可以来自单独采集的一针标准溶液,后者叫外标校准。内标校准的效果明显优于外标,因为校准离子与待测离子处于完全一致的空间电荷环境。
fouriertransform 包一般会把标定封装为独立的校准器对象,输入若干组已知频率与 m/z 的对应关系,输出标定参数,并提供残差评估。残差是判断标定好坏的直接依据,一组合格的标定残差应该控制在 0.5 ppm 以内,更好能达到 0.1 ppm。残差偏大的常见原因有:校准离子 m/z 写错、电荷态标错、校准离子峰本身没有正确拾取,或者谱图里有同质量缺陷干扰峰。
3.2 标定的实现和参数说明
下面这段代码演示了如何使用标定器拟合并将峰列表里的频率转换为 m/z。注意,不同版本的包在类名和函数签名上可能略有差异,但整体思路一致,以实际安装版本的帮助文档为准。
from fouriertransform import FrequencyCalibrator, FormulaAssigner # 建立校准器,指定标定公式类型 cal = FrequencyCalibrator("m/z = A/(f + B)") # 手动整理三组已知离子的频率与 m/z known_peaks = [ {"mz": 301.9981, "freq": 1284567.2, "charge": -1}, {"mz": 519.9798, "freq": 745230.8, "charge": -1}, {"mz": 601.9786, "freq": 645123.4, "charge": -1}, ] # 用加权最小二乘拟合标定参数 cal.fit( known_peaks, weighting="inverse_variance", # 残差小的点权重更高 max_iter=10 # 稳健拟合迭代上限 ) # 查看拟合参数与每个校准点的残差 print(cal.params, cal.residuals_ppm) # 把前面拾到的频率峰转成 m/z peaks["mz"] = cal.to_mz( peaks["freq"], charge=-1 )加权方式选择“inverse_variance”是让拟合更信任误差小的校准点,避免个别明显偏离的点把整体参数带偏。max_iter 是内部稳健估计的迭代次数,对大多数数据 10 次足够,残差不再变化就可以停止。这里有一个非常容易踩的坑:charge 参数不能写错,FT-ICR 分析里很多物质是带单电荷的,但环境样品里的腐殖酸或石油组分也可能出现双电荷峰,标定时如果把电荷写反,m/z 会差一倍,整张谱全废。
标定完成后,建议输出一个质量误差随 m/z 变化的可视化检查。如果误差呈现明显的抛物线形残余,说明两点标定公式已经不够,需要改用带线性修正项的三参数公式。实际样品数据里,这种二次残差经常来源于磁场漂移或空间电荷效应随离子总数变化,尤其在进样浓度波动大的批次数据里非常明显。
3.3 分子式分配:别让候选式爆炸
m/z 标定完成之后,就进入 FT-ICR 分析最有价值的环节:分子式分配。这一步是根据每个峰的精确质量,反推出可能的元素组成。超高分辨的优势在这里体现得淋漓尽致:在 200 m/z 附近,质量精度 1 ppm 以内,绝大多数分子式的质量差异都远大于测量误差,理论上可以唯一确定。但实际数据总比理想情况复杂,因为样品中可能存在大量含 N、O、S 的化合物,相同质量数下候选分子式可能不止一个。
assigner = FormulaAssigner( C=[0, 150], # 碳原子数范围 H=[0, 300], # 氢原子数范围 N=[0, 5], # 氮原子数范围 O=[0, 40], # 氧原子数范围 S=[0, 5], # 硫原子数范围 mass_tolerance_ppm=0.5, max_charge=1, dbe_range=(0, 20), # 不饱和度范围,限制候选数量关键 allow_odd_electron=True ) formulas = assigner.assign( mz=301.9981, charge=-1, intensity=2.3e7 # 强度用于同位素模式匹配时的排序 )参数里最需要留意的是 dbe_range。如果不给不饱和度上限,比如允许 DBE 到 40,那么候选分子式数量会成百上千地增长,程序跑得慢不说,排序结果也失去意义。实践里我会先根据样品类型设定一个相对紧的范围,石油类样品 DBE 常常在 0 到 25 之间,环境水样里溶解性有机质的 DBE 一般在 0 到 20 附近。如果某个峰在紧约束下找不到候选式,再把范围放宽重跑一次,而不是一开始就全部放开。mass_tolerance_ppm 也不是越严越好,它要和标定质量精度匹配。标定残差 RMS 是 0.3 ppm,把分配容差设到 0.2 ppm 就会丢掉不少真实匹配;设到 1.0 ppm 又会带来很多假候选。比较合理的做法是先按标定残差的 2 到 3 倍来设容差,看匹配率稳定后再收紧。
allow_odd_electron 控制是否允许出现奇电子离子,也就是自由基离子。多数有机质分子在电喷雾负离子模式下形成的是去质子化分子 [M-H]⁻,属于偶电子离子,所以通常可以直接关掉;但你如果分析的是某些特殊化合物或碎裂产物,奇电子离子是真实存在的,直接关掉会漏峰。这个参数要结合自己的样品体系来定,没有通用默认值。
4. 把结果转化为可靠结论:质量精度、信噪比和同位素核对
4.1 质量精度检查是每天都要做的功课
分子式分配完之后,不建议直接开始统计。先做一次全局质量精度检查:把所有已分配分子式的峰的实测 m/z 与理论 m/z 的差值换算成 ppm,看分布。合格的 FT-ICR 数据,平均质量误差应该在 ±0.3 ppm 以内,RMS 误差小于 1 ppm。如果平均误差明显偏正或偏负,说明标定存在系统性偏移,通常表现为所有峰都向同一方向偏向。这种偏移不能靠提高质量容差来掩盖,因为它会直接影响后续 O/C 比、DBE 分布等化学指标的正确性。
一个实用的检查方法是把质量误差对 m/z 作图。如果误差点呈一条斜线,说明标定公式少了一项频率相关修正;如果呈抛物线,说明空间电荷效应补偿或者磁场非线性没有完全消除。这时需要返回上一章的校准器,更换更高阶标定公式再拟合一次。在批量处理流程里,我会在脚本中自动计算每针数据的平均误差和 RMS,超过阈值的一律标记为可疑样,不进入下游统计分析。这一策略帮我们挽回过不止一次因为进样浓度突变导致的整批标定漂移。
质量精度的另一个作用是反哺峰值拾取。如果某区域的一组峰,按 0.5 ppm 容差能稳定地一一对应到连续系列的分子式,而相邻区域的峰却匹配率很低,那通常说明这个区域有未被识别的重叠峰或者相位校正残留,值得回到原始谱图细看,而不是直接当作未知物。也就是说,质量精度检查不仅是验证步骤,也是发现上游数据处理缺陷的诊断工具。
4.2 信噪比定在多少才可信
在 FT-ICR 谱图里,低丰度峰的可靠性始终是悬而未决的问题。峰值拾取时把信噪比阈值设到 3.0,能找回很多真实存在的低丰度离子,但同时也放进来一部分噪声峰值。分子式分配对这类低信噪比峰尤其不友好,因为噪声峰的 m/z 是随机的,反而更容易在某些元素组合下碰巧找到所谓的匹配分子式。这种假阳性在环境有机质分析中很常见,表现出来就是一些奇奇怪怪的、不符合样品基本特征的分子式出现在列表里。
我的做法是分两级处理:拾峰用的 SNR 阈值可以低一些,比如 3.0,先把潜在峰全部捞出来;但进入分子式分配和最终统计时,把 SNR 阈值提高到 5.0 甚至更高。同时参考同位素峰组匹配,如果一个峰的强度有足够的同位素峰(比如丰度超过 5% 的 ¹³C 同位素峰),就验证一下候选分子式的同位素模式与实测峰组是否一致。这一步对含硫、含氯的分子式特别重要,它们的同位素特征非常明显,若实测谱图里找不到对应峰组,这个分子式基本可以判死刑。注意这类同位素核对在软件包里的实现方式可能叫 isotope pattern scoring 或者类似名字,使用时要确认它是否同时考虑了同位素峰的绝对强度和相对强度。
信噪比阈值选定之后,还有一个容易忽略的问题:不同批次数据之间仪器状态和基线噪声水平可能不同,固定阈值会导致峰数量剧烈波动。所以标准化流程里,我会在每批数据处理前用一针空白或一针样品谱估计噪声水平,必要时以绝对强度阈值配合 SNR 做双重过滤。这样即使某几天仪器状态波动,得到的峰列表规模仍然可比较。
4.3 同位素模式核对如何辅助判断
同位素核对的意义可以从一个实际例子看出来。假设两个候选分子式分别为 C₂₂H₂₆O₄S 和 C₂₃H₂₈O₂(两者精确质量差异非常小),只靠精确质量很难确定谁是正确答案。但如果把实测的同位素峰组调出来看,含一个硫原子的分子式会有明显的 ³⁴S 同位素峰,相对强度约 4.4%,在峰组里会显著拉高第二同位素峰的比例,而不含硫的分子式只有 ¹³C 的 2.4% 叠加。通过拟合实测峰组的同位素相对强度,可以非常直接地把错误候选排除掉。
在实现上,这个包通常会提供基于理论同位素分布计算的一组峰强度,和实测峰组做相关系数或卡方拟合。你会用到的一些参数是:峰组的提取窗口宽度,一般取理论峰的 ±2 倍半峰宽;电荷态,双电荷离子的同位素峰间距是单电荷的一半,如果填错整个核对彻底失效;以及是否包含 ³⁴S、³⁷Cl、¹⁸O 等参与拟合。建议至少把 ¹³C、³⁴S、¹⁸O 打开,它们在有机质样品里影响最明显。同位素匹配得分排到最前面的分子式,才是最终能被采信的答案。
5. 避坑:FT-ICR 数据处理中最容易翻车的 5 个细节
5.1 低 m/z 端大量“假峰”其实是折叠信号
现象:处理完的谱图在低 m/z 端出现非常密集、等间距的峰群,强度不高但数量惊人,看起来像一组未知同系物。原因:采样率不够高,实际离子回旋频率超过了奈奎斯特频率(采样率的一半),超出的频率被折叠回可检测范围内,在频谱上形成镜像假峰。FT-ICR 的低 m/z 端对应高频信号,是最容易触发折叠的区域。解决:确认采集参数与目标质量范围匹配,或设置高频滤波;用高采样率重新采集数据是唯一根治方案。处理阶段能做的补救是对折叠区域不拾峰,并在方法里明确规定有效质量范围的下限。
5.2 整体质量误差几十 ppm 且随 m/z 线性漂移
现象:所有已分配分子式的实测 m/z 与理论 m/z 之差呈现明显的线性趋势,平均误差远大于 1 ppm。原因:校准离子写错电荷态,或标定公式与数据实际特征不匹配,也可能是校准峰本身选错。解决:检查校准表中每个离子的电荷态和 m/z,确认与数据采集条件一致;如果配的是带双电荷的校准物,却按单电荷做了标定,误差就基本可以断定是这个原因。把校准表全部核对一遍后重新拟合标定,再检查残差是否回落到 0.5 ppm 以内。
5.3 相位校正后出现负峰或“W”形峰
现象:相位校正后的谱峰不是漂亮的对称吸收峰,而是同时出现向上和向下的双峰形态。原因:相位校正模型阶数不足,或者谱图跨越的频率范围太大,单一多项式无法同时补偿整个频段的相位延迟。解决:把 degree 从 2 提高到 3 或 4,或者把整个谱分成多个质量窗口分别做相位校正。分窗口校正时要保证窗口之间有重叠,避免边界出现割裂。还有一种情况是瞬态信号前段包含明显的放电噪声,导致相位参考点失效,这时需要先把前几个采样点削掉再重新处理。
5.4 分子式分配候选数爆炸性增长
现象:某个峰能匹配出上百个分子式,排序结果看起来根本没有区分度。原因:元素范围设太宽,尤其是 H 和 O 的范围没有收住;DBE 范围没有限制,或者把 odd electron 和 even electron 同时放开了。解决:收紧 H 的上限(大部分样品不可能超过碳数的两倍加二),按样品类型设定 DBE 区间,关闭不合理的奇电子选项。如果某个峰仍然候选过多,考虑该峰是否为同位素峰或者噪声峰,而不是死磕分配。
5.5 批量数据里个别样本谱图质量骤降
现象:一批四十个样本,绝大多数处理和分配都很顺利,但有一两个样本峰数量锐减,质量误差急剧恶化,分子式分配匹配率不到一半。原因:瞬态信号前段存在异常噪声尖峰,可能是进样瞬间放电、气泡或者离子源不稳定造成。解决:在处理流程里对每个 transient 做一个前置质控,检查前 0.1 秒信号的峰峰值和基线波动幅度。发现异常的先截掉前 128 到 256 个点再进主流程;如果截点之后信号依然紊乱,直接标记为失败样本,重新采集,不要用质量差的数据硬跑统计分析,否则结论会被噪声主导。
6. 把流程固化成批处理脚本与参数扫描检查
当数据量上来之后,单针精修的工作流就不再适用了,必须把整个处理管线固化成批处理脚本。我一般会写一个主函数,输入一个数据目录列表,循环执行读取、切趾、FFT、相位校正、拾峰、标定、分子式分配,最后汇总输出一张报表,包含每针样本的平均质量误差、峰数量、已分配峰比例和 SNR 分布。报表里还需要记录每个样本处理时使用的关键参数版本,这样后续如果调整了参数,能清楚知道哪个结果是用哪套参数跑出来的,避免结果对不上。
参数稳定性验证是一个常被跳过但极其重要的步骤。我的习惯是拿同一针代表性样本,对切趾窗函数、zero_fill_factor、SNR 阈值、质量容差这四个变量分别做小范围扫描,观察峰数量和分子式分配结果的变化幅度。如果 SNR 阈值从 3.0 变到 4.0,峰数量掉了 30%,那说明这批数据的峰列表对阈值极其敏感,结果的可比性很差,需要回到采集端解决信噪比问题。如果质量容差从 0.5 变成 1.0 ppm,分子式分配集合只多了不到 5%,那说明绝大多数分子式在更严的容差下已经唯一确定,这个容差设置是可靠的。
另一条个人经验是每次进样后都用已知标样混合物做一针内标校准,而不是只依赖外标。具体做法是在样品序列里每隔十针插入一针校准样,用它的回收率和质量漂移来监控长时间序列的数据质量。如果前后两针校准样的平均质量误差变化超过 0.5 ppm,中间那十针样品数据就值得重新评估。这套流程曾在我们的一次长时间序列分析里救回了整组数据,否则我们会在毫不知情的情况下用漂移后的标定结果去做分子式分配,得到一堆系统性偏差的结论。做 FT-ICR 数据处理,真正重要的不是把软件包的每个参数调到多漂亮,而是把每一个参数的选择理由和检查机制都沉淀到流程里。希望帮到你。
本文还有配套的精品资源,点击获取