R語言分布滯后非線性模型(DLNM)空氣污染研究溫度對死亡率影響建模應(yīng)用|附代碼數(shù)據(jù)
全文下載鏈接:http://tecdat.cn/?p=23564
最近我們被客戶要求撰寫關(guān)于DLNM的研究報告,包括一些圖形和統(tǒng)計輸出。
在本文中,環(huán)境應(yīng)激源往往表現(xiàn)出時間上的滯后效應(yīng),這就要求使用足夠靈活的統(tǒng)計模型來描述暴露-反應(yīng)關(guān)系的時間維度。在此,我們開發(fā)了分布式滯后非線性模型(DLNM),這是一個可以同時代表非線性暴露-反應(yīng)依賴性和滯后效應(yīng)的建模框架。這種方法是基于 "交叉基準(zhǔn) "的定義,這是一個雙維的函數(shù)空間,它同時描述了沿預(yù)測空間和其發(fā)生的滯后維度的關(guān)系形狀
通過這種方式,該方法為以前用于該環(huán)境的一系列模型提供了一個統(tǒng)一的框架。為了說明這個方法,我們用DLNMs的例子來表示溫度和死亡率之間的關(guān)系,使用1987-2000年期間國家發(fā)病率、死亡率和空氣污染研究中的數(shù)據(jù)。
簡介
有時特定暴露事件的影響并不局限于觀察到的那段時間,而是在時間上有所滯后。這就帶來了一個問題,即對暴露事件與未來一系列結(jié)果之間的關(guān)系進(jìn)行建模,指定事件發(fā)生后不同時間的影響分布(定義的滯后期)。最終,這一步需要定義暴露-反應(yīng)關(guān)系的額外滯后維度,描述影響的時間結(jié)構(gòu)。
在評估環(huán)境應(yīng)激源的短期影響時,這種情況經(jīng)常發(fā)生:一些時間序列研究報告稱,暴露在高水平的空氣污染或極端溫度下,會在發(fā)生后的幾天內(nèi)影響健康。此外,當(dāng)一個應(yīng)激源主要影響一批脆弱的個體時,就會出現(xiàn)這樣的現(xiàn)象,這些個體的事件只因暴露的影響而提前了短暫的時間。
在已經(jīng)提出的處理之后效應(yīng)的各種方法中,分布式滯后模型(DLM)發(fā)揮了主要作用,最近在空氣污染和溫度研究中被用來量化健康效應(yīng)。這種方法的主要優(yōu)點(diǎn)是,它允許模型包含暴露-反應(yīng)關(guān)系的時間過程的詳細(xì)表述,這反過來又提供了對存在滯后貢獻(xiàn)或收獲的總體效應(yīng)的估計。
雖然傳統(tǒng)的DLMs適合于描述線性效應(yīng)的滯后結(jié)構(gòu),但在用于表示非線性關(guān)系時,它們顯示出一些局限性。我們提出了一個解決方案,進(jìn)一步放寬對關(guān)系的假設(shè),并將這種方法擴(kuò)展到分布式滯后非線性模型(DLNM),這是一個模型家族,可以以靈活的方式描述沿預(yù)測器空間和其發(fā)生的滯后維度同時變化的效應(yīng)。通過這種方式,DLNM類也為現(xiàn)有的較簡單的方法提供了一個統(tǒng)一的框架。
DLNMs以前只在流行病學(xué)方面進(jìn)行過簡單的描述:本文的目的是嚴(yán)格地發(fā)展這種方法,并描述在統(tǒng)計軟件R中專門編寫的軟件包dlnm中的實(shí)現(xiàn),提供一個使用真實(shí)數(shù)據(jù)集的應(yīng)用實(shí)例。我們簡要描述了時間序列分析中使用的基本模型,并介紹了基礎(chǔ)的概念,作為描述變量和因變量之間非線性關(guān)系的一般方法。我們概述了在時間上滯后效應(yīng)的復(fù)雜性,并提供了一個簡單的DLMs的一般表示。然后說明了這種方法在溫度對死亡率影響的建模中的應(yīng)用。最后我們提供了一些討論并提出了可能的進(jìn)一步發(fā)展。
基本模型
一般的表示法
描述結(jié)果Yt的時間序列(t=1,...,n)的一般模型表示方法為

其中,≡E(Y ),g是一個單調(diào)的函數(shù),Y被假定來自屬于指數(shù)族的分布。函數(shù)sj表示變量x j和線性預(yù)測器之間的平滑關(guān)系,由參數(shù)向量bj定義。變量uk包括其他預(yù)測因子,其線性效應(yīng)由相關(guān)系數(shù)k指定。函數(shù)sj也可以通過基于廣義加性模型的非參數(shù)方法來指定。然而,在目前的發(fā)展中,我們依靠的是一種完全的參數(shù)化方法。
在環(huán)境因素的時間序列分析中,結(jié)果Yt通常是每日計數(shù),假定來自所謂的過度分散泊松分布。這些研究利用了過去幾年中統(tǒng)計方法的重大改進(jìn),來量化空氣污染的短期影響。通常,這些方法包括一個平滑的時間函數(shù),以識別隨時間緩慢變化的混雜因素的影響,表現(xiàn)為季節(jié)性或長期趨勢。也包括溫度和濕度等氣象因素的非線性影響。分類變量,如一周中的幾天或年齡組被作為因素進(jìn)行模擬。盡管空氣污染通常用線性關(guān)系來描述,但為了評估非線性效應(yīng),這一假設(shè)可以放寬。
在這里,我們關(guān)注的是一個一般的函數(shù)s,它指定了預(yù)測因子x的潛在非線性和滯后效應(yīng),通常指的是空氣污染或溫度,但不失一般性。
基函數(shù)
x和g()之間的關(guān)系由s(x)表示,它作為一個線性項的總和包含在廣義線性模型的線性預(yù)測器中。這可以通過選擇一個基數(shù)來實(shí)現(xiàn),基數(shù)是一個函數(shù)空間,我們認(rèn)為s是其中的一個元素。相關(guān)的基函數(shù)包括一組完全已知的原始變量x的變換,產(chǎn)生一組新的變量,稱為基變量。估計關(guān)系的復(fù)雜性取決于基數(shù)的類型和它的維度。幾個不同的基礎(chǔ)函數(shù)被用來描述環(huán)境因素對健康的潛在非線性影響,其選擇取決于對關(guān)系形狀的假設(shè)、調(diào)查的具體目的所要求的近似程度以及解釋問題。在完全參數(shù)化的方法中,主要的選擇通常依賴于描述平滑曲線的函數(shù),如多項式或樣條函數(shù),或使用線性閾值參數(shù)化,由截斷的線性函數(shù)(x-)+表示,當(dāng)x>時等于(x-),否則等于0。上述簡單模型的一般表示方法為

滯后效應(yīng)
額外維度
在存在滯后效應(yīng)的情況下,在給定時間t的結(jié)果可以用過去的暴露xt-來解釋,滯后代表暴露和反應(yīng)之間所經(jīng)過的時間。一個相對簡單的方法是對有序暴露的原始向量x進(jìn)行轉(zhuǎn)換,得出n×(L+1)矩陣Q,如
這一步規(guī)定了暴露-反應(yīng)關(guān)系的額外滯后維度。最終,這里提出的建??蚣艿哪康氖峭瑫r描述兩個維度的依賴關(guān)系:通常的預(yù)測器空間和新的滯后維度。
分布滯后模型
當(dāng)假設(shè)有線性關(guān)系時,滯后效應(yīng)可以自然地用分布式滯后模型(DLM)來描述。這種方法允許將單一暴露事件的影響分布在一個特定的時間段內(nèi),用幾個參數(shù)來解釋不同滯后期的貢獻(xiàn)。這些模型已被廣泛用于評估環(huán)境因素的滯后效應(yīng)。最簡單的表述是無約束的DLM,通過為每個滯后期加入一個參數(shù)來指定。不幸的是,由于相鄰天數(shù)的暴露之間的高度相關(guān)性以及由此產(chǎn)生的模型中的串聯(lián)性,對特定滯后期效應(yīng)的估計精度往往非常差。
為了使分布式滯后曲線的估計更加精確,可以施加一些限制條件,例如假設(shè)滯后區(qū)間內(nèi)的效應(yīng)不變,或者使用連續(xù)函數(shù)如多項式或樣條來描述平滑曲線。一個以前L天暴露量的移動平均數(shù)為預(yù)測因子的簡單模型可以被視為DLM的一個特例:這樣的模型已被廣泛用于空氣污染流行病學(xué)領(lǐng)域,有時也被用于量化溫度的影響。這類模型以前只給過多項式 DLMs。有可能制定一個更簡單和通用的DLM定義,其中沿滯后期的分布效應(yīng)的形狀由一個適當(dāng)?shù)幕A(chǔ)指定。在矩陣記號中

我們可以定義

通過構(gòu)建每個滯后期的隱含線性效應(yīng)b,可以幫助解釋估計的參數(shù)g?,具體如下。

分布式滯后非線性模型?
有完善的方法來描述簡單滯后模型的靈活暴露-反應(yīng)關(guān)系,或者是簡單線性效應(yīng)的靈活DLM,但很少同時對這兩個部分進(jìn)行建模。已經(jīng)提出了描述非線性效應(yīng)的擴(kuò)展方法,通過對閾值或分段函數(shù)的每個項或?qū)€性和二次項分別應(yīng)用約束矩陣C,可以構(gòu)建一個DLM。盡管如此,這些方法在描述這種復(fù)雜的依賴關(guān)系的能力方面仍然有些局限。通過產(chǎn)生一個新的模型框架,可以描述預(yù)測器空間和滯后期的非線性關(guān)系,從而實(shí)現(xiàn)一個有用的概括,這就是DLNM家族。
交叉基的概念
雖然DLNM的代數(shù)符號可能相當(dāng)復(fù)雜,涉及到三維數(shù)組,但基本概念是建立在交叉基數(shù)的定義上的,是很簡單的。交叉基點(diǎn)可以被描繪成一個雙維的函數(shù)空間,同時描述沿x的關(guān)系的形狀及其分布的滯后效應(yīng)。選擇交叉基點(diǎn)相當(dāng)于選擇兩組基函數(shù),它們將被組合起來產(chǎn)生交叉基函數(shù)。
DLNM
為了對我們所考慮的兩個空間的關(guān)系形狀進(jìn)行建模,我們需要同時應(yīng)用描述的兩個轉(zhuǎn)換。首先,如(2)所述,我們?yōu)閤選擇一個基礎(chǔ)來定義預(yù)測器空間中的依賴關(guān)系,指定Z。然后,如(3)所述,我們?yōu)榇鎯υ赯中的x的每個派生基變量創(chuàng)建額外的滯后維度。該結(jié)構(gòu)是對稱的,即兩個轉(zhuǎn)換的順序可以顛倒,將基函數(shù)直接應(yīng)用于矩陣Q的每一列。

解釋DLNM
盡管參數(shù)化很復(fù)雜,但對DLNM參數(shù)的估計和推斷并不比任何其他廣義線性模型產(chǎn)生更多的問題,而且在指定交叉基變量后,可以用普通的統(tǒng)計軟件進(jìn)行。然而,雖然(4)中較簡單的DLM的解釋是直接的,包括報告(6)中每個滯后的估計線性效應(yīng)b?,但更復(fù)雜的DLNM的結(jié)果與平滑的非線性依賴關(guān)系很難總結(jié)。一個解決方案是為每個滯后期和適當(dāng)?shù)谋┞吨到⒁粋€預(yù)測網(wǎng)格,使用三維圖來提供沿兩個維度變化的影響的總體情況。
預(yù)測網(wǎng)格,用預(yù)測效果E的m×(L+1)矩陣和相關(guān)的標(biāo)準(zhǔn)誤差Esd矩陣表示,可以使用估計系數(shù)的向量g?,從包括交叉基函數(shù)矩陣W的擬合模型中計算得出。

并且,給定V(g?)為估計系數(shù)的方差矩陣

這個網(wǎng)格對于計算滯后p的暴露效果或滯后x p的暴露效果的估計很有用,只需分別取e-p和ex p-。最后,通過將不同滯后期的所有貢獻(xiàn)相加,可以計算出總體效應(yīng)的估計值。矢量etot和相關(guān)的標(biāo)準(zhǔn)誤差esd tot,由每個滯后期的貢獻(xiàn)相加得到,說明整個滯后期的暴露效果。


應(yīng)用
數(shù)據(jù)和模型選擇
我們應(yīng)用DLNMs來研究1987-2000年期間溫度對總體死亡率的影響。數(shù)據(jù)集來自國家發(fā)病率、死亡率和空氣污染研究。
它包括5114個總體和特定病因的死亡率、天氣和污染數(shù)據(jù)的每日觀測。
分析基于(1)中的模型,通過準(zhǔn)泊松族的廣義線性模型進(jìn)行擬合,在控制混雜因素方面有以下選擇:每年有7個自由度(df)的時間自然立體樣條,以描述長期趨勢和季節(jié)性;每周一天的指標(biāo)變量;滯后0-1的露點(diǎn)溫度平均值有3個自由度的自然立體樣條;滯后0-1的臭氧和CO的平均值的線性項。
glm(death?~?ns.basis?+?ns(dp01,df=3)?+?dow?+?o301?+?co01?+ns(date,df=14*7),family=quasipoisson(),?data)
這些選擇是根據(jù)幾篇關(guān)于時間序列分析的方法學(xué)和實(shí)質(zhì)性論文。通過選擇兩個基點(diǎn)來描述溫度和滯后期空間的關(guān)系,研究了平均溫度的影響;我們說明了一個靈活的模型,用自然立體樣條來描述每個維度的關(guān)系。結(jié)點(diǎn)被放置在溫度范圍內(nèi)等距的數(shù)值上,以便在尾部有足夠的靈活性,而在滯后期的對數(shù)尺度上等距放置,以便在分布式滯后期曲線的第一部分有更多的靈活性,因?yàn)樵谀抢镱A(yù)計會有更多的變化。最大的滯后期L被設(shè)定為30天。為了比較,我們用前幾天溫度的移動平均數(shù)擬合了比較簡單的模型。
我們根據(jù)修改后的赤池和貝葉斯信息標(biāo)準(zhǔn)來選擇結(jié)的數(shù)量,它定義了每個維度上的df,用于通過準(zhǔn)似然法擬合的具有過度分散反應(yīng)的模型,具體內(nèi)容如下。

所有的分析都是用R軟件進(jìn)行的。
#?3-D?圖?plot(ns.pred,label="Temperature")

點(diǎn)擊標(biāo)題查閱往期內(nèi)容

R語言分布滯后線性和非線性模型(DLMs和DLNMs)分析時間序列數(shù)據(jù)

左右滑動查看更多

01

02

03

04

結(jié)果
當(dāng)用于比較不同的建模選擇時,QAIC導(dǎo)致了一個相對復(fù)雜的模型,預(yù)測器空間有11df,滯后維度有5df,總共有55個參數(shù)用于定義關(guān)系。相比之下,QBIC表明是一個5×5df的模型,用25df來描述總體效果。由于對DLNM框架內(nèi)這些標(biāo)準(zhǔn)的表現(xiàn)沒有任何了解,我們選擇了后者作為我們的最終模型。
圖1提供了溫度對死亡率影響的總體情況,顯示了與參考值21?C(總體最低死亡率點(diǎn))相比,沿著溫度和滯后的相對風(fēng)險(RR)的三維圖。該圖顯示了熱的非常強(qiáng)烈和直接的影響,并表明對極熱的溫度有更多的滯后影響。寒冷溫度的最大影響大約在滯后2-3年達(dá)到。
盡管3-D圖是總結(jié)兩個維度的總體關(guān)系的有用工具,但不能包括估計的不確定性。為了對這種關(guān)系進(jìn)行更具體的評估,我們可以繪制特定溫度或滯后期的影響。圖2顯示了特定滯后期(0、5、15和28)的溫度和特定溫度(-10.8、-2.4、26.5和31.3?C)的滯后期的RR,大約對應(yīng)于溫度分布的第0.1、5、95和99.9百分位數(shù)(稱為中度和極端寒冷和炎熱)。溫度的總體影響,將分析中考慮的30天滯后期的貢獻(xiàn)相加,包括在下面。溫度-死亡率關(guān)系似乎隨著滯后期而變化,滯后期0和5的最低死亡率點(diǎn)不同(左上角的前兩個圖)。該圖證實(shí),如果與中度高溫相比,極端高溫的影響更為滯后,其顯著風(fēng)險分別持續(xù)10天和3天(右上角第三和第四張圖)。盡管如此,只有極端高溫表明可能存在收獲效應(yīng),在滯后15天后開始。相對于21?C的總體估計RR是1.24(95%CI:1.13-1.36)和1.07(95%CI:1.03-1.11),對于極端和中度高溫來說。寒冷的溫度顯示出完全不同的模式,中度寒冷的影響持續(xù)到滯后25天(右上角的前兩個圖)。此外,寒冷的影響似乎趨于平緩,中度寒冷的總體RR略高,為1.30(95%CI:1.20-1.40),而極度寒冷的RR為1.20(95%CI:1.04-1.39)(如下圖)。
plot(ns.pred,"overall"

為了將這一DLNM與更簡單的替代方法進(jìn)行比較,對滯后0-1和滯后0-30的移動平均和溫度空間的相同樣條函數(shù)的模型進(jìn)行了擬合。前者對高溫的影響提供了類似的估計,但顯示低溫的影響較弱,中度寒冷的估計RR為1.06(95%CI:1.03-1.09)。這一差異可能是由于低估了,因?yàn)榈蜏禺a(chǎn)生的影響持續(xù)時間超過2天。相反,滯后0-30的移動平均模型對寒冷的影響相似,但對高溫的估計較低,對中度和極端高溫的RR分別為1.01(95%CI:0.97-1.04)和1.06(95%CI:0.97-1.17)??紤]到滯后期內(nèi)的每一個先前的暴露都被假定為對每一天的影響提供了相同的貢獻(xiàn),平均31天的估計值可能會造成一些偏差,這是可信的。上述標(biāo)準(zhǔn)表明DLNM的擬合效果更好,如果與滯后0-1和0-30移動平均模型相比,QAIC的差異為571和517,QBIC為468和445。
已經(jīng)進(jìn)行了敏感性分析,以評估模型選擇的影響。特別是,我們評估了與改變用于指定交叉基函數(shù)(沿兩個維度)以及季節(jié)性和長期趨勢部分的df有關(guān)的估計總體效果的變化。增加溫度空間的結(jié)數(shù),產(chǎn)生的平滑曲線要少得多,可能是由于過度擬合,而在滯后維度上選擇不同的樣條,沒有明顯的變化。使用更多的df來控制季節(jié)和長期趨勢并不影響估計值,除了在非常低的溫度下溫度-死亡率曲線有不太明顯的下降。
此外,對滯后和特定溫度曲線的檢查顯示,當(dāng)增加季節(jié)性控制時,在長滯后期的負(fù)面效應(yīng)完全消失了。因?yàn)榫哂休^長滯后期的模型的效果對季節(jié)性成分更敏感。
討論
在本文中,我們描述了DLNMs的類別,可以用來模擬同時顯示非線性依賴和滯后效應(yīng)的因素的影響。DLNM在概念上是簡單的,但又足夠靈活,允許有廣泛的模型,包括以前使用的簡單模型和更復(fù)雜的新變體。
概念上的簡單性允許構(gòu)建一個R包來擬合這種廣泛的模型。這種豐富的選擇(基礎(chǔ)類型、結(jié)的數(shù)量和位置、最大滯后)所強(qiáng)調(diào)的一個困難是,可以用什么標(biāo)準(zhǔn)來選擇替代品。
在上面的例子中,我們用信息標(biāo)準(zhǔn)來指導(dǎo)結(jié)點(diǎn)數(shù)量的選擇,但在選擇基類型和最大滯后時,我們用的是先驗(yàn)論證。以前從流行病學(xué)的角度對DLNM的選擇進(jìn)行了討論,由于對什么是 "最佳 "模型沒有共識,敏感性分析特別重要,可以評估關(guān)鍵結(jié)論對模型選擇的依賴性。
DLMN的范圍很廣,這有助于實(shí)現(xiàn)這一目標(biāo)。回歸診斷,如殘差和部分自相關(guān)圖,也可能有幫助。此外,我們已經(jīng)討論了DLNM的選擇,假設(shè)它集中在感興趣的變量上(在我們的例子中是溫度)。還有一個協(xié)變量的模型選擇問題,其中的一些部分也可能是DLNMs。
這個問題,有時被稱為調(diào)整的不確定性。同樣,在什么方法是最佳的問題上還沒有形成共識,對模型選擇的這一部分的敏感性分析也很重要。
參考文獻(xiàn)
Zanobetti A, Schwartz J, Samoli E, Gryparis A, Touloumi G, Atkinson R, Le Tertre A, Bobros J, Celko M, Goren A, Forsberg B, Michelozzi P, Rabczenko D, Aranguez Ruiz E, Katsouyanni K. The temporal pattern of mortality responses to air pollution: a multicity assessment of mortality displacement. Epidemiology 2002; 13(1):87--93.
Braga AL, Zanobetti A, Schwartz J. The time course of weather-related deaths. Epidemiology 2001; 12(6):662--667.
Schwartz J. Is there harvesting in the association of airborne particles with daily deaths and hospital admissions? Epidemiology 2001; 12(1):55--61.
本文摘選?《?R語言分布滯后非線性模型(DLNM)空氣污染研究溫度對死亡率影響建模應(yīng)用?》?,點(diǎn)擊“閱讀原文”獲取全文完整資料。
點(diǎn)擊標(biāo)題查閱往期內(nèi)容
R語言分布滯后線性和非線性模型(DLM和DLNM)建模
分布滯后線性和非線性模型(DLNM)分析空氣污染(臭氧)、溫度對死亡率時間序列數(shù)據(jù)的影響
R語言中的分布滯后非線性模型DLNM與發(fā)病率和空氣污染示例
【視頻】R語言中的分布滯后非線性模型(DLNM)與發(fā)病率,死亡率和空氣污染示例
R語言分布滯后線性和非線性模型(DLNM)分析空氣污染(臭氧)、溫度對死亡率時間序列數(shù)據(jù)的影響
R語言分布滯后線性和非線性模型(DLMs和DLNMs)分析時間序列數(shù)據(jù)
R語言分布滯后非線性模型(DLNM)空氣污染研究溫度對死亡率影響建模應(yīng)用R語言分布滯后非線性模型(DLNM)研究發(fā)病率,死亡率和空氣污染示例
R語言分布滯后線性和非線性模型(DLM和DLNM)建模
R語言廣義相加模型 (GAMs)分析預(yù)測CO2時間序列數(shù)據(jù)
Python | ARIMA時間序列模型預(yù)測航空公司的乘客數(shù)量
R語言中生存分析模型的時間依賴性ROC曲線可視化
R語言ARIMA,SARIMA預(yù)測道路交通流量時間序列分析:季節(jié)性、周期性
ARIMA模型預(yù)測CO2濃度時間序列-python實(shí)現(xiàn)
R語言基于遞歸神經(jīng)網(wǎng)絡(luò)RNN的溫度時間序列預(yù)測
R語言用多元ARMA,GARCH ,EWMA, ETS,隨機(jī)波動率SV模型對金融時間序列數(shù)據(jù)建模
R語言神經(jīng)網(wǎng)絡(luò)模型預(yù)測車輛數(shù)量時間序列
卡爾曼濾波器:用R語言中的KFAS建模時間序列
在Python中使用LSTM和PyTorch進(jìn)行時間序列預(yù)測
R語言從經(jīng)濟(jì)時間序列中用HP濾波器,小波濾波和經(jīng)驗(yàn)?zāi)B(tài)分解等提取周期性成分分析
使用PYTHON中KERAS的LSTM遞歸神經(jīng)網(wǎng)絡(luò)進(jìn)行時間序列預(yù)測
Python中的ARIMA模型、SARIMA模型和SARIMAX模型對時間序列預(yù)測
R語言k-Shape時間序列聚類方法對股票價格時間序列聚類
R語言多元Copula GARCH 模型時間序列預(yù)測