Chad Milando, PhD, MS

cmilando@bu.edu

Interested in similar spreadsheets? see here:

# R code to reproduce this spreadsheet

https://chadmilando.com/teaching/

library(dlnm)

DLNM in Excel, with no macros (Solver used for Poisson)

library(splines)

Recommended view % is 75%

library(tidyverse)

library(lubridate)

https://onlinelibrary.wiley.com/doi/full/10.1002/sim.5963

https://github.com/gasparrini/dlnm

df <- dlnm::chicagoNMMAPS

df2 <- df %>% filter(month %in% 5:9, year == 1990)

df2$dow_fct <- paste0('x', wday(df2$date, week_start = 1))

Step 1: getting model coefficients

Chicago data sample

k50 <- quantile(df2$temp, probs = 0.5)

bs example

an example spreadsheet for the bs() function in R

using bs() on the temperature variable

cb1.temp <- dlnm::crossbasis(

creating a polynomial basis of the 30-lag

df2$temp, argvar=list(knots = k50, fun = "bs", degree = 3, intercept = F),

offset basisvar

offsetting the temperature basis to match the 30day lag

lag = 30, arglag = list(fun="poly", degree = 2, intercept = F))

crossbasis

taking row products of the two bases to create the cross basis

model3 <- glm(death ~ cb1.temp + dow_fct, family=poisson(), df2)

Poisson predict

getting model coefficients via unconditional poisson regression

Step 2: using model coefficients to make the graphs

pred3.temp <- crosspred(cb1.temp,

using the bs() from above, convert the center point to basis form

model3, bylag = 1, cumul = T, at = seq(6.5, 30, by = 0.5), cen = 20.5)

for the newdata you want to pred on, make basis, center on ^

combine varvec bs with varlag bs, which is poly^ in this case

plot(pred3.temp)

crosspred output

take different slices to view RR by lag or by temp overall

# for temp = 10, what is the RR by lag

plot(pred3.temp, var = 10, cumul = T,

The major first hurdle was replicating the cross-basis function

ylab="Cumulative RR", col = 2)

* basically, you use a crossbasis when you want to:

# what is the overall relationship

(1) model the relationship between exposure and outcome (non-linearly, independent of time)

plot(pred3.temp, 'overall', col = 2,

so, exposure = 50 C-> probabiliy of outcome 0.12, exposure = 40 C -> P(outcome) = 0.4

ylab="Cumulative RR", xlab="Temperature")

(2) model the impact of time (indepedent of exposure level)

your cumulative exposure today is 1 * todays exposure + 0.5 * yesterday's exposure

using cross basis allows you to model both together. you multiply the basis of each together

Note: I ignore the group component for simplicity, so this just assumes 1 strata etc.

but otherwise it’s the same, lags are just put inside each group

On the choice of centering value in crosspred

(from Gasparrini 2013)

Therefore the cross-basis in (7) should always be defined without an intercept in the basis functions for
x. Also, these basis functions can be centered on a specific exposure value x0, which will represent the
reference for the risk summaries computed by (8)–(10).

The centering value for crosspred doesn't have a statistical interpretation, other than as a reference

Just choose a comparison value that make the interpretation make sense

(mean daily max temperature)

Questions I'm still working on answering

> what is the rationale behind the formula for SE in varvec