5.1 原始数据预处理：从谱图到峰表

5.1 原始数据预处理：从谱图到峰表

预处理把每个样品上百 MB 的三维数据（保留时间 × m/z × 强度）压缩成"行=样品、列=特征峰"的矩阵。四步：**峰检测、对齐归并、缺失填补、归一化**。这一步的质量决定后面一切分析的上限。

峰检测：在噪声里找信号

以 R 的 xcms（LC-MS 事实标准）为例，centWave 算法在 m/z 维找高斯形色谱峰：

library(xcms)
raw <- readMSData(filenames, mode = "onDisk")
cwp <- CentWaveParam(peakwidth = c(5, 30),   # 峰宽秒数，宽了窄了都丢峰
                     ppm = 25,               # 质量精度，按仪器规格设
                     snthresh = 10)          # 信噪比门槛
xset <- findChromPeaks(raw, param = cwp)
xset <- adjustRtime(xset, ObiwarpParam())    # 保留时间校正
xset <- groupChromPeaks(xset,                # 跨样品对齐归并
        PeakDensityParam(sampleGroups = groups, minFraction = 0.7))
​

minFraction = 0.7 的含义：一个特征至少要在 70% 的样品里出现，否则视为噪声丢弃。调低它，矩阵变宽、假峰变多；调高它，真实但低丰度的代谢物被误杀。没有万能值，要看 QC 中该特征的重复性再定。

填补：缺失不是随机事件

代谢组学数据的缺失有明确机制：浓度低于检出限。这意味着不能用"零"或随机小数胡填——常用 kNN 或最小值一半填补，并给填补比例设上限：

import pandas as pd
from sklearn.impute import KNNImputer

feat = pd.read_csv("peak_table.csv", index_col=0)
missing = feat.isna().mean(axis=0)
keep = missing < 0.30            # 缺失超30%的特征直接不要
feat = feat.loc[:, keep[keep].index]
feat[:] = KNNImputer(n_neighbors=5).fit_transform(feat)
​

归一化与缩放：两个不同的问题

样品间归一化处理的是进样量/提取效率差异：概率商归一化（PQN）、中位数归一化、或用内标。变量间缩放处理的是代谢物浓度量级差异（葡萄糖 mmol 级 vs 激素 pmol 级）：pareto 或自 autoscaling。别混用两个词，论文方法学部分写清楚用的哪个。

常见缩放方式对比

⚠️ 常见坑：把归一化和缩放做反了顺序，或在 QC 校正（3.3 节 LOESS）之前就缩放——流程顺序是漂移校正 → 缺失填补 → 样品归一化 → 变量缩放，颠倒会互相污染。

本节要点回顾

• 四步流水线：找峰 → 对齐 → 填补 → 归一化
• 参数没有默认正确值：峰宽、ppm、minFraction 要按仪器与 QC 表现调
• 缺失机制是检出限，kNN 填补 + 30% 缺失上限
• 归一化（样品间）与缩放（变量间）是两件事，顺序不能乱

下一节给这些匿名特征峰起名字。

峰表生成的五步流水线与参数敏感性

格式转换（ProteoWizard 转 mzML）→ 峰提取（XCMS 的 centWave）→ 比对归组（grouping）→ 保留时间校正（retcor）→ 填补缺失（fillPeaks），这五步的参数对结果的影响常大于生物学差异本身。centWave 的两个关键参数 ppm 与 peakwidth 的设置不当，要么把噪声当峰（假峰倍增），要么把宽峰切碎（一个代谢物变三个特征）。参数评估的客观方法是 pooled QC 的重复性：调参后重算 QC 样本 CV 中位数，取使其最小的参数组，而不是沿用文献默认值。

缺失值处理：不是统计小节而是化学问题

# 三种缺失机制对应三种处理：删、填、分层建模
import numpy as np
# 模拟某代谢物在 60 个样本中的强度（含 15% 缺失）
raw = np.random.default_rng(3).normal(1e5, 1e4, 60)
miss = np.random.default_rng(3).random(60) < 0.15
x = raw.copy(); x[miss] = np.nan
print(f"缺失率 {np.isnan(x).mean()*100:.0f}%")
# 方案A：最小值/5 填充（假定缺失=低于检测限）
x_a = np.where(np.isnan(x), np.nanmin(x)/5, x)
# 方案B：kNN 填充（假定缺失=随机机制）
from math import isnan
obs = [v for v in x if not isnan(v)]
x_b = np.where(np.isnan(x), np.median(obs), x)   # 单变量简化版
print(f"方案A 均值 {np.mean(x_a):.3e}（偏低，适合左删失）")
print(f"方案B 均值 {np.mean(x_b):.3e}（保分布，适合随机缺失）")
​

代码里两种填充的差异看似微小，但在"组间缺失率不同"的场景（比如病例组整体浓度更低、更多样本低于检测限）会直接制造假阳性——病例组被系统性填低。判断缺失机制的第一步永远是看"缺失率是否与组别相关"，相关的缺失必须回到仪器层面解决（降低稀释倍数、换更灵敏方法），而不是统计层面掩盖。

标度化方法的选择同样有纪律：autoscating（均值中心化除以标准差）让低浓度高方差的噪声代谢物获得与大代谢物同等的话语权，可能放大噪声；pareto（除以标准差的平方根）是代谢组学折中主流；log 变换先把乘性误差转为加性、再做标度，适合浓度跨数量级的全谱数据。顺序不可颠倒：先缺失处理、再标度化、最后建模。审稿中最常见的低级错误是先标度化再填补缺失——填补值在标度后的坐标系里失去物理含义，全部下游分析被污染。