企业官网建设流程全解析

第一次接手FAERS数据项目时，我以为最难的是R语言的代码，真正动手之后才发现，最耗时间的不是写模型，而是摸清这个数据库的脾气。几百万行原始报告、七个文件拆开存放、同一份报告可能存在多条重复记录、药品名称五花八门——如果前期没把这些理顺，后面不管你用ROR还是PRR算信号，结果都站不住脚。

FAERS全称FDA Adverse Event Reporting System，也就是美国FDA的不良事件报告系统，收录的是药品上市后自发提交的不良事件报告。对于我们做药物警戒、临床药理、药物流行病学的人来说，它是一个绕不开的真实世界数据源。而R语言的强项正好覆盖了FAERS分析的全流程：批量下载、百万行级清洗、列联表统计、信号挖掘、可视化输出。这篇文章就当一份完整实操记录来写，适合刚接触FAERS、想用R语言跑药物不良事件信号挖掘的读者。数据规模大、字段杂、历史版本多，但你只要把下面这几关走通，后面所有分析都能自己掌控。

1. FAERS到底存了什么：七个季度文件和一个去重规则

FAERS数据按季度发布，每个季度一个压缩包，解压后是若干个以$作为字段分隔符的文本文件。这个数据库最让人头大的地方在于：它不是一张大宽表，而是拆成多个逻辑模块，靠主键串起来。

1.1 七个核心文件分别记录了什么

我按惯例处理FAERS时，默认会先认识这几个文件。虽然历史上AERS（FAERS的前身）和不同年份的文件字段有出入，但近几年的季度发布基本稳定为下面这套结构：

• DEMO：报告的人口学信息，包括PRIMARYID、CASEID、报告日期FDA_DT、性别、年龄、报告类型、报告国家等。它是整个分析的主表。
• DRUG：报告涉及的药品信息，包括药品名称、给药途径、剂量、疗程，以及最重要的角色代码role_cod。
• REAC：不良事件结果表，核心字段是PT（MedDRA首选术语），有时候还带SOC（系统器官分类）。
• OUTC：患者结局，比如死亡、住院、残疾等，一个报告可以对应多个结局。
• RPSR：报告来源，比如是医生报告、药师报告、消费者报告还是申办方报告。
• THER：药物治疗方案，包括用药起止时间。做信号挖掘时用得少，但做药物相互作用或疗程分析时会用到。
• INDI：适应症，也就是用药是针对什么疾病开的。

每个文件之间靠什么关联？主要是PRIMARYID。DEMO表里PRIMARYID是报告主键，DRUG、REAC、OUTC这些表里每条记录都带有对应的PRIMARYID。注意，DRUG表和REAC表不是一对一的，一份报告可能有多条药品记录、多条不良事件记录，所以要避免把它们直接粗暴横向合并造成笛卡尔积。

1.2 文件格式里那些容易翻车的细节

FAERS发布的数据文件后缀大多是.txt，但你如果用默认参数去读，会死得很难看。字段分隔符不是制表符也不是逗号，而是$。我第一次读的时候，打开文件看了眼觉得像文本，直接read.table(file, header=TRUE)，结果整列数据全串位，连表头都对不上，排查了半天才发现分隔符的问题。

另一个坑是字段内容中可能包含引号、换行符或空字段，尤其REAC文件里$分隔字段时，有的字段是空白。稳妥做法是：

raw <- read.delim("REAC20Q1.txt", header = TRUE, sep = "$", quote = "", fill = TRUE, stringsAsFactors = FALSE)

对大文件来说，read.delim其实不够快，几百万行的REAC文件会读到人发困。我更推荐用data.table::fread：

library(data.table) reac <- fread("REAC20Q1.txt", sep = "$", quote = "", fill = TRUE, header = TRUE)

fread读取速度比基础函数快一个数量级，而且它会自动识别列类型，处理百万行没有任何压力。还有一点需要注意：不同年份发布的FAERS文件，字段顺序和列名可能做过调整，尤其是DEMO文件。比较稳妥的做法是读取后先遍历一下列名，再根据自己需要的字段进行选择，而不是硬编码第几列。

1.3 为什么信号挖掘一般要合并多个季度数据

FAERS单季度报告量看起来很多，但分摊到具体某一个药物和某一种不良事件上，数字往往很小。你要是只用单季度数据算ROR，可能目标事件记录只有2条，统计功效完全不够，算出来的置信区间宽得没意义。

国内外发表在期刊上的FAERS信号挖掘研究，普遍采用连续几年的数据，常见的是合并8个季度（2年）甚至更长窗口。拉长时间段能缓解自发报告漏报、报告滞后的问题，也能让罕见不良事件有足够的观察频次。合并的思路很简单，把多个季度的CSV或TXT文件读进来，用rbindlist堆叠，或者先按季度读成list再合并。

library(data.table) files <- list.files("faers_raw", pattern = "REAC.*\\.txt$", full.names = TRUE) reac_all <- rbindlist(lapply(files, fread, sep = "$", quote = "", fill = TRUE))

不过合并不是简单的叠加，紧接着就是下一个关键步骤：去重。

2. 从官网到R环境：数据下载与导入的完整链路

标准的数据获取路径是去FDA官网的FAERS Quarterly Data Extract Files页面，手动选择季度下载压缩包。但如果要追溯多年数据，手动点击特别累。你可以用R脚本半自动化完成，也可以借助openFDA API快速取数。

2.1 三种取数方式的适用场景

我常用的取数方式有三种，按效率和适用场景区别很大：

openFDA的drug/event.json接口适合先快速看某个药的报告量，比如你还不确定某药物是否值得做全量信号挖掘，可以先发一个查询确认有没有数据基础。但真做全量挖掘，还是建议用季度文件，自己掌握全部细节，灵活度最高。

2.2 用R脚本批量下载并解压季度文件

官网季度压缩包下载链接有规律可循，适合写循环批量下载。下面这段代码是我自己的常用写法，把需要下载的季度拼接好，下载到本地并解压：

library(data.table) quarters <- c("20q1", "20q2", "20q3", "20q4", "21q1", "21q2", "21q3", "21q4") for (q in quarters) { url <- paste0("https://example.fda.fis.faers/", q, ".zip") # 替换为实际页面链接 dest <- paste0("faers_raw/", q, ".zip") download.file(url, dest, mode = "wb", quiet = TRUE) unzip(dest, exdir = paste0("faers_raw/", q)) }

下载之前建议先用curl或浏览器确认链接模式。压缩包内文件名一般带有季度标识，比如DEMO20Q1.txt、DRUG20Q1.txt。解压后按文件前缀分别读取即可。这个过程看着简单，但有几个细节值得注意：下载模式一定要设mode = "wb"，否则Windows下压缩包容易损坏；官网有时候会调整链接，定期访问页面确认没有失效比写一套“永久脚本”更实际。

2.3 读取文件的兼容处理

不同季度的FAERS文件列名不完全一致，所以读取之后要做的第一件事不是直接分析，而是统一列名。以DRUG文件为例，历史上有的季度叫DRUGNAME，有的可能略有差异，需要先检查：

drug_files <- list.files("faers_raw", pattern = "DRUG.*\\.txt$", full.names = TRUE) drug_list <- lapply(drug_files, fread, sep = "$", quote = "", fill = TRUE) # 查看所有文件的列名，确认差异 lapply(drug_list, names) # 统一为小写，方便后续处理 drug_list <- lapply(drug_list, function(dt) { setnames(dt, names(dt), tolower(names(dt))) dt }) drug_all <- rbindlist(drug_list, fill = TRUE)

rbindlist的fill=TRUE参数很关键，它允许不同季度的文件字段不完全对齐，多出来的列保留为NA，缺失的列也能正确堆叠起来。这一步能省掉大量手工对齐工作。

3. 清洗阶段真正决定分析质量：角色筛选、去重与药品名归一化

数据清洗是FAERS分析中最不性感却最关键的部分。我在项目里反复跟人强调：信号检测的代码人人都能写，但数据清洗质量直接决定你算出来的信号是真信号还是垃圾。

3.1 同一份报告为什么会出现多条记录

FAERS的重复问题非常突出。原因在于，同一份病例可能被不同途径重复上报，比如医院报告了一份，药企又因为同一病例递交了一份；或者FDA在后续跟进中更新了信息，生成了新的记录。如果不去重，这些重复记录会放大目标药物-事件组合的报告数，最终高估信号强度。

FDA官方给出了明确的去重规则：以PRIMARYID为键，如果同一PRIMARYID对应多条记录且FDA_DT（FDA收到日期）不同，保留FDA_DT最新的一条；如果FDA_DT相同，则保留CASEID最大的一条。翻译成R代码就是：

demo_all <- demo_all[order(FDA_DT, CASEID, decreasing = TRUE)] demo_all <- demo_all[!duplicated(demo_all$PRIMARYID), ]

注意排序的方向一定不能搞反。有个同事当时图省事，按CASEID降序排列去重，结果FDA_DT老的数据反而被留下了，整个分析全部重跑。这个去重结果还要同步过滤DRUG和REAC表，只保留去重后的PRIMARYID集合，保证后续关联是在同一份报告基础上进行。

3.2 DRUG角色代码到底该不该全保留

FAERS的DRUG文件里有个role_cod字段，常见取值有：

• PS：Primary Suspect，首要怀疑药
• SS：Secondary Suspect，次要怀疑药
• C：Concomitant，伴随用药
• I：Interacting，相互作用药

常规信号挖掘通常只保留PS和SS，理由很直观：不良事件之所以发生，最有解释力的嫌疑药是PS，SS次之。如果一股脑把所有伴随药都算进目标药物，噪声会非常大。比如患者同时吃着七八种药，真正引起肝损伤的可能只有一种，但其他伴随药也会出现在该不良事件的报告里，计算ROR时全部拉进来，会严重稀释真实信号。

我个人的经验是在做正式分析时默认只取PS和SS，但建议把包含C和I的全量数据也保留一份。如果后续发现某个信号的ROR很异常，你可以回去看一下加入伴随药之后信号是否消失，这能帮判断是否混杂了适应症因素或联合用药因素，这种敏感性分析放在论文里也很加分。

drug_ps <- drug_all[role_cod %in% c("PS", "SS")]

3.3 药品名称归一化：别名问题绕不开

这是FAERS清洗里最让人头疼的一步。DRUGNAME是自由文本，同一个活性成分在不同报告里可能有几十种写法。比如某个经典药物，报告里可能同时出现商品名、通用名、带剂型的写法、带盐基的写法、甚至大小写混搭。

我能给出的实用建议是分两步走。第一步是文本粗清洗，统一大写、去掉空格和常见剂型后缀；第二步是构建一个映射字典，把常见的商品名和别名映射到标准活性成分名。粗清洗代码长这样：

library(stringr) drug_ps[, drug_clean := toupper(drugname)] drug_ps[, drug_clean := str_squish(drug_clean)] drug_ps[, drug_clean := str_remove_all(drug_clean, "\\b(TABLET|CAPSULE|INJECTION|SOLUTION|CREAM)\\b")] drug_ps[, drug_clean := str_remove_all(drug_clean, "[0-9]+\\s*(MG|ML|MCG|G|IU)")]

这样清洗完，还能剩不少无法归一的别名，那就必须建人工映射表。有条件的单位会去买WHO Drug词典做标准化映射，但对于个人研究或预算有限的团队，自己维护一个高频名称字典足够了。原则是：宁可保留少量无法归一的名称，也不要为了合并而把两个不同活性成分归成同一个。

4. 信号检测的R实现：ROR、PRR与贝叶斯方法的取舍

清洗结束后就可以进入正题了：计算药物-不良事件组合的失衡信号。这部分说白了就是把每对“药物-不良事件”放到一个2×2列联表里，看目标组合的观测频数是否显著高于背景频率。

4.1 2×2列联表的构建逻辑

以目标药物drug_i和目标不良事件event_j为例，列联表的四个格子分别是：

• a：同时报告了drug_i和event_j的报告数
• b：报告了drug_i但没报告event_j的报告数
• c：没报告drug_i但报告了event_j的报告数
• d：既没报告drug_i也没报告event_j的报告数

很多教程会直接把这个表的构建跳过，但实际上，对全量药物-事件对批量构建列联表，是计算性能上最大的瓶颈。高效的实现方法是用data.table做条件计数。先把药品表按PRIMARYID展开成每个报告是否使用目标药物，事件表也展开成每个报告是否发生目标事件，然后做merge，再对每个药物-事件组合分组计数。

# 假设 drug_long 是每个报告对应的药物 # 假设 reac_long 是每个报告对应的不良事件 report_drug <- drug_ps[, .(PRIMARYID, drug_clean)] %>% unique() report_reac <- reac_all[, .(PRIMARYID, pt)] %>% unique() # 全体报告的药物-事件对组合 combo <- merge(report_drug, report_reac, by = "PRIMARYID", allow.cartesian = TRUE) combo_counts <- combo[, .N, by = .(drug_clean, pt)] # 每个报告的数量，用于计算背景 n_report <- demo_dedup[, .N]

真正细化的列联表计算要分别统计含目标药物的报告总数、含目标事件的报告总数、同时含两者的报告数，再反推b、c、d。代码可以封装成函数，对每一个目标药物跑一遍。如果你分析的药物数量不多，循环足够用；如果是全量筛选，建议用矩阵方式预先计算所有药物和事件的边缘频数。

4.2 ROR和PRR的计算公式与95%置信区间

ROR，Reporting Odds Ratio，报告比值比，是FAERS信号挖掘最主流的方法之一。它的原理与此前分析“暴露-疾病关联”的思路基本一致，只是这里的“病例”变成了“目标不良事件的报告”，“对照”变成了“目标药物不相关的不良事件报告”。

ROR的公式：

ROR = (a / b) / (c / d) = (a × d) / (b × c)

95%置信区间则取对数后再变换回来：

SE(ln(ROR)) = sqrt(1/a + 1/b + 1/c + 1/d) 95% CI = exp(ln(ROR) ± 1.96 × SE(ln(ROR)))

判断信号通常用下限大于1，同时要求报告数a不少于3。有的顶刊文章会要求更严格，比如a≥5，这个看你投稿期刊的习惯。

PRR，Proportional Reporting Ratio，比例报告比值，是另一种经典频率派方法。它的出发点是：在目标药物所有不良事件报告中，目标事件所占的比例，与所有其他药物报告里目标事件所占比例之间的比值：

PRR = (a / (a + b)) / (c / (c + d)) = (a × (c + d)) / (c × (a + b))

PRR的置信区间公式里有个坑，它和ROR不一样，加减号里前面是加，后面是相减：

SE(ln(PRR)) = sqrt(1/a - 1/(a+b) + 1/c - 1/(c+d)) 95% CI = exp(ln(PRR) ± 1.96 × SE(ln(PRR)))

实践中行业常采用的PRR阳性判断标准是三个条件同时满足：PRR≥2、卡方≥4、a≥3。这个标准来自Evans等人在2001年发表的经典文章，至今仍被许多药物警戒文章引用。

signal_table <- combo_counts[, { a <- N b <- drug_total - N # drug_total 为该药物在报告集中出现的报告数 c <- event_total - N # event_total 为该事件在报告集中出现的报告数 d <- n_report - a - b - c ror <- (a * d) / (b * c) ror_ci_low <- exp(log(ror) - 1.96 * sqrt(1/a + 1/b + 1/c + 1/d)) ror_ci_up <- exp(log(ror) + 1.96 * sqrt(1/a + 1/b + 1/c + 1/d)) prr <- (a / (a + b)) / (c / (c + d)) prr_ci_low <- exp(log(prr) - 1.96 * sqrt(1/a - 1/(a+b) + 1/c - 1/(c+d))) prr_ci_up <- exp(log(prr) + 1.96 * sqrt(1/a - 1/(a+b) + 1/c - 1/(c+d))) list(ror = ror, ror_ci_low = ror_ci_low, ror_ci_up = ror_ci_up, prr = prr, prr_ci_low = prr_ci_low, prr_ci_up = prr_ci_up, a = a) }, by = .(drug_clean, pt)]

这个代码看起来短，跑起来可能会出现很多极端值。比如b为0时ROR是无穷大，或者分母太小导致置信区间爆炸。我的建议是过滤掉b或c为0的组合，或者对ROR设一个合理的截断值，避免画图时个别点把整个坐标轴都拉伸变形。

4.3 BCPNN和EBGM：贝叶斯方法的应用场景

频率派方法ROR和PRR计算简单、解释直观，但在报告数极少时极不稳定。你在小样本场景下会看到ROR被算成一个巨大的数，置信区间从0.2跨到几百，这种结果基本没有应用价值。这时候贝叶斯收缩方法更靠谱。

BCPNN，贝叶斯置信传播神经网络，是WHO乌普萨拉监测中心常用的方法，它计算信息成分IC（Information Component）。IC的本质是量化目标药物-事件组合的共现概率相对独立假设的偏离程度。R中有PhViD包可以直接实现，核心思想是通过贝叶斯后验分布对罕见事件的估计值进行收缩，减少假阳性。

另一个常见方法是EBGM（Empirical Bayesian Geometric Mean），FDA自己的信号检测系统用的就是MGPS（Multi-item Gamma Poisson Shrinker）框架，EBGM就是它输出的核心指标。R的openEBGM包提供了一套从数据预处理、EM算法估计超参数到输出EBGM及置信区间的完整流程。

我个人的取舍习惯是：初步筛选用ROR和PRR这种轻量方法，快速得到候选信号列表；对候选里比较重要或者报告数偏少的组合，再用BCPNN或EBGM做二次确认。三板斧全过，才值得进入后续人工评估。

4.4 阈值不能生搬硬套

每次我帮人看分析结果，都会问一句：你的阈值是哪里来的？ROR下限大于1是一个常用标准，但不代表它适用于所有分析场景。不同药物类别、不同适应症人群的信号分布完全不一样，用一种固定阈值筛出来的结果，可能遗漏真实信号，也可能堆满假阳性。

我更推荐的做法是多方法交叉验证。比如把同时满足ROR下限>1和PRR≥2、a≥3的组合作为候选信号，再对候选信号做BCPNN的IC置信区间检查。这样虽然会扩大初步计算量，但至少在“宁可漏报”和“宁可错报”之间做了平衡。最终落到论文里的阳性标准，建议参考同领域已发表文献的做法，不要自己临时发明一个组合阈值，否则审稿人问起依据，你答不上来。

另外一个重要的提醒：FAERS是自发报告数据库，存在报告偏倚、媒介影响、适应症混杂、上市时间效应等问题。信号再强，也只是“信号”，不是“因果关系”。我见过有人把PRR=20的结果直接当成药物导致不良反应的实锤，这是对药物警戒基本逻辑的误解。信号挖掘的作用是生成假设，后续还要做因果关系评估，参考Bradford Hill准则或WHO-UMC因果关系分级，结合临床前数据、临床试验证据和其他真实世界数据库复核，才能下结论。

5. 从数字到结论：可视化输出与信号解读的边界

算完信号，下一步就是把结果变成能让人一眼看懂的东西。FAERS信号挖掘中最常画的是火山图、气泡图和森林图，R语言在这一步的优势非常明显。

5.1 火山图快速定位强信号

火山图的横轴放效应值，纵轴放显著性指标，是信号筛选的标准视图。用ggplot2可以轻松实现：

library(ggplot2) library(ggrepel) signal_table[, log_ror := log2(ror)] signal_table[, neg_log_p := -log10(p_value)] # p_value 由卡方或Fisher精确检验得到 ggplot(signal_table, aes(x = log_ror, y = neg_log_p)) + geom_point(aes(size = log10(a)), alpha = 0.6, color = "#2c6fbb") + geom_vline(xintercept = 0, linetype = "dashed") + geom_hline(yintercept = -log10(0.05), linetype = "dashed") + labs(x = "log2(ROR)", y = "-log10(P)") + theme_minimal(base_size = 14)

点越靠右上，说明目标药物-目标事件组合的信号强度和统计显著性越突出。点的尺寸映射报告数之后，也能快速分辨出一个信号是建立在几十条报告上还是几千条报告上。

5.2 气泡图和森林图各有用武之地

气泡图很适合做多个药物和多个事件的横向对比。横轴放log2(PRR)，纵轴放log2(ROR)，点大小映射报告数，颜色映射信号等级。如果某种药物在两种方法下都显示强信号，点会集中在右上角，一眼就能识别。

森林图更适合展示单个目标药物下多个不良事件的风险谱。比如你要看某降糖药到底和哪些不良事件关联最强，把每个PT的ROR点估计和置信区间横着画出来，能非常清楚地呈现风险排序：

target_drug <- signal_table[drug_clean == "SITAGLIPTIN" & a >= 3] target_drug[order(ror), pt := factor(pt, levels = pt)] ggplot(target_drug, aes(x = ror, y = pt)) + geom_point(size = 2.5, color = "#c0392b") + geom_errorbarh(aes(xmin = ror_ci_low, xmax = ror_ci_up), height = 0.3) + geom_vline(xintercept = 1, linetype = "dashed") + labs(x = "ROR (95% CI)", y = NULL) + theme_minimal(base_size = 13)

画这种图的时候，过滤条件一定写清楚。我是建议至少保留a≥3且置信区间不跨越1的事件，不然森林图里塞几十个无意义事件，整张图不仅看不清楚，还会让人觉得你完全没有做信号筛选的章法。

5.3 信号解读的边界：一次分析能说明什么

做FAERS分析久了，我越来越觉得，中间的计算过程反而不是最难的，最难的是知道结果的“边界”。曾经有一次我分析某个DPP-4抑制剂品种，膀胱癌的ROR显著升高，朋友看到结果很兴奋，觉得这是重大发现。但仔细看报告源，相当一部分报告集中在该品种上市后媒体高度关注的几年，患者和医生在媒体报道期间对这类事件的上报意愿明显提高，这就是典型的报告偏倚。

自发报告数据库的固有缺陷决定了FAERS结果必须谨慎解读。疗程长、媒体关注度高、适应症人群本身风险较高的药物，天然会产生更多信号。所以我在做结论时，一定会交叉核对三件事：该事件在同类药物中是否也出现、目标药物是否有相关的临床前证据或临床试验信号、报告的临床描述是否提供了合理的时间关联和机制解释。只看ROR数字就下结论，建议直接放弃。

6. 一点个人经验：先建缓存数据，再反复喂给模型

最后分享一个实际项目里帮我节省了大量时间的小习惯：处理完FAERS原始文件后，不要反复重读那些几百MB的TXT文件，而是把清洗后的长表保存成RDS或parquet格式，后续所有分析都从这个缓存文件读取。

saveRDS(drug_clean_final, "faers_clean/drug_clean.rds") saveRDS(reac_clean_final, "faers_clean/reac_clean.rds")

下次写脚本时直接readRDS，读取时间从几十秒缩短到一两秒。分析又是一个不断迭代的过程，今天你筛了300个信号，明天可能想换阈值再看一遍，每次都从原始TXT重新解析，时间成本完全不可接受。我自己最初就是没有这个意识，每次调参都要等十几分钟读数据，喝了两杯咖啡才反应过来该做个缓存。

FAERS整条链路跑通之后，你会发现这个库其实没有想的那么高不可攀。它的门槛不在于统计学或R语言，而在于你是否愿意花时间把数据结构、清洗规则、方法取舍这些基本功做扎实。希望这篇实操记录能帮你在做药物警戒数据分析时少踩几个坑，把时间花在真正有价值的问题上。