時間序列之外的應用

分布滯後模型首先是在很久以前的計量經濟時間序列分析中提出的[Almon，1965]，然後在環境流行病學Schwartz [2000]的時間序列數據中重新提出。DLNM的擴展是由Armstrong [2006]構想的。Gasparrini等人對時間序列數據的建模框架進行了重新評估。[2010]。有趣的是，已經在不同的研究領域中提出了這種暴露-滯後-反應關聯的模型。一般的想法是通過特定函數加權過去的暴露，這些函數的參數由數據估算。在癌症流行病學[Hauptmann等，2000；Langholz等，1999；Richardson，2009；Thomas，1983；Vacek，1997]和藥物流行病學[Abrahamowicz等]中，說明了類似於DLM的線性-暴露-反應關係模型。

基本函數

指定標準暴露反應和滯後反應關係的基本函數，例如多項式，分層或閾值函數。例如，樣條線由推薦的包樣條線中包含的函數ns（）和bs（）指定。多項式是通過函數poly（）獲得的。這是一個簡單向量的轉換示例：

poly(1:5,degree=3)

1 2 3

[1,] 0.2 0.04 0.008

[2,] 0.4 0.16 0.064

[3,] 0.6 0.36 0.216

[4,] 0.8 0.64 0.512

[5,] 1.0 1.00 1.000

attr(,"degree")

[1] 3

attr(,"scale")

[1] 5

attr(,"intercept")

[1] FALSE

attr(,"class")

[1] "poly" "matrix"

第一個未命名的參數x指定要轉換的向量，而參數度設置多項式的度。定義分層函數是通過strata（）指定的。

strata(1:5,breaks=c(2,4))[,]

1 2

[1,] 0 0

[2,] 1 0

[3,] 1 0

[4,] 0 1

[5,] 0 1

結果是帶有附加類別「層」的基礎矩陣。轉換是定義對比的虛擬參數化。參數break定義了層的右開放區間的下邊界。

閾值函數通過thr（）指定。一個例子：

thr(1:5,thr.value=3,side="d")[,]

1 2

[1,] 2 0

[2,] 1 0

[3,] 0 0

[4,] 0 1

[5,] 0 2

結果是具有附加類別「 thr」的基礎矩陣。參數thr.value定義一個帶有一個或兩個閾值的向量，而side用於指定高（「 h」，默認值），低（「 l」）或雙精度（「 d」）閾值參數化。

點擊標題查閱往期內容

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

左右滑動查看更多

01

02

03

04

基本轉換

此函數代表以dlnm為單位進行基本轉換的主要函數，適用於指定暴露-反應和滯後-反應關係。它的作用是應用選定的轉換並以適用於其他函數（例如crossbasis（）和crosspred（））的格式生成基本矩陣。以下示例複製了該部分中顯示的多項式變換：

onebasis(1:5,fun="poly",degree=3)

b1 b2 b3

[1,] 0.2 0.04 0.008

[2,] 0.4 0.16 0.064

[3,] 0.6 0.36 0.216

[4,] 0.8 0.64 0.512

[5,] 1.0 1.00 1.000

attr(,"fun")

[1] "poly"

attr(,"degree")

[1] 3

attr(,"scale")

[1] 5

attr(,"intercept")

[1] FALSE

attr(,"class")

[1] "onebasis" "matrix"

attr(,"range")

[1] 1 5

結果是帶有附加類「 onebasis」的基礎矩陣。同樣，第一個未命名參數x指定要轉換的向量，而第二個參數fun將字符轉換定義為應用轉換而調用的函數的名稱。具體來說，基本矩陣包括fun和range屬性，以及定義轉換的被調用函數的參數。如前所述，onebasis（）還可以根據特定要求調用用戶定義的函數。一個簡單的例子：

> mylog <- function(x) log(x)

> onebasis(1:5,"mylog")

b1

[1,] 0.0000000

[2,] 0.6931472

[3,] 1.0986123

[4,] 1.3862944

[5,] 1.6094379

attr(,"fun")

[1] "mylog"

attr(,"range")

[1] 1 5

attr(,"class")

[1] "onebasis" "matrix"

交叉基

這是dlnm軟體包中的主要函數。它在內部調用onebasis（）來生成暴露-反應和滯後-反應關係的基矩陣，並通過特殊的張量積將它們組合起來，以創建交叉基，該交叉基在模型中同時指定了暴露-滯後-反應關聯性。它的第一個參數x的類定義如何解釋數據。可以使用第二個變量lag修改滯後期。

作為一個簡單的示例，我模擬了2-5個滯後期內3個對象的暴露歷史矩陣：它們中的每一個都將傳遞給onebasis（）來分別構建暴露-反應和滯後-反應關係的矩陣。僅用於時間序列數據的附加參數組定義了被視為單獨無關序列的觀察組，例如在季節性分析中可能有用。作為一個簡單的示例，我模擬了2-5個滯後期內3個對象的暴露歷史矩陣：它們中的每一個都將傳遞給onebasis（）來分別構建暴露-反應和滯後-反應關係的矩陣。作為一個簡單的示例，我模擬了2-5個滯後期內3個對象的暴露歷史矩陣：

> hist

lag2 lag3 lag4 lag5

sub1 1 3 5 6

sub2 7 8 9 4

sub3 10 2 11 12

然後，我應用交叉基參數化，將二次多項式作為暴露反應函數，並將分層函數2-3和4-5定義為滯後反應函數的分層函數：

lag=c(2,5),argvar=list(fun="poly",degree=2),

arglag=list(fun="strata",breaks=4))[,]

v1.l1 v1.l2 v2.l1 v2.l2

sub1 1.250000 0.9166667 0.4930556 0.4236111

sub2 2.333333 1.0833333 1.4583333 0.6736111

sub3 2.916667 1.9166667 2.5625000 1.8402778

該函數返回「 crossbasis」類的矩陣對象。它首先使用argvar和arglag列表中的參數調用onebasis（），以建立暴露反應空間和滯後反應空間的矩陣基礎。在另一個示例中，我將crossbasis（）應用於數據集中的變量temp，該數據集表示1987-2000年期間日平均溫度序列：

> summary(cb)

CROSSBASIS FUNCTIONS

observations: 5114

range: -26.66667 to 33.33333

lag period: 0 30

total df: 10

BASIS FOR VAR:

fun: thr

thr.value: 10 20

side: d

intercept: FALSE

BASIS FOR LAG:

fun: ns

knots: 1 4 12

intercept: TRUE

Boundary.knots: 0 30

此處，將暴露反應建模為閾值為10和20的雙閾值函數。滯後時間設置為0到30。滯後反應函數留給默認的自然三次樣條（fun =「 ns」），其滯後值為1、4和12。

預測

crossbasis（）生成的交叉基矩陣需要包含在回歸模型公式中才能擬合模型。關聯通過函數crosspred（）進行匯總，該函數針對默認值或用戶直接選擇的預測值和滯後值的組合的網格進行預測。例如，我使用創建的交叉基矩陣cb，使用數據集時間序列數據來研究溫度與心血管疾病死亡率之間的關聯。首先，我將一個簡單的線性模型與模型公式中包含的交叉基矩陣擬合。然後，我通過使用cross-basis和回歸模型對象作為前兩個參數調用crosspred（）來獲得預測：

crosspred(cb,model,at=-20:30)

結果是「 crosspred」類的列表對象，其中的存儲預測和有關模型的其他信息，例如係數和與交叉基參數有關的關聯（協）方差矩陣的一部分。可以為特定的預測器-滯後組合選擇預測的網格。例如，我提取溫度為-10°C且滯後5的預測和置信區間，然後提取25°C的整體累積預測：

> pred$allfit["25"]

25

1.108262

第一個結果表明，在給定的一天中，-20°C的溫度會在五天後導致0.95例心血管死亡的增加，或者在給定的一天中，溫度為-6攝氏度時，心血管死亡的數目增加0.95。其他類型的預測可以通過crosspred（）獲得。特別是，如果模型連結等於log或logit，則將自動返回取冪的預測。如果參數cum設置為TRUE，則是累積預測的矩陣cum。

crosspred（）的另一種用法是預測特定的暴露歷史記錄集的影響。這可以通過輸入暴露歷史矩陣作為參數來實現。例如，我們可以從擬合模型中預測出，在過去10天暴露於30°C和在滯後期的其餘時間暴露於22°C之後，心血管死亡的總體累積增加：如果參數cum設置為TRUE，則包括增量累積預測的矩陣cum，並將其存儲在組件cum中。crosspred（）的另一種用法是預測特定的暴露歷史記錄集的影響。這可以通過輸入暴露歷史矩陣作為參數來實現。例如，我們可以從擬合模型中預測出，在過去10天暴露於30°C和在滯後期的其餘時間暴露於22°C之後，心血管死亡的總體累積增加：

> crosspred(cb,model,at=histpred)$allfit

1

5.934992

dlnm軟體包的主要優點之一是，用戶可以使用標準回歸函數執行DLNM，只需在模型公式中包括交叉基矩陣即可。函數crosspred（）自動處理來自回歸函數lm（）和glm（），gam（）（軟體包mgcv），coxph（）的模型。

降維

DLNM的擬合度可以降低到預測變量或滯後的一個維度，僅以諸如總累積暴露反應表達。該計算通過函數crossreduce（）進行，該函數的工作原理與crosspred（）非常相似。前兩個自變量base和model指定交叉基矩陣和需要對其執行計算的模型對象。減少的類型由類型定義，帶有選項「 overall」-「 lag」-「 var」，用於匯總總體累積暴露反應，滯後特異性暴露反應或預測變量特異性滯後反應。

繪圖

一維或二維關聯的解釋通過圖形表示來輔助。通過方法函數plot（），lines（）和points（）為類「 crosspred」和「 crossreduce」提供高級和低級繪圖功能。例如，我使用對象pred中的預測。plot（）方法可以通過參數ptype為「 crosspred」對象生成不同類型的圖。具體來說，它會生成整個二維暴露-滯後-反應關聯的圖形。二維關聯可以繪製為3-D或等高線圖，例如：

> plot(pred,ptype="3d",main="3D plot"

可以通過選擇不同的ptype獲得定義的關聯的摘要。

> plot(pred,"overall"

在這種情況下，方法函數plot（）在內部調用函數plot.default（），如上面的示例所示，可以將其他特定參數添加到函數調用中。通過設置ptype =「 slices」，可以將滯後特異性和預測因子特異性關聯分別繪製為暴露-反應和滯後-反應關係，因為它們是在3-D曲面中沿特定維度切割的切片。例如：

> plot(pred,"slices",lag=5

這兩個圖分別代表了滯後5的暴露反應和特定於25°C溫度的滯後反應。參數lag和var指定必須分別繪製lag和特定於預測變量的關聯的值。

點擊文末 「閱讀原文」

獲取全文完整代碼數據資料。

本文選自《R語言分布滯後非線性模型（DLNM）研究發病率，死亡率和空氣污染示例》。

點擊標題查閱往期內容

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

分布滯後線性和非線性模型（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 模型時間序列預測