简介:这是一套面向R语言用户的代谢组学工具包,专用于处理UK Biobank中Nightingale NMR生物标志物数据,适合从事大规模队列生物标志物分析的研究人员。包内提供 extract_biomarkers() 和 compute_nightingale_ratios() 两个核心函数,前者可从原始数据集中提取NMR字段并简化列名,按基线与首次重复评估拆分度量;后者能基于168种生物标志物计算临床常用比率,并附带 nmr_info 数据框供查询标志物信息。压缩包共含19个文件,以R源码、Rd帮助文档、MD说明文档及一个Rda数据文件为主,整体仅18KB,涵盖完整函数实现、使用说明和GitHub安装指引,可直接安装到R环境中使用。已有909人学习,对希望快速在UK Biobank中应用Nightingale数据的R用户而言,可直接获得工具源码与文档,省去自行整理字段和计算比率的时间,也可根据源码进一步扩展分析流程,便于二次开发与学习参考。
1. 拿到 UK Biobank 的 Nightingale 数据包之后:为什么我第一时间装 ukbnmr
申请 UK Biobank 数据时,很多人会在代谢组学那一栏勾上 Nightingale NMR 数据,下载解压后却在满屏字段号里发懵:几百列全是 23400 开头的数字,列名和 SHOWCASE 手册对不上,默认查看器打开文件就卡死,想跑 GWAS 或孟德尔随机化,连哪几列是 QC 标记都要翻半天文档。ukbnmr 正是为这个场景设计的 R 包:它内置了 Nightingale 生物标志物与 UK Biobank 字段的完整映射、样本级 QC 标记和检测限阈值,一条命令就能把原始数据包变成可以直接做分析的对象。对做代谢流行病学、全基因组关联分析和 MR 的研究者来说,这几乎是目前处理这批数据最省事的入口。
2. 先看懂数据包:Nightingale NMR 生物标志物在 UK Biobank 里怎么组织
在做任何分析之前,先花十分钟搞清楚你从 UK Biobank 下载回来的数据包里到底有什么。磨刀不误砍柴工,这一步省了,后面全是玄学报错。
2.1 解压后的三类配套文件:sbin、csv、sample 各管什么
UK Biobank 的标准数据包通常不是单个文件,而是一组配套文件放在同一个目录里。最常见的三类是 .sbin、.csv 和 .sample。
.sbin 是二进制数据,按行存储所有参与者的所有字段,体积大但读取效率高;.csv 其实是制表符分隔的文本,是同一份数据的文本形态,方便抽样查看;.sample 是样本清单,每行对应一个参与者,带 FID、IID 这类标识列,也包含性别、缺失率等元信息。这三者描述的是同一批样本,但行顺序不一定完全一致——这是后续最容易翻车的地方。
ukbnmr 的 build 函数按「整个解压目录」来读数据。你只需要把路径指到最外层文件夹,它会自动找到配套文件并完成关联。即使你申请数据时只拿到了文本版,没有 .sbin,也能读,只是加载速度会慢一些。我自己的习惯是:拿到任何 UK Biobank 数据包,先不急着打开 csv,而是确认目录下有这几类文件,再交给 ukbnmr 处理。
注意:无论文件用什么后缀,核心是把同一批样本的标识列(eid)找出来。ukbnmr 依赖 eid 列做关联。如果你手里的数据是从其他途径导出、没有标准 eid 列,先自己补一列再交给你后面要用的处理逻辑。
2.2 Nightingale 的 250+ 指标不是乱放的:六大类与字段编号规律
Nightingale NMR 平台测的是血浆或血清中能被核磁共振稳定定量的小分子和脂蛋白。UK Biobank 把结果按字段编号存进数据包,字段集中在 23400 开头的一段里。很多第一次接触的人都会问:几百列全是数字,怎么知道谁是谁?
按测量对象分,这批指标大致归为六大类。第一类是脂蛋白亚类,按密度和粒径拆成 VLDL、LDL、HDL 的子类,每个子类又分别报告颗粒浓度、胆固醇、甘油三酯、磷脂等不同表示;第二类是脂肪酸,包括总量、饱和度(SFA、MUFA、PUFA)和具体脂肪酸如 DHA 的占比;第三类是氨基酸,包括丙氨酸、谷氨酰胺、甘氨酸、组氨酸、异亮氨酸、亮氨酸、缬氨酸、苯丙氨酸、酪氨酸;第四类是糖酵解相关代谢物,如葡萄糖、乳酸、丙酮酸、柠檬酸;第五类是酮体,如乙酰乙酸、β-羟基丁酸、丙酮;第六类是炎症相关指标,最常用的是糖蛋白乙酰化,也就是 GlycA。
| 类别 | 测量对象 | 指标示例 |
|---|---|---|
| 脂蛋白亚类 | VLDL/LDL/HDL 浓度、胆固醇、甘油三酯、磷脂、粒径 | LDL_P、HDL_C、VLDL_TG |
| 脂肪酸 | 饱和、单不饱和、多不饱和及具体脂肪酸占比 | SFA、MUFA、PUFA、DHA |
| 氨基酸 | 9 种氨基酸 | Ala、Gln、Gly、His、Ile、Leu、Val、Phe、Tyr |
| 糖酵解 | 葡萄糖、乳酸、丙酮酸、柠檬酸 | Glc、Lac、Pyr、Cit |
| 酮体 | 乙酰乙酸、β-羟基丁酸、丙酮 | AcAce、bOHBut |
| 炎症 | 糖蛋白乙酰化 | GlycA |
这六大类是 Nightingale 自己定义的面板结构,不是 UK Biobank 临时拼的。明白这个结构,你才不会被变量名带偏。比如 HDL 胆固醇和 HDL 颗粒数是两个不同测量维度,LDL 粒径是一个描述符而不是浓度。很多新手把同一子类的所有表示一股脑塞进模型,结果共线性爆炸,回归系数方向都是乱的。
2.3 ukbnmr 内置的三张参考表:biomarker、sample、gp 是处理逻辑的地基
ukbnmr 包里自带参考数据,装好包就能用,不需要联网下载。第一张 biomarker 表记录每个指标的短命名、描述、对应的 UK Biobank 字段编号、类别、单位,以及检测限和空白阈值。第二张 sample 表记录每个样本的 eid、QC 标记和稀释因子、浓度因子等信息。第三张 gp 表把指标归成组,比如按 VLDL、LDL、HDL 分,或者按脂蛋白类、氨基酸类分。
这三张表是包里所有处理逻辑的基准。build 函数读入你的数据包后,本质上做的是三件事:把你的字段编号翻译成短命名、把你的样本 QC 标记整理成统一列、把你的指标按类别归档。理解了这个地基,后面手动改代码时就不容易迷路。
library(ukbnmr) # 查看内置参考表结构 str(ukbnmr$biomarker) head(ukbnmr$biomarker[, 1:6])这段代码用来确认你装的包里 biomarker 表有哪些列。不同版本的列名可能有差异,但至少会包含短命名和字段编号两列。先看一眼,后面写过滤代码时心里有数。
3. 用 ukbnmr 把原始数据包拉成可分析的对象
这一章是干活的核心。我会给你一套从数据包到整洁数据框的最小流程,所有代码都能直接复制运行,但有几处需要根据你装的包版本微调——我会在代码后面说明怎么调。
3.1 最小可用代码:把路径指到数据包目录,一行读入
library(ukbnmr) # 数据包解压后的目录,里面应有 sbin/csv/sample 等配套文件 data_dir <- "D:/UKB/nightingale" ukb_obj <- build_ukbnmr(data_dir)逻辑说明:build_ukbnmr 的输入参数只有一个,就是数据包目录的路径。函数会自动识别目录里的配套文件,把文本或二进制数据读进来,再按内部映射表整理成标准结构。读入时间取决于文件大小,几百 MB 的 .sbin 需要等一会儿,这是正常的。
参数说明:如果你只有文本版数据,没有 .sbin,路径指向包含 csv 的目录即可;如果版本支持 verbose 参数,可以设成 TRUE 查看读取进度,用来判断程序是卡住了还是在慢慢读。
3.2 读入之后先看结构:str、names、head 三连
# 看一眼整个对象的层级结构 str(ukb_obj, max.level = 1) # 列出对象内部的各个部件 names(ukb_obj) # 抽样查看样本信息和生物标志物数据 head(ukb_obj$sample) head(ukb_obj$biomarker)逻辑说明:ukb_obj 是一个列表对象,里面至少包含样本级信息和生物标志物数据两块。str 能看到每个部件的维度,names 给出部件名。因为包版本之间的部件命名略有差异,跑这三行是为了确认你手上版本的实际结构,避免后续写错对象名。
这里有个容易忽略的点:head 看的是前几行,不代表整体顺序。真正做数据合并前,一定要用 eid 列做关联,不能假设两个对象行号一一对应。
3.3 把生物标志物矩阵取出来:长表转宽表的最稳写法
不同版本的 ukbnmr 输出的 biomarker 部件可能是宽表(一行一个样本,一列一个指标),也可能是长表(一行一个样本-指标组合)。我习惯先跑一遍 str 判断,再决定要不要转换。
library(dplyr) library(tidyr) # 如果输出是长表:每行一个 eid + 字段号 + 值 bio_wide <- ukb_obj$biomarker |> select(eid, ukb_field, value) |> pivot_wider(names_from = ukb_field, values_from = value) # 如果输出已经是宽表,直接赋值即可 # bio_wide <- ukb_obj$biomarker逻辑说明:这段代码把长表转成宽表,行是样本,列是字段编号。pivot_wider 的关键参数是 names_from 和 values_from,前者决定新列名来源,后者决定单元格值来源。转换后列名是数字字段号,后面通过映射表翻译成短命名。
参数说明:如果你发现自己的长表里除了 ukb_field 还有 biomarker 短命名列,也可以把 names_from 换成短命名列,一步到位。但我不建议这么做——保留字段编号至少能让你在复盘时回查 SHOWCASE 和审稿材料。
3.4 输出对象里到底有什么:三块内容一次看清
build 完成后,你手里的有效内容可以概括成三块。第一块是样本级信息,包含每个 eid 的 QC 标记和样本处理描述变量;第二块是生物标志物数值,行数是样本数,列数是 Nightingale 面板里的有效指标;第三块是映射与阈值信息,用于把字段编号翻译成短命名、把数值跟检测限做比较。
拿到这三块之后,我的建议是立刻存一份中间结果:把 bio_wide 和样本 QC 表按 eid 合并,输出成一个 RDS 文件。这样下次分析不用重新跑 build,大量节省时间。RDS 的保存用 saveRDS,读回用 readRDS,这是 R 里最稳的中间格式。
3.5 不要绕过工具手工硬读 csv:三个你迟早会踩的坑
有人会问:我自己用 fread 读 csv、再按字段号筛列,岂不是更简单?坦白说,短期的确可以,但你会先后撞上三个坑。
第一个坑是 QC 标记分散在多个字段段里,你不知道哪些字段是 QC、哪些是测量值,手工筛容易漏;第二个坑是低于检测限的数值在交付表里可能已经被替换成普通数字,原始「测不到」的信息被抹掉了,你不知道阈值是多少;第三个坑是 csv、sbin、sample 三个文件的行顺序不一定一致,按行号合并必出错。ukbnmr 把这三个问题做成了固定逻辑,省掉的不只是时间,还有排查错误的痛苦。
4. 样本能不能进分析,取决于 QC:Nightingale QC 标记与过滤逻辑
Nightingale 数据不是「测出来就能用」。仪器自动整合峰面积时,会遇到谱图质量差、样本浓度异常、样本处理出错等情况。UK Biobank 把这些问题用 QC 标记记录了下来,而 ukbnmr 的 sample 表把这些标记统一成了几列。这一章讲清楚怎么读懂、怎么过滤、怎么不滤错。
4.1 三个 QC 标记的分工:核磁 QC、平台 QC、随访 QC
sample 表里通常会有三组 QC 标记:核磁 QC、平台 QC,以及随访数据对应的 QC 标记。核磁 QC 判断原始谱图是否合格,基线漂移、信噪比过低都会被标为失败;平台 QC 是 Nightingale 平台对样本整体质量的综合判断,涵盖浓度异常、峰识别失败等;随访 QC 对应的是随访时采的血液样本,如果你同时分析基线和随访数据,各自的 QC 要分开用。
# 查看各组 QC 标记的分布 table(ukb_obj$sample$qc_nmr) table(ukb_obj$sample$qc_nigh)逻辑说明:table 输出会告诉你每个 QC 标记下有多少样本通过、多少失败。通过率如果低得离谱,先怀疑是不是列名选错了,或者 1/0 方向和你想的相反。
参数说明:1/0 的编码方向在不同数据版本里有差异。不要凭记忆假设 1 一定是通过,先跑 table 看数量级——UK Biobank 的 NMR 数据 QC 通过率通常很高,失败只是少数,如果 table 结果显示失败占大多数,多半是方向理解反了。
4.2 推荐的过滤顺序:先样本级,再变量级
QC 过滤的正确顺序是先样本级、后变量级。样本级 QC 失败意味着这一整条谱都不可信,所有指标都不能用;变量级检查则是看某个指标在所有样本里有多少低于检测限、多少缺失,这决定这个指标是否值得进入后续模型。
library(dplyr) # 第一步:只保留通过核磁 QC 和平台 QC 的样本 valid_eid <- ukb_obj$sample |> filter(qc_nmr == 1, qc_nigh == 1) |> pull(eid) # 第二步:把生物标志物矩阵筛到这些样本上 bio_clean <- bio_wide |> filter(eid %in% valid_eid)逻辑说明:filter 是按行筛选,qc_nmr == 1 和 qc_nigh == 1 两个条件同时满足才保留。pull 把过滤后的 eid 列提取成向量,供下一步使用。bio_clean 行数会比 bio_wide 少,这是预期结果,不是 bug。
参数说明:如果你还要分析随访数据,需要把随访 QC 标记也加进 filter 条件,并且确认 valid_eid 里包含的是随访样本而非基线样本。基线和随访样本不要混在同一个 bio_clean 里做后续的线性模型,原因在 4.4 讲。
4.3 低于检测限不是 0:blanking 阈值与稳健性分析
NMR 定量有一个检测限概念。低于检测限的样本,仪器给出的数值不可靠,甚至可能只是背景噪声,直接当成真实浓度参与均值和回归,会把人体的真实分布拉偏。ukbnmr 的 biomarker 表里保存了每个指标的检测限或空白阈值,这是它比手工读 csv 强很多的地方。
# 查看每个指标的检测限 ukb_obj$biomarker |> select(shortname, ukb_field, limit_of_detection)逻辑说明:把这表打出来,你就能知道每个指标在什么数值以下不可信。拿到阈值后,常见做法有两种:一是把低于检测限的值替换为缺失,二是替换成检测限的一半,作为保守估计。两种策略各跑一遍分析,看结论稳不稳定,这是审稿人认的稳健性分析。
# 示例:把某个指标的低于检测限值替换为缺失 threshold <- 0.2 # 从上一张表里查到的具体阈值 bio_clean <- bio_clean |> mutate(biomarker_x = ifelse(biomarker_x < threshold, NA, biomarker_x))参数说明:ifelse 的三个参数分别对应判断条件、真值、假值。这里的 biomarker_x 要替换成你实际要处理的指标列名。批量处理时可以用 across 加自定义函数,但一次处理一个指标更能看清过滤对样本量的影响。
4.4 批次与随访:另一个必须写进协变量里的因素
Nightingale 数据是分批处理的。UK Biobank 在几十万人的规模上分多个批次完成测量,批与批之间的样本储存时间、试剂批次、仪器校准都可能存在细微差异。ukbnmr 的 sample 表里通常包含批次相关信息,比如板号、孔位、稀释因子等。
我的建议是:不管这些变量最后用不用,先保留在分析数据框里。做关联分析时,在模型里加上批次协变量,或者至少在 QC 过滤后按批次画一遍指标均值和方差,观察有没有系统性偏移。基线和随访数据尤其要小心:随访样本往往在几年后统一检测,和基线的批次差异可能比真实生物学变化还大。把基线和随访混在一起直接比较,很容易得到假的时间趋势。
5. 避坑:我用 ukbnmr 处理 NMR 数据时踩过的五个坑
工具能解决大部分格式问题,但分析设计上的坑还得自己一个个踩。下面五条是我在这类数据处理流程里真实遇到过的案例,现象、原因、解决都列出来,你可以对照排查。
5.1 现象:行数对不上,合并后缺失一大片
把生物标志物矩阵和样本 QC 表合并后,eid 匹配率只有六成;或者直接用 cbind 拼数据,结果变量错位、数值对不上。
原因:sbin、csv、sample 三个文件的行顺序不完全一致,靠行号合并必然错位。
解决:一律以 eid 为主键做合并,绝不按行号。读入后第一件事就是验证两边 ID 集合一致:
all(sort(ukb_obj$sample$eid) == sort(bio_wide$eid))返回 FALSE 就说明读入过程丢了行或存在重复 eid,先排查再继续。这一步能省下后面好几个小时。
5.2 现象:QC 过滤后样本量骤减,跟文档对不上
filter(qc_nmr == 1) 之后样本剩下不到十分之一,与 UK Biobank 官方公布的 QC 通过率完全不符。
原因:想当然认为 1 是通过。某些字段或某些导出版本里,1 和 0 的含义可能和你预期相反。
解决:过滤前先跑 table() 看分布,再对照数据字典确认方向。更稳妥的办法是找一个已知有问题的样本,比如数据字典里注明处理失败的样本,看它在你用的 QC 列里是不是被标成了预期方向。方向确认无误后再写过滤代码。
5.3 现象:变量太多,被描述符指标带跑偏
把所有指标标准化后丢进 PCA,第一个主成分全是粒径和浓度的线性组合;跑 LASSO 选出来的变量也全是同一个脂蛋白子类的不同表示。
原因:Nightingale 面板里,「浓度」「粒径」「占比」不是相互独立的测量,而是一套物理量导出的不同表示。描述符类指标携带大量重复信息,把它们全部放进模型会造成共线性爆炸。
解决:做关联分析前剔除描述符类指标。ukbnmr 提供了剔除描述符的便捷函数,具体名字不同版本有差异,装好后看帮助文档即可;或者自己写代码,在变量名里筛掉 descriptor 关键词。实在不确定时,先看相关性矩阵,同一子类内部相关系数超过 0.99 的,只保留生物学上最可解释的那一个。
5.4 现象:某指标分布直方图在 0 附近有一根异常高柱,关联分析结果还特别「显著」
低丰度指标在 0 附近聚集了大量样本,均值和方差都被拉偏,跑出来的关联结果显著得可疑。
原因:低于检测限的值被当成了真实的低浓度,甚至被当成 0 参与计算。这批「0」其实是测不到,不是真的没有。
解决:用 ukbnmr 的检测限阈值把低于阈值的值标记为缺失,再做一个把低于阈值替换为阈值一半的敏感性分析。两种处理结论不一致时,说明这个指标在该人群里信噪比不足,不适合当主要结局,建议降级为次要探索。
5.5 现象:基线和随访放一起分析,代谢物随年龄「显著升高」,换个年龄段又变成降低
把基线、随访样本放一起算均值,得到的时间趋势在不同年龄段方向不一致,甚至完全相反。
原因:基线和随访数据往往是不同批次、不同储存时间下检测的,批次效应比真实生物学变化大。
解决:基线和随访数据拆开分析,或者把批次变量作为协变量放进模型。ukbnmr 的样本 QC 表对随访样本有专门标记,先按标记拆开,分别做 QC 过滤,再合并做纵向模型。批次效应的处理不是什么高深技巧,但很多人下意识忽略了它。
6. 进阶:用 QC 通过率和检测限给生物标志物排序,快速锁定值得做的指标
掌握了加载和过滤之后,下一步是让 QC 信息服务于研究设计。Nightingale 面板有 250+ 指标,不可能每个都进主分析,怎么快速筛掉不靠谱的指标、给后续建模划定范围,才是进阶的分水岭。
6.1 一个实用的变量筛选流程:非缺失率与低于检测限比例双排序
# 每个指标在 QC 通过样本中的非缺失率 non_miss <- colMeans(!is.na(bio_clean)) # 每个指标低于检测限的比例(阈值从 ukbnmr 的 biomarker 表取) lod_vec <- setNames(ukb_obj$biomarker$limit_of_detection, ukb_obj$biomarker$shortname) below_lod <- sapply(names(non_miss), function(v) { thresh <- lod_vec[v] if (is.na(thresh)) return(NA_real_) mean(bio_clean[[v]] < thresh, na.rm = TRUE) }) # 汇总并排序 tibble( biomarker = names(non_miss), non_missing_rate = non_miss, below_lod_rate = below_lod ) |> arrange(desc(non_missing_rate), below_lod_rate)逻辑说明:non_miss 反映每个指标有多少样本能给出数值,below_lod 反映这些数值里有多少低于检测限。两个维度结合,能把「测不到」和「缺失」区分开。非缺失率高、低于检测限比例低的指标优先进入分析;两者冲突时,优先看低于检测限——非缺失率只能说明有值,不能说明值可信。
参数说明:setNames 的作用是把检测限向量命名为指标短名,方便 sapply 按名字取阈值。如果你的版本里阈值列名不叫 limit_of_detection,先跑 colnames(ukb_obj$biomarker) 确认再替换。
6.2 把映射表迁移到自有数据上
如果以后拿到别的 Nightingale 检测数据,不限于 UK Biobank 来源,ukbnmr 的 biomarker 表依然有用。只要字段编号沿用同一套标准,直接拿映射表翻译列名,你的分析脚本就能整套复用。我的习惯是把这份映射表单独存一份 CSV 放在项目目录里,所有脚本统一引用,不每次现查 SHOWCASE。这样不仅代码干净,审稿时需要回查字段定义也方便。
6.3 一个教训
最早我图快,把 Nightingale 数据用通用工具一次性读成宽表,跳过 QC 直接跑相关性,结果冒出一堆假信号,花了整整两周排查才发现是低于检测限的样本在作怪。后来老老实实按「映射表 → 样本 QC → 检测限」三步走,半小时就能产出一张可靠的分析矩阵。现在不管多急,这三件事我都会先做完再谈建模。希望帮到你。
本文还有配套的精品资源,点击获取