关于算法:R语言分布滞后线性和非线性模型DLNM分析空气污染臭氧温度对死亡率时间序列数据的影响

6次阅读

共计 8513 个字符,预计需要花费 22 分钟才能阅读完成。

原文链接 http://tecdat.cn/?p=23947 

摘要

散布滞后非线性模型(DLNM)示意一个建模框架,能够灵便地形容在工夫序列数据中显示潜在非线性和滞后影响的关联。该方法论基于穿插基的定义,穿插基是由两组根底函数的组合示意的二维函数空间,它们别离指定了预测变量和滞后变量的关系。本文在 R 软件实现 DLNM,而后帮忙解释后果,并着重于图形示意。本文提供指定和解释 DLNM 的概念和实际步骤,并举例说明了对理论数据的利用。

关键字:散布滞后模型,工夫序列,平滑,滞后效应,R。

1. 简介

统计回归模型的次要目标是定义一组预测变量与后果之间的关系,而后预计相干影响。当依赖项显示某些滞后影响时,会进一步减少复杂性:在这种状况下,预测变量的产生(咱们称其为裸露事件)会在远远超出事件周期的工夫范畴内影响后果。此步骤须要定义更简单的模型以表征关联,并指定依赖项的工夫构造。

1.1 概念框架

对滞后效应的适当统计模型的阐明及其后果的解释,有助于建设适当的概念框架。这个框架的次要特点是定义了一个额定的维度来形容关联,它指定了裸露和后果之间在滞后维度上的工夫依赖性。这个术语,借用了工夫序列剖析的文献,代表了评估影响滞后时裸露事件和后果之间的工夫距离。在长时间裸露的状况下,数据能够通过等距时间段的划分来结构,定义一系列裸露事件和后果实现。这种划分也定义了滞后单位。在这个工夫构造中,裸露 - 反馈关系能够用两种相同的观点中的任何一种来形容:咱们能够说一个特定的裸露事件对将来的多个后果产生影响,或者说一个特定的后果能够用过来多个裸露事件的奉献来解释。而后,能够应用滞后的概念来形容向前(从固定后果到将来后果)或向后(从固定后果到过来的后果)的关系。

最终,滞后效应统计模型的次要特色是它们的二维构造:该关系同时在预测变量的通常空间和滞后的维度上进行形容。

1.2 散布滞后模型

最近,在评估环境压力因素的短期影响的钻研中曾经解决了滞后影响的问题:一些工夫序列钻研报告说,裸露于高水平的净化或极其温度会在其产生后的几天内继续影响衰弱(Braga 等,2001;Goodman 等,2004;Samoli 等,2009;Zanobetti 和 Schwartz,2008)。

给定定义的数据工夫构造和简略的滞后维度定义,工夫序列钻研设计可提供多种劣势来解决滞后影响,其中工夫划分是由等距离和有序的工夫点间接指定的。在这种状况下,滞后效应能够用散布滞后模型(DLM)来优雅地形容,该模型最后是在计量经济学中开发的(Almon 1965),最近在环境因素钻研中用于量化衰弱效应(Schwartz 2000; Zanobetti et al。2000; 2007)。Muggeo 和 Hajat,2009 年)。通过这种办法,能够应用多个参数来解释在不同时滞下的影响,从而将单个裸露事件的影响散布在特定的时间段内,

1.3 本文目标

统计环境 R 提供了一组用于指定和解释 DLNM 后果的工具。本文的目标是提供该程序包函数的全面概述,包含函数的具体摘要以及以理论数据为例的示例。该示例波及 1987-2000 年期间两个环境因素(空气污染(臭氧)和温度)对死亡率的影响。在本文中,我重新考虑了定义 DLNM,预测成果并借助图形函数解释后果的次要概念和实际步骤。

2. 非线性和滞后效应

在本节中,我介绍了工夫序列模型的根本公式,而后介绍了形容非线性效应和滞后效应的办法,后者通过简略 DLM 的模型来形容。

2.1 根本模型

工夫序列数据的模型通常能够示意为:

其中 µt≡E(Yt),Yt 是 t = 1 时的一系列后果 …,n,假如来自指数族的散布。函数 sj 指定变量 xj 和线性预测变量之间的关系,该变量由参数向量 βj 定义。变量 uk 蕴含具备由相关系数 γk 指定的线性效应的其余预测变量

之前形容的数据说明性示例中,后果 Yt 是每日死亡计数,假设是泊松散布,其中 E(Y)= µ,V(Y)= φµ。

臭氧和温度的非线性和滞后影响通过函数 sj 建模,该函数定义了预测变量和滞后变量两个维度之间的关系

2.2 非线性裸露 - 反馈关系

DLNM 开发的第一步是定义预测变量空间中的关系。通常,非线性裸露 - 反馈依赖性通过适当的函数 s 在回归模型中示意。在齐全参数化的办法中,提出了几种不同的函数,每个函数都具备不同的假如和灵活性。次要抉择通常依赖于形容润滑曲线的函数,例如多项式或样条函数(Braga 等,2001;Dominici 等,2004)。对于线性阈值参数化的应用(Muggeo 2010; Daniels et al。2000); 或通过虚构参数化进行简略分层。

所有这些函数都对原始预测变量进行了转换,以生成蕴含在模型中作为线性项的一组转换变量。相干的根底函数包含原始变量 x 的一组齐全已知的转换,这些转换生成一组称为根底变量的新变量。代数示意能够通过以下形式给出:

定义 DLNM 的第一步是在函数 mkbasis()中执行的,该函数用于创立根底矩阵 Z。此函数的目标是提供一种通用的形式来蕴含 x 的非线性效应。举例来说,我建设了一个将所选基函数利用于向量 的基矩阵:

R> mkais(1:5, tpe = "s", df = 4, egree = 2, cenvlue = 3)

后果是一个列表对象,存储根底矩阵和定义该矩阵的自变量。在这种状况下,所选基准是具备 4 个自由度的二次样条,由参数类型 df 和度定义。

能够通过第二个参数类型抉择不同类型的根底。可用的选项是天然三次方或简略的 B 样条(类型 =“ns”或“bs”);虚构变量层;多项式(“poly”);阈值类型的函数和简略的线性(“lin”)。参数 df 定义了根底的维数(根底的列数,基本上是转换后的变量的数目)。该值可能取决于参数“结点”。如果未定义,则默认状况下将结搁置在等距的分位数上。自变量度数抉择“bs”和“poly”的多项式度数。

参数 cen 和 cenvalue 用于使连续函数(类型“ns”,“bs”,“poly”和“lin”)的基准居中,如果未提供 cenvalue,则默认为原始变量的均值。

2.3 滞后效应

定义 DLNM 的第二步是指定函数,以对附加滞后维度中的关系进行建模,以实现滞后成果。在这种状况下,给定工夫 t 的后果 Yt 能够用过来的裸露量 xt- L 来解释。给定最大滞后 L 时,附加滞后维度能够由 n×(L +1)矩阵 Q 示意,例如:

简略的 DLM 应用形容后果与滞后危险之间的依赖关系的函数来容许线性关系的滞后效应。

第二步通过函数 mklagbasis()进行,该函数调用 mkbasis()来构建根底矩阵 C。例如:

R> mkgbais(mxlag =5,type ="strta", kots = c(2, 4))

在此示例中,在通过第一个参数 maxlag 将最大滞后固定为 5 之后,滞后向量 0:maxlag 对应于,将主动创立并利用所选函数。

3. 定义 DLNM

DLNM 标准的最初一步波及同时定义预测器和滞后两个维度中的关系。只管非线性和滞后效应的术语不同,但这两个过程在概念上是类似的:定义示意相干空间中关系的根底。

而后,通过穿插基的定义来指定 DLNM,穿插基是二维函数空间,同时形容了沿预测变量范畴及其滞后维度的依存关系。首先,抉择 x 的基函数得出 Z,而后为 x 的每个基变量创立附加的滞后维度,从而生成一个 数组 R˙。通过定义的 C,DLNM 能够示意为:

抉择穿插基等于如上所述抉择两组基函数,将其组合以生成穿插基函数。这是通过函数 crossbasis()执行的,该函数调用函数 mkbasis()和 mklagbasis()别离生成两个根本矩阵 Z 和 C,而不是通过张量积将它们组合起来以产生 W。能够应用此函数指定臭氧和温度的两个穿插基。相干代码为:

basi.o3 <- crossbasis(o3 varype= "hthr"
+ vnots = 40, laty = "sata", lanot = c(2,6), mag= 10)
bai.te <- crossbasis(tmp varype = "bs",
+ vrgre  3, vad = 6 cevalu = 25 ladf = 5, malag = 30)

在此示例中,臭氧的穿插基包含一个预测空间的阈值函数,线性关系超过 40.3 µgr / m3,并且虚构参数化假如沿滞后 0 -1、2- 5 和 6 -10 的层具备恒定的散布滞后效应。相比之下,温度的选项是:以 25 摄氏度为核心的 6 自由度的立方样条(默认为等距的结点),以及以 5 自由度的立方样条(默认为 lagtype =“ns”)(结为 25℃)。默认状况下,最多 30 个滞后。

如果未设置核心值,则默认的中心点是预测变量的平均值(例如,对于上述温度的穿插基,温度为 25℃)。该值代表来自 DLNM 的预期成果的参考。参考值的抉择不影响模型的拟合,并且能够依据解释问题抉择不同的值。

这些抉择能够通过函数 summary()进行查看。例如:

R> summary(basis.temp)

为了预计相应参数 η,能够在通用回归函数的模型公式中包含穿插基矩阵。在该示例中,最终模型还包含一个天然立方样条,以模仿季节性趋势和长期趋势重量,代码是:

odel <- glmdeath ~ bais.temp+ basis.o +ns(tim 7 * 14)  dw,
+ fmily = quasiposson())

4. 依据 DLNM 进行预测

如第 3 节所示,DLNM 的标准波及裸露序列的简单参数化,然而参数 η 的估算是应用常见的回归命令进行的。然而,定义沿两个维度的关系的此类参数的含意并不简略。能够通过预测在具备适当裸露值和 L + 1 滞后的网格上的滞后特定成果来辅助解释。此外,能够通过将滞后特定奉献相加来计算从滞后 L 到 0 继续裸露所预测的总体成果。预测的成果通过函数 crosspred()在 dlnm 中计算。以下代码在示例中计算了对臭氧和温度的预测:

pre.o <- crosspred(basis, odel at = c(0:6,0., .3))

传递给 crosspred()的前两个参数是“crossbasis”类的对象和用于预计的模型对象。像下面的第一个示例一样,能够通过 at 参数间接指定必须为其预测成果的裸露值向量。在这里,我抉择了臭氧中从 0 到 65 µgr / m3 的整数,再加上所选阈值的值和 10 个单位以上的值(别离为 40.3 和 50.3 µgr / m3)。而后,该函数调用 crossbasis()来构建预测基准,并依据模型中的参数生成预测成果和标准误差。后果是“crosspred”类的列表对象,该对象存储了预测的成果。它包含滞后效应矩阵和总体效应向量,以及相应的标准误差矩阵和向量。如第 5 节所示。例如,臭氧减少 10 个单位的总体成果示意为 RR 和 95%置信区间,能够通过以下公式得出:

R> pred.o3$allRRfit\["50.3"\]

R> cbind(lRlow,alRigh)\["50.3",\]

 

5. 形容 DLNM

由 DLNM 估算的二维裸露 - 反馈关系可能难以概括。关联的图形示意提供了个别形容。调用高级函数 plot.default(),persp()和 filled.contour()来生成散点图,3- D 和等高线图。例如,臭氧和死亡率之间的关系能够通过 RR 进行总结,即每次滞后会比阈值高出 10 µgr / m3。该图如图 1(左)所示,可通过以下形式取得:

图 1:在阈值(40.3 µgr / m3)以上的臭氧减少 10 个单位时,滞后效应(左)和总体效应(右)对死亡率的影响。

R> plot(re.o3)

参数 ptype =“slices”指定图的类型,在这种状况下,沿着滞后空间在预测值 var = 50.3 处的预测成果矩阵的切片,对应于在 40.3 µgr / m3 的阈值之上减少了 10 个单位。自变量 ci 示意置信区间的图类型。如果应用 cumul = TRUE,则绘制累积成果。

依据概念定义,能够应用两种不同的观点来读取图 1 中的左图:它示意在第 t 天以 50.3 µgr / m3 的臭氧进行单次裸露后,将来每一天的危险减少。

或者,能够绘制总体成果,该总体成果是通过应用参数 ptype =“overall”将滞后效应相加得出的:

R> plot(pred)

图 2:温度和全因死亡率之间的裸露 - 反馈关系的三维图,以 25°C 为参考。

一种更具体的办法来示意温度与死亡率之间的平滑关系,其中样条函数已用于定义这两个维度的相关性。能够应用 3 - D 和等高线图对这种简单的依赖关系进行个别形容,该图阐明了由预测成果的整个网格给出的成果外表。所示的图是通过以下形式取得的:

R> plot(pred.temp, "contour")

参考点(此处为 25℃)是 crossbasis 函数在 crossbasis()中核心的值。

三维图或等高线图提供了关系的全面摘要,但在示意特定预测值或滞后值的影响方面的能力无限。上面给出了更全面的图,该图片通过以下形式取得:

R> plot(pred.temp, "slices
+ ci.g , ltensity =20 colr(0)))

图 3(左)显示了由 plot()和 lines()中的参数 var 抉择的温度值的预测滞后效应影响。另外,图 3(右)显示了针对特定滞后的沿温度的预测效应的多重曲线图(左),以及图 3(右)中绘制的雷同滞后效应,以及 99%的置信区间。

这些图表显示了低温和高温影响的不同模式,低温的影响十分强烈且迅速,高温影响更为提早,在最后的滞后中为负。

6. 建模策略

DLNM 框架提供了机会,能够通过为预测变量和滞后变量两个维度中的每个维度抉择根本函数来指定宽泛的模型抉择。后面各节中阐明的示例代表了一种潜在的建模代替办法。为了探讨该办法的灵活性以及模型抉择的相干问题,上面显示了与不同模型的比拟,以预计与温度的关联。具体来说,为预测变量的空间抉择多项式和档次函数,同时放弃雷同的天然三次样条,以模仿长达 30 天的滞后散布的滞后曲线。指定穿插根底,运行模型并预测成果的代码为:

R> basis.temp2 <- crossbasis(emp, vrtpe = "poly",
R> model2 <- update(mdel, .~. - bsis.emp + baiste2)
R> model3 <- updat(model .~. -bais.tmp + bass.mp3)

对于预测变量,第一种办法倡议应用与第 5 节中的原始三次样条雷同的自由度的多项式函数。第二种模型基于一个更简略的双阈值函数,将单个阈值置于 25°C,之前确定为最低死亡率。此抉择还便于模型比拟,因为这是其余两个连续函数的中心点。这三个模型预计的总体成果显示在由代码产生的图 4(左)中:

R> plot(pred.temp, "overall", ylim = c(0.5, 2.5), ci = "n", lwd = 1.5,
+ main = "Overall effect")
R> lines(pred.temp2, "overall", col = 3, lty = 2, lwd = 2)
R> lines(pred.temp3, "overall", col = 4, lty = 4, lwd = 2)
R> legend("top", c("natural spline", "polynomial", "double threshold"),
+ col = 2:4, lty = c(1:2, 4), lwd = 1.5, inset = 0.1, cex = 0.8)

正如预期的那样,代替模型会产生不同的后果。特地是,如果与具备等距结点的三次样条进行比拟,则多项式模型会预计出高温的“摆动”关系。取而代之的是,这两个函数提供了十分靠近的低温影响估算值。相同,尽管双阈值模型的线性假如仿佛足以模仿高温的依赖性,但有一些证据表明,这种办法往往会低估热的影响。预计的散布滞后曲线的第二次比拟如图 4 所示(右),如下所示:

R> plot(pred, slices", va =32, im =95 .2="n"

只管在所有三个模型中都为滞后空间抉择了完全相同的函数,但对预测变量的不同抉择提供了散布滞后曲线的不同估计值,与 32°C 的参考点相比,代表了 32°C 的影响。

图 4:温度为 32°C 时的总体效应(左)和滞后特异性效应(右)对 3 种代替模型的全因死亡率的影响(以 25°C 为参考)。芝加哥 1987-2000。

特地是,样条曲线和多项式模型会产生十分类似的成果(正如预期的那样,思考到高温度尾部曲线在其余维度上的拟合简直雷同),而双阈值模型的曲线显示出截然不同的形态。具体而言,因为不足此模型的灵活性,因而暗示播种成果(较长滞后的负预计)可能示意伪像。

不足通用规范,无奈在可用的抉择中抉择总结关联的最佳模型,从而加重了对各种代替产品的规格要求的这种丰富性。在下面的示例中,我对样条线模型体现出了显著的偏爱。这种抉择既基于对函数属性的理解,例如灵活性和稳定性,又基于给出图 4 所示后果的正当论据。然而,该论断是有问题的,而不是基于牢靠的和个别的统计抉择规范。此外,论断是基于几个先验的抉择,就像阈值地位或结数或多项式次数一样。

通常,在 DLNM 中,能够形容两个不同的抉择级别。第一个波及不同函数的标准。如上所示,该抉择应既基于假如的裸露反馈形态的合理性,又基于复杂性,可概括性和易于解释之间的折衷。第二级重点关注特定函数内的不同抉择,例如用于定义样条曲线基的结的数量和地位。后者更难解决,只管不是 DLNM 开发所固有的。一些钻研人员在工夫序列剖析中钻研了这个问题,提出了基于信息准则(Akaike,Bayesian 和其余变体),偏自相干或(狭义)穿插验证的办法(Peng 等,2006;Baccini 等,2006)。2007)。用户能够在 DLNM 中利用雷同的办法,然而他应该记住,这些模型的二维性质带来了额定的复杂性,例如最大滞后的定义。此外,对于执行不同准则的根据还不是结论性的(Dominici 等人,2008 年)。须要进一步钻研以提供无关 DLNM 中模型抉择的一些领导。

能够倡议应用其余办法。Muggeo(2008)提出了一个模型,该模型具备对预测变量空间进行束缚的分段参数化,以及基于惩罚性样条的双重惩办基于散布滞后的参数化。此办法包含主动抉择阈值和散布滞后曲线的平滑度,并且已在 R(Muggeo 2010)中齐全实现。这种办法与灵便的 DLNM 的比拟能够放宽对预测变量维度上形态的假如,从而能够提供无关此关系的其余一些见解。

7. 数据要求

本文介绍的 DLNMs 框架是为工夫序列数据开发的。(1)中根本模型的个别表达式容许将此办法利用于(狭义)线性模型(GLM)中的任何族散布和链接函数,并扩大到狭义加法模型(GAM)或基于狭义预计方程的模型(GEE)。然而,DLNM 的以后实现须要一系列等距,残缺和有序的数据。

还应用选定滞后时间段中蕴含的先前察看值来计算一系列转换变量中的每个值。因而,将转换变量中的第一个最大滞后观测值设置为 NA。容许在 x 中短少值,然而因为雷同的起因,将雷同且下一个 maxlag 转换后的值设置为 NA。只管正确,但对于零散的缺失观测值存在的较长滞后工夫的 DLNM,这可能会产生计算问题。在这种状况下,能够思考一些插补办法。

dlnm 的次要长处之一是,用户能够应用规范回归函数执行 DLNM,只需在模型公式中包含穿插基矩阵即可。通过函数 lm(),glm()或 gam(),能够间接应用它。然而,用户能够与数据的工夫序列构造兼容地利用不同的回归函数。这些函数应该具备针对 coef()和 vcov()的办法,或者用户必须提取参数并将其蕴含在 crosspred()的参数 coef 和 vcov 中(请参见第 4 节)。

8. 最终论断

DLNM 类代表形容形容非线性效应和滞后效应的景象的对立框架。该模型系列的次要长处是在一个独特的框架中对立了许多以前的办法来解决滞后效应,还为关系提供了更灵便的抉择。DLNM 的标准仅波及抉择两个基数以生成(5)中的穿插基函数,例如,包含线性阈值,档次,多项式和样条变换。

穿插基和参数估计的拆散提供了多个长处。首先,如示例中所示,能够通过穿插基函数转换多个显示滞后成果的变量,并将其蕴含在模型中。其次,能够应用规范回归命令进行预计,并应用默认的诊断工具和相干函数集。更重要的是,此实现提供了一个开放平台,能够在其中实现应用不同回归命令指定的其余模型,来帮忙在其余状况下或钻研设计中开发方法。


最受欢迎的见解

1. 在 python 中应用 lstm 和 pytorch 进行工夫序列预测

2.python 中利用长短期记忆模型 lstm 进行工夫序列预测剖析

3. 应用 r 语言进行工夫序列(arima,指数平滑)剖析

4.r 语言多元 copula-garch- 模型工夫序列预测

5.r 语言 copulas 和金融工夫序列案例

6. 应用 r 语言随机稳定模型 sv 解决工夫序列中的随机稳定

7.r 语言工夫序列 tar 阈值自回归模型

8.r 语言 k -shape 工夫序列聚类办法对股票价格工夫序列聚类

9.python3 用 arima 模型进行工夫序列预测

正文完
 0