![]()
系列简介
这是我们一系列原创技术贴,从易到难,每天学习一点。所有内容均为疾控数据分析、科研论文相关,或者说很多和现在的热门监测预警相关,所以我们这个系列就叫“监测预警基础”。
之前很多老师反馈这种小专题的形式特别好,可以一口气彻底学透一种方法,好的,我们又来了!
这次是分布滞后非线性模型,我们初步计划用8小节左右的内容讲通它的前世今生,以及如何使用。所以我们这个小专题的名字就叫做《八小篇搞懂DLNM》
今天是第7小篇——R语言实操以后结果的解读,好的,可能也是我们最后一篇,当然,如果后面还有时间,我们还会补一篇,以一个实际例子走完DLNM整个过程。
上一篇我们带大家从零跑通了第一个 DLNM 模型,很多同行跟着复现后,对着输出结果提出了一系列非常典型的疑问:为什么交叉基是 12 个变量?时间样条为什么有 21 个系数?为什么只拿 28℃和 20℃比?数据里的相对湿度怎么没放进模型?什么时候不用泊松而要换成负二项回归?
今天这篇,我们就对着真实输出结果,把跑完模型后最常见的疑问一次性讲透,带你从「能跑通代码」进阶到「看得懂原理、说得清结论」。
![]()
很多新手打开summary(model)的第一反应是头晕:几十行系数,到底看哪个?其实不用逐行死磕,把输出分成 4 个区域,按优先级看就够了。
1. 模型调用区(Call):先校验有没有写错
作用是核对模型配置,避免低级错误看一眼公式里的变量、分布类型、连接函数是不是和你预期的一致。比如有没有漏加星期效应、是不是选错了分布,这一步就能发现。你的结果里,交叉基、时间样条、星期效应都在,分布是泊松对数连接,完全符合设定。
![]()
2. 系数表区(Coefficients):分三类区别对待
这是内容最多的区域,但不是每个系数都有业务意义,按变量类型分三类:
第一是交叉基系数(cb 开头):单个系数无直接解读价值,它们是基函数的组合系数,必须通过crosspred()整合后才能算出相对危险度 RR。
第二是混杂控制系数(时间样条、星期效应):可以看整体显著性,验证混杂控制是否有效;其中星期效应的系数可以直接解读,验证是否符合业务常识。
第三截距项:仅数学意义,无业务解读价值。
![]()
以上结果如果我们看一下的话就是
交叉基系数(cbv1.l1cbv2.l3这类 cb 开头的)不用单独解读,单个系数无直接流行病学意义。交叉基是「暴露维度 4 个样条基 × 滞后维度 3 个样条基」= 12 个基变量,所以你会看到 12 个cb开头的系数;它们是组合生效的,必须通过crosspred()函数整合起来,才能算出具体温度、具体滞后天数对应的相对危险度(RR);单个系数的正负、大小没有业务含义,不用纠结它显不显著,这也是新手最容易踩的坑。
时间趋势系数(ns(time, df = 7 * 3)开头的)也不用逐个解读,看整体显著性即可。这 21 个系数是时间自然样条的基变量,共同拟合了长期趋势和季节波动;输出里大部分系数都高度显著(***),说明数据存在明显的季节 / 长期趋势,我们把它作为混杂因素控制住是正确的,避免把季节波动误判成温度的效应。
星期效应系数(factor(dow)开头的,这部分可以直接解读,用来验证模型是否符合常识。参考组是周一(dow=1),系数代表相对于周一,当天发病数的对数变化量,取指数就是相对倍数。举几个例子:factor(dow)5系数 = 0.413,p<0.001,高度显著:exp(0.413)≈1.51,说明周五的发病数比周一高约 51%,符合疾控 “周中报告数高” 的普遍规律;factor(dow)7系数 = 0.003,p=0.59,不显著:说明周日和周一的发病数没有统计学差异。这部分结果符合业务常识,侧面说明模型的混杂控制是有效的。
3. 拟合优度区:判断模型整体质量
核心看三点:
第一是残差离差远小于零离差 → 模型加入自变量后解释了大部分变异,模型整体有效。
第二是自由度变化是否合理 → 验证自变量数量是否和预期一致。
第三是AIC 值 → 用于横向对比不同参数的模型,越小越好。
![]()
我们看看上面的结果
第一零离差 vs 残差离差。零离差:只有截距项的模型误差,代表 “完全不加入任何自变量时的总变异”;残差离差:加入所有自变量后的剩余误差。我们结果里残差离差(27 万)远小于零离差(138 万),说明加入温度、时间趋势、星期效应后,模型解释了绝大部分数据变异,模型整体是有效的。
第二自由度变化。自由度从 1080 降到 1041,减少了 39 个,刚好对应自变量总数:12 个交叉基变量 + 21 个时间样条变量 + 6 个星期效应变量 = 39 个,参数数量完全对应,没有异常。
第三AIC(赤池信息准则)。数值 279105,单独看没有绝对意义;作用是横向比较不同模型:比如调整自由度、更换滞后期时,AIC 越小,说明模型在 “拟合效果” 和 “简洁性” 之间平衡得越好。
第四离散参数(重点提醒):泊松回归默认离散参数 = 1,但疾控真实数据普遍存在过度离散(方差远大于均值)。简单判断:残差离差 / 残差自由度 = 271557 / 1041 ≈ 260.9,这个值远大于 1,说明存在严重过度离散。
本次是模拟数据,效应设置得比较强,所以离散度很高;真实研究中如果比值 > 1.5,建议换成负二项回归,否则标准误偏小,容易出现假阳性,会被质疑。
4. 效应结果区:这才是研究结论
也就是crosspred()输出的相对危险度(RR)和置信区间。这是 DLNM 分析的最终产出,也是论文里要汇报的核心结果。所有系数、模型配置,最终都是为了算出这个数值。
![]()
解读一下上面的结果。在 14 天的累积效应下,日平均气温 28℃时的传染病发病风险,是参照温度 20℃时的2.71 倍;95% 置信区间为 (2.38, 3.08),区间不包含 1,说明该效应具有统计学意义。
这和我们模拟数据的设定完全吻合:20℃是风险最低点,温度偏离 20℃越远,发病风险越高,28℃属于高温区间,风险显著上升。
4.绘图
最后运行出来的图片包括
![]()
这两个图实质上是一致的,说明的是风速多大以及滞后久对发病的影响如何,也是分布滞后非线性模型主要说明的内容。
![]()
累计效应图也非常重要,主要说明风速大小对其发病的影响。
![]()
这两个个图实质上就是对3D图的切片。
![]()
疑问 1:为什么交叉基是 4×3=12 个变量?
这是 DLNM 最核心的原理问题,对应你输出里cbv1.l1cbv4.l3一共 12 个系数。
原理:基函数的交叉组合。我们构建交叉基时设置了两个关键参数:暴露维度argvar = list(fun = "ns", df = 4) → 用自然样条,4 个自由度,滞后维度arglag = list(fun = "ns", df = 3) → 用自然样条,3 个自由度。
自由度是什么?大白话讲,自由度决定了用多少个基础积木(基函数)来拼接曲线。暴露维度 4 个自由度 → 生成 4 个基变量,用来拼接非线性的暴露 - 反应曲线,滞后维度 3 个自由度 → 生成 3 个基变量,用来拼接平滑的滞后效应曲线。
交叉基的本质,是让两个维度的基函数两两配对组合:暴露维度的每一个基变量,都和滞后维度的每一个基变量相乘,生成一个新的二维基变量。所以总变量数 = 暴露自由度 × 滞后自由度 = 4 × 3 = 12 个。
这就是为什么你会看到 12 个cb开头的系数。它们是一个整体,共同描述了「温度 - 滞后 - 发病」的三维曲面,单独拿出任何一个系数都没有流行病学意义。
疑问 2:时间样条为什么是 21 个系数?
对应你输出里ns(time, df = 7 * 3)1ns(time, df = 7 * 3)21共 21 个系数,原理非常简单:自然样条函数ns()的自由度是多少,就会生成多少个基变量,对应多少个系数。我们代码里写的是ns(time, df = 7*3),也就是总自由度 21。
那么为什么是 7 / 年?这是时间序列研究的行业惯例,用每年 6~8 个自由度,就能很好地控制住季节波动和长期趋势,同时不会过度拟合。3 年数据,每年 7 个自由度 → 7×3=21 个总自由度,对应 21 个基变量、21 个系数。
这 21 个系数共同拟合了一条平滑的时间趋势曲线,用来剥离季节、长期变化对发病数的影响,避免把季节波动误算成温度的效应。同样,单个系数没有意义,看整体显著即可。
疑问 3:为什么拿 28℃和 20℃比?参照温度是随便选的吗?
很多人纳闷:为什么偏偏是 28 和 20,不是 30 和 25?答案很简单:这是我们代码里手动指定的
正式研究不能随便选参照。入门演示我们随便选了 20℃,但真实科研中,必须选择「累积风险最低的最适温度」作为参照。原因很简单:以最低点为参照,U 型曲线的两端高低温 RR 都会大于 1,符合「偏离最适温度风险升高」的逻辑,解读最清晰;如果随便选中位数、均值当参照,可能会出现一半大于 1、一半小于 1 的情况,曲线解读非常别扭,也不符合行业惯例。
疑问 4:为什么会删除 14 个观测?是数据有问题吗?
这是分布滞后模型的正常机制,不是错误。
因为我们设置了最大滞后期 14 天,意味着计算第t天的发病风险,需要用到当天、前 1 天…… 一直到前 14 天,共 15 天的温度数据。数据集的前 14 天,没有足够的历史数据,交叉基会生成缺失值 NA;而glm()默认会自动删除所有含缺失值的行,所以就少了 14 个观测。
数值验证:3 年共 1095 天,删除 14 天后剩 1081 个样本。输出里零模型自由度是 1080 = 1081 - 1,完全吻合,说明一切正常。滞后越长,删除的开头数据越多。只要总样本量足够(几年的日数据),十几天的损失对结果几乎没有影响,无需特殊处理。
疑问 5:残差离差那么大,是不是模型不好?
很多人看到残差离差二十多万,会觉得模型很差。其实不是,离差的绝对值和样本量、发病基数有关,关键看离散程度
泊松回归默认离散参数 = 1,判断标准很简单:
![]()
你的结果:271557 / 1041 ≈ 260.9,远大于 1,说明存在严重的过度离散。
本次是模拟数据,我们设定的温度效应很强,所以离散度很高;真实疾控数据也普遍存在过度离散,只是幅度没这么大。过度离散的危害:会导致标准误偏小,p 值偏低,容易出现假阳性,论文中会被审稿人质疑。一般建议:离散度 > 1.5 时,就换成负二项回归。
最后总结:看 DLNM 结果的正确顺序,给大家整理了一个新手检查清单,照着一步步来,就不会乱:
- 看 Call:确认变量、分布没写错,避免低级失误;
- 看样本量:确认删除的滞后数据数量和设定的最大滞后期匹配;
- 看控制变量:星期效应、时间趋势整体显著,说明混杂控制有效;
- 算离散度:残差离差 / 自由度,大于 1.5 就换成负二项回归;
- 算 RR 和置信区间:用crosspred()提取核心效应,这才是最终结论;
- 不要死磕交叉基单个系数:它们是组合积木,单独看没有意义。
![]()
![]()
编辑:普通疾控人 | 审核:诗酒趁年华
文章来源 | 原创
说明 | 转载只为分享,如有侵权联系删除
©版权声明 | 部分信息和图片来自公开网络
转载请注明
再次转载请注明出处
![]()
科普健康 | 宣传疾控
本号为多位疾控机构从业者运营
重点关注国内外健康事件
致力于疾控科普
在做好科普服务大众的同时
做好疾控机构的宣传
让更多的人了解疾控,拥抱健康
欢迎加「小编」微信(cdcjkr126com)
本文具体说明
本文为原创内容,文章为个人理解所学,不涉及疫情信息及内部保密数据,发表的目的为自我总结及给有需求的人士学习使用。如有不妥之处,欢迎联系小编修改、删除。
更多精彩视频,尽在“CDC疾控人”视频号
![]()
特别声明:以上内容(如有图片或视频亦包括在内)为自媒体平台“网易号”用户上传并发布,本平台仅提供信息存储服务。
Notice: The content above (including the pictures and videos if any) is uploaded and posted by a user of NetEase Hao, which is a social media platform and only provides information storage services.