R语言数据科学在现场流行病学调查分析中的应用
以An unusal outbreak of food poisoning, Thailand, 1990为例
ZSCDC
2021-01-13
1990年8月25日在泰国某大学召开国家残障运动会后20小时，485人餐后2-3小时出现了上消化道症状。当时大约1300人参会，包括老师、学生、教练和残障人氏。次日，流行病学专家组进驻开展调查。
1 材料与方法
1.1 数据探索性分析
1.2 描述性流行病学研究
三间分布
环境调查
1.3 分析流行病学研究
1.4 病例定义
病例：居住在事发大学内且食用1990年8月25日提供晚餐，出现恶心、呕吐、腹痛、腹泻等症状之一者；
非病例：居住在事发大学内且且食用1990年8月25日提供晚餐，未出现上述症状者。
1.5 实验室检测
2 结果
2.1 探索性分析
2.1.1 数据概览
## Rows: 1,094
## Columns: 13
## $ id <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17...
## $ sex <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ age <dbl> 13, 14, 13, 15, 14, 11, 19, 17, 15, 20, 18, 17, 21, 16, 1...
## $ exptime <I<chr>> 25330825180000, 25330825180000, 25330825180000, 253308...
## $ beefcurry <dbl> 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ saltegg <dbl> 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ eclair <dbl> 1.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.0, 0.0, 1....
## $ water <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ onset <I<chr>> 25330825220000, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA...
## $ nausea <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, ...
## $ vomiting <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, ...
## $ abdpain <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, ...
## $ diarrhea <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, ...## id sex age exptime
## Min. : 1.0 Min. :0.000 Min. : 1.00 Length:1094
## 1st Qu.: 274.2 1st Qu.:0.000 1st Qu.:14.00 Class :AsIs
## Median : 547.5 Median :1.000 Median :18.00 Mode :character
## Mean : 547.5 Mean :0.659 Mean :23.69
## 3rd Qu.: 820.8 3rd Qu.:1.000 3rd Qu.:24.00
## Max. :1094.0 Max. :1.000 Max. :99.00
## beefcurry saltegg eclair water
## Min. :0.0000 Min. :0.000 Min. : 0.00 Min. :0.000
## 1st Qu.:1.0000 1st Qu.:1.000 1st Qu.: 0.00 1st Qu.:1.000
## Median :1.0000 Median :1.000 Median : 2.00 Median :1.000
## Mean :0.9534 Mean :0.957 Mean :11.48 Mean :1.021
## 3rd Qu.:1.0000 3rd Qu.:1.000 3rd Qu.: 2.00 3rd Qu.:1.000
## Max. :9.0000 Max. :9.000 Max. :90.00 Max. :9.000
## onset nausea vomiting abdpain
## Length:1094 Min. :0.0000 Min. :0.0000 Min. :0.0000
## Class :AsIs 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000
## Mode :character Median :0.0000 Median :0.0000 Median :0.0000
## Mean :0.3985 Mean :0.3793 Mean :0.3501
## 3rd Qu.:1.0000 3rd Qu.:1.0000 3rd Qu.:1.0000
## Max. :1.0000 Max. :1.0000 Max. :1.0000
## diarrhea
## Min. :0.0000
## 1st Qu.:0.0000
## Median :0.0000
## Mean :0.2148
## 3rd Qu.:0.0000
## Max. :1.0000## [1] 668| id | sex | age | exptime | beefcurry | saltegg | eclair | water | onset | nausea | vomiting | abdpain | diarrhea |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 0 | 39 | 0 | 0 | 0 | 0 | 629 | 0 | 0 | 0 | 0 |
经概览，结合数据集说明，发现此数据集以下几个问题：
本数据集缺省值有668个，其中暴露时间缺失39个，发病时间缺失629个，其余变量均有赋值；
性别为双精度数值型，取值为0，1，65.9%为男性；
年龄为双精度数值型，99或为缺失值代码；
暴露时间、发病时间均为字符型，时间格式为佛历；
beefcurry-water等4种食物（饮品）均为双精度数值型，9或为缺失值代码；其中eclair中取值较多，90为缺失值，80为食用但不记不清数量；
症状均为双精度数值型，取值为0，1.
必要时，可制作数据集的探索性分析报告，方法如下：
2.1.2 源数据可视化探索
2.1.2.1 变量相关矩阵可视化
上图为经规整后数据集中各变量间的相关关系图，可以重点关注一下是否发病case变量与其他变量之间的相关关系，不难发现是否发病与4种食物暴露因素间仅有是否食用eclair的eclair1变量相关关系较大（\(r = 0.48\)）。后续分析中值得关注。
2.1.2.2 变量分布点阵图
将发病病例筛选出，然后以病例ID为横坐标，各疑似暴露因素为纵坐标，以各变量值区分散点颜色，绘制病例的各暴露因素间的散点图。病例较多时，可以随机抽样100例以内的病例作展示。一般情况下，可以发现发病点与暴露源因素间的伴随关系。本起疫情中由于供应的菜品较少，只有4种，目前还看不出明确的伴随关系。也就是说各种暴露因素暴露率均高时，此图不易发现征兆，但在暴露众多的情况下，此图极有价值。
以上两图极有价值，多数案例经过上述两图展示，基本可以锁定暴露源因素。
2.2 描述性流行病学研究
2.2.1 时间分布
约1300人参加运动会，调查组调查了1210人，其中113人未食用晚餐，3人发病但未调查到信息。 最终取得了1094人的调查信息。其中病例469人，非病例625人，罹患率为42.87%。首例病例出现在1990-08-25 15:00:00，末例病例出现在1990-08-26 09:00:00，约87.63%的病例出现在1990年8月25日20时至26日零时这4个小时内。本起疫情最短潜伏期为1小时，最长潜伏期为14.5小时，平均潜伏期为3.5小时。
疫情发病曲线如下图。 
暴露与发病时序关系图如下：
潜伏期分布如下图：
该起疫情发病曲线呈现发病单峰形态提示应考虑本次疫情由一次性同源暴露产生，疑似暴露餐次与发病时序图可以提供某种印证。
2.2.2 人群分布
469病例中男性332例，男性罹患率为46.05%(332/721)；女性137例，女性罹患率为36.73%(137/373)。
最高罹患率出现在< 10岁组，罹患率为70%(7/10)。年龄组罹患率情况见下图。
年龄组罹患率情况与临床症状分布如下图。
2.2.3 空间分布
因数据集缺乏相关数据，略。
2.3 分析性研究
2.3.1 食物特异罹患率
rowname | df.Tot.Exp. | df.Exp.Cases | df.AR. | df.Tot.Unex. | df.Unex.Cases | df.AR..1 | df.RR | df.CI.ll | df.CI.ul | df.p.Fisher. |
eclair1 | 683 | 383 | 56 | 294 | 15 | 5.1 | 11.0 | 6.68 | 18.1 | 0.00 |
beefcurry | 998 | 447 | 45 | 91 | 22 | 24.2 | 1.8 | 1.28 | 2.7 | 0.00 |
saltegg | 1,002 | 448 | 45 | 87 | 21 | 24.1 | 1.8 | 1.27 | 2.7 | 0.00 |
water | 1,063 | 458 | 43 | 25 | 10 | 40.0 | 1.1 | 0.66 | 1.8 | 0.84 |
由上表可见，食用eclair者的罹患率为56.08%(383/683)，其患病风险是未食用者的10.99(95%CI 6.68-18.07)倍。同样的食用beefcurry和seategg者的患病风险是分别是未食用者的1.85(95%CI 1.28-2.68)倍和1.85(95%CI 1.27-2.71)倍。饮水者与未饮水者的患病风险差异无统计学意义。
2.3.2 调整相对危险度
items | RR | lowerCI | upperCI | p.value |
(Intercept) | 0.058 | 0.0092 | 0.37 | 2.6e-03 |
sex1 | 1.949 | 1.4080 | 2.70 | 5.7e-05 |
agegroup10 - 29 | 0.653 | 0.1111 | 3.83 | 6.4e-01 |
agegroup> 30 | 0.334 | 0.0547 | 2.04 | 2.3e-01 |
beefcurry | 1.006 | 0.4857 | 2.08 | 9.9e-01 |
saltegg | 0.761 | 0.3671 | 1.58 | 4.6e-01 |
water | 1.302 | 0.4305 | 3.94 | 6.4e-01 |
eclair1 | 23.978 | 13.7768 | 41.73 | 2.7e-29 |
| Characteristic | OR1 | 95% CI1 | p-value |
|---|---|---|---|
| sex | |||
| 0 | — | — | |
| 1 | 1.95 | 1.41, 2.70 | <0.001 |
| agegroup | |||
| < 10 | — | — | |
| 10 - 29 | 0.65 | 0.09, 3.67 | 0.6 |
| > 30 | 0.33 | 0.05, 1.95 | 0.2 |
| beefcurry | 1.01 | 0.49, 2.11 | >0.9 |
| saltegg | 0.76 | 0.36, 1.58 | 0.5 |
| water | 1.30 | 0.43, 4.08 | 0.6 |
| eclair1 | 24.0 | 14.2, 43.4 | <0.001 |
|
1
OR = Odds Ratio, CI = Confidence Interval
|
|||
单暴露因素分析中发现3种食物食用后发病风险显著增加。为了控制可能混杂因素，采用非条件logistic回归法进行分析发现：食用eclai者的发病风险经调整后增大至23.98(95%CI 13.78-41.73)。
2.3.3 剂量反应关系
剂量 | 发病 | 未发病 | RR | CI下限 | CI上限 | p值 |
0- | 15 | 279 | 1.0 | |||
0.5-1 | 51 | 54 | 9.5 | 5.6 | 16 | 0 |
2-4 | 298 | 226 | 11.2 | 6.8 | 18 | 0 |
4- | 34 | 20 | 12.3 | 7.2 | 21 | 0 |
为进一步评估因果关系，计算食用eclairs的剂量反应关系可见，食用越多eclairs，发病风险越高(\(p< 0.01\))。
3 调查结论
从流行病学调查分析来看，本起疫情是参加1990年8月24-26日全国残障运动会的人员食用了25日大会供应的晚餐中被污染的eclair后引发的同源一次暴露消化道聚集性疫情。
4 参考文献
Thaikruea, L. , Pataraarechachai, J. , Savanpunyalert, P. , & Naluponjiragul, U. . (1995). An unusual outbreak of food poisoning. Southeast Asian J Trop Med Public Health, 26(1), 78-85.
Data source: package{epiDisplay}