第2章 质谱数据读取与解析：用Python搞定那些烦人的格式

说起质谱数据格式，你想想看，各家仪器厂商都有自己的“方言”。Thermo的.raw，Bruker的.d，Agilent的.d（对，又是.d，但结构完全不同），还有通用的mzML、mzXML、mgf……当年我第一次接手多平台数据时，差点被这些格式搞疯掉。

这一章，我来聊聊怎么用pyOpenMS和pymzML这两个库，把这些格式统一读进Python。说白了，就是给机器学习模型当个靠谱的“翻译官”。

2.1 质谱文件格式的“家族图谱”

搞懂格式是第一步。我个人习惯把它们分成三类：

• 厂商原生格式：Thermo .raw、Bruker .d、Agilent .d、Waters .raw。封闭、二进制、压缩，只有厂商自己的工具能读。
• 开放交换格式：mzML、mzXML、mgf。基于XML或文本，社区标准，靠得住。
• 过渡/分析格式：mzTab、mzIdentML、featureXML。存的是分析结果，不是原始谱图。

你肯定会问：为什么不能统一？嗯，商业原因嘛，谁都不想轻易把自己的核心算法暴露出去。但好在社区给力，pyOpenMS能直接怼开大部分厂商格式，pymzML更偏轻量级，适合快速解析。

2.2 pyOpenMS：重武器，能打硬仗

pyOpenMS基于OpenMS C++库，性能好，功能全。我一般在处理大型队列、需要做图谱对齐或定量分析时用它。安装稍微有点门槛：

pip install pyopenms

# 如果Conda更方便
conda install -c bioconda pyopenms

2.2.1 读取mzML文件

mzML是HUPO-PSI推荐的标准格式。读取流程很直白：

import pyopenms as oms

# 加载文件
exp = oms.MSExperiment()
oms.MzMLFile().load("sample.mzML", exp)

print(f"谱图数量: {exp.size()}")
print(f"光谱范围: {exp.getMinRT():.2f} - {exp.getMaxRT():.2f} 秒")
print(f"质荷比范围: {exp.getMinMZ():.2f} - {exp.getMaxMZ():.2f} Th")

我在项目中遇到过一个问题：仪器采集时可能记录了大量低质量谱图（比如空白、冲洗段），直接读进来会拖慢后续分析。我的办法是先筛一遍保留时间：

import numpy as np

# 获取所有谱图
spectra = exp.getSpectra()

# 只保留RT在60-600秒之间的（假设这是有效分析窗口）
valid_spectra = [s for s in spectra if 60 < s.getRT() < 600]

print(f"剔除掉 {len(spectra) - len(valid_spectra)} 张无效谱图")

2.2.2 从谱图里取出元数据和峰值列表

每张谱图都有两个核心东西：元数据（保留时间、MS级别、极性等）和峰值列表（m/z, 强度对）。提取方法如下：

def extract_spectrum_info(spectrum):
    """提取单个谱图的核心信息"""
    info = {
        "rt": spectrum.getRT(),
        "ms_level": spectrum.getMSLevel(),
        "polarity": spectrum.getPolarity(),
        "tic": sum(spectrum.getIntensityArray()),
        "peaks": list(zip(
            spectrum.getMZArray(),
            spectrum.getIntensityArray()
        ))
    }
    return info

# 遍历所有谱图
for i, spec in enumerate(spectra[:5]):  # 先看前5个
    info = extract_spectrum_info(spec)
    print(f"谱图 #{i}: MS{info['ms_level']}, RT={info['rt']:.2f}, 肽段数={len(info['peaks'])}")

2.2.3 读mzXML和mgf

mzXML和mzML结构类似，只是标签名不同。pyOpenMS处理起来几乎一样：

# 读mzXML
exp = oms.MSExperiment()
oms.MzXMLFile().load("sample.mzXML", exp)

# 读mgf（最常见的是PMF或MS/MS谱库）
exp = oms.MSExperiment()
oms.MascotGenericFile().load("sample.mgf", exp)

要注意的是，mgf文件不一定包含保留时间信息。有些来源的mgf只存m/z和强度，需要额外从文件名或元数据里补全。

2.3 pymzML：轻量快跑，适合做特征提取

pymzML是纯Python写的，安装零依赖，解析原理直观。我更喜欢在快速原型验证、小批量数据探索时用它。安装极简：

pip install pymzml

2.3.1 读mzML并访问谱图

pymzML把谱图当迭代器，内存友好：

import pymzml

run = pymzml.run.Reader("sample.mzML")

# 遍历所有谱图
for i, spectrum in enumerate(run):
    if i > 5:  # 只处理前6张
        break
    
    print(f"谱图 #{i}: MS{spectrum.ms_level}, RT={spectrum.time:.2f}")
    print(f"  峰值数量: {spectrum.peaks.shape[0]}")
    
    # 提取前10个最高峰
    top_idx = spectrum.peaks[:, 1].argsort()[-10:][::-1]
    top_peaks = spectrum.peaks[top_idx]
    print(f"  最高峰: m/z={top_peaks[:3, 0]}, 强度={top_peaks[:3, 1]}")

你想想看，pymzML的peaks属性直接是numpy数组，省去了从列表转数组的麻烦，写后续特征计算时顺手得多。

2.3.2 获取元数据和离子谱

pymzML把元数据存成字典，按图索骥就行：

for i, spectrum in enumerate(run):
    if i > 5:
        break
    
    # 元数据
    meta = spectrum.extracted_meta
    print(f"极性: {meta.get('polarity', 'N/A')}")
    print(f"前体离子m/z: {spectrum.selected_precursors}")
    print(f"碰撞能量: {meta.get('collision_energy', 'N/A')}")
    
    # 如果是MS2，可以取出碎片谱
    if spectrum.ms_level == 2:
        precursors = spectrum.selected_precursors
        if len(precursors) > 0:
            precursor_mz = precursors[0][0]  # 第一个前体离子
            print(f"Fragments for {precursor_mz:.4f}")

2.4 常见格式对比：选哪个库？

我个人是混合着用——pymzML做特征提取的探索，pyOpenMS做最终的批量处理管道。各取所长嘛。

2.5 一次完整的读取与解析流程

我把上面讲的东西串起来，画成一个可复用的流程框架。后续特征构建里，你都会用到这个骨架。

这个流程你肯定要反复用。第3章讲具体的特征计算时，我会基于这个框架，在每个阶段给出代码实现。

2.6 避坑指南：我踩过的三个坑

既然分享了这些，顺便把踩过的坑也列出来，算是帮你省点时间。

1. 内存爆炸：别一次性把整个大文件读到内存里。用过pymzML的迭代器机制，加上分块处理，能省一半内存。
2. 乱序谱图：有些仪器采集时保留时间不严格递增。pyOpenMS的exp.sortSpectra(True)能自动补排序。不过我更建议你在读取时就按RT排好，免得后续搞错时间轴。
3. 缺失元数据：厂商文件转成mzML时，常有字段丢失。比如极性字段在Thermo .raw里是明码，但转成mzML后可能变成枚举值。我建议做一步元数据完整性检查，缺的字段手动补全。

好了，质谱数据读取这件事，说白了就这么回事——选对库，读对结构，别漏元数据。你想想看，后面要做的特征工程、模型训练，全都依赖这一步。把它理得清清爽爽，后面就省心多了。