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的臭氧進行單次暴露後，未來每一天的風險增加。

點擊標題查閱往期內容

R語言分布滯後線性和非線性模型（DLMs和DLNMs）分析時間序列數據

左右滑動查看更多

01

02

03

04

或者，可以繪製總體效果，該總體效果是通過使用參數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(pre.temp, "overall", ylim = c(0.5, 2.5), ci = "n", lwd = 1.5,

+ main = "Overall effect")

R> lines(pretemp2, "overall", col = 3, lty = 2, lwd = 2)

R> lines(pretemp3, "overall", col = 4, lty = 4, lwd = 2

+ p, 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）中的交叉基函數，例如，包括線性閾值，層次，多項式和樣條變換。

交叉基和參數估計的分離提供了多個優點。首先，如示例中所示，可以通過交叉基函數轉換多個顯示滯後效果的變量，並將其包含在模型中。其次，可以使用標準回歸命令進行估計，並使用默認的診斷工具和相關函數集。更重要的是，此實現提供了一個開放平台，可以在其中實現使用不同回歸命令指定的其他模型，來幫助在其他情況下或研究設計中開發方法。

本文摘選 《 R語言分布滯後線性和非線性模型（DLNM）分析空氣污染（臭氧）、溫度對死亡率時間序列數據的影響 》 ，點擊「閱讀原文」獲取全文完整資料。

點擊標題查閱往期內容

R語言中的分布滯後非線性模型DLNM與發病率和空氣污染示例

【視頻】R語言中的分布滯後非線性模型（DLNM）與發病率，死亡率和空氣污染示例

R語言分布滯後線性和非線性模型（DLNM）分析空氣污染（臭氧）、溫度對死亡率時間序列數據的影響

R語言分布滯後線性和非線性模型（DLMs和DLNMs）分析時間序列數據

R語言分布滯後非線性模型（DLNM）空氣污染研究溫度對死亡率影響建模應用R語言分布滯後非線性模型（DLNM）研究發病率，死亡率和空氣污染示例

R語言分布滯後線性和非線性模型（DLM和DLNM）建模

R語言廣義相加模型 (GAMs)分析預測CO2時間序列數據

Python | ARIMA時間序列模型預測航空公司的乘客數量

R語言中生存分析模型的時間依賴性ROC曲線可視化

R語言ARIMA，SARIMA預測道路交通流量時間序列分析:季節性、周期性

ARIMA模型預測CO2濃度時間序列-python實現

R語言基於遞歸神經網絡RNN的溫度時間序列預測

R語言用多元ARMA,GARCH ,EWMA, ETS,隨機波動率SV模型對金融時間序列數據建模

R語言神經網絡模型預測車輛數量時間序列

卡爾曼濾波器：用R語言中的KFAS建模時間序列

在Python中使用LSTM和PyTorch進行時間序列預測

R語言從經濟時間序列中用HP濾波器，小波濾波和經驗模態分解等提取周期性成分分析

使用PYTHON中KERAS的LSTM遞歸神經網絡進行時間序列預測

Python中的ARIMA模型、SARIMA模型和SARIMAX模型對時間序列預測

R語言k-Shape時間序列聚類方法對股票價格時間序列聚類

R語言多元Copula GARCH 模型時間序列預測