# CREATE THE CROSS-BASIS: DOUBLE THRESHOLD AND NATURAL SPLINE
cb4 <- crossbasis(chicagoNMMAPS$temp, lag=30, argvar=list(type="dthr",
knots=c(10,25)), arglag=list(df=5))
# RUN THE MODEL AND GET THE PREDICTION FOR TEMPERATURE
library(splines)
model4 <- glm(death ~ cb4 + ns(time, 7*14) + dow,
family=quasipoisson(), chicagoNMMAPS)
pred4 <- crosspred(cb4, model4, by=1)
# REDUCE TO OVERALL ASSOCIATION
redall <- crossreduce(cb4, model4)
summary(redall)
# REDUCE TO LAG-SPECIFIC ASSOCIATION FOR LAG 5
redlag <- crossreduce(cb4, model4, type="lag", value=5)
# REDUCE TO PREDICTOR-SPECIFIC ASSOCIATION AT VALUE 33
redvar <- crossreduce(cb4, model4, type="var", value=33)
# NUMBER OF PARAMETERS OF THE ORIGINAL MODEL
length(coef(pred4))
# REDUCED NUMBER OF PARAMETERS FOR OVERALL AND LAG-SPECIFIC ASSOCIATIONS
length(coef(redall)) ; length(coef(redlag))
# REDUCED NUMBER OF PARAMETERS FOR PREDICTOR-SPECIFIC ASSOCIATIONS
length(coef(redvar))
# TEST: IDENTICAL FIT BETWEEN ORIGINAL AND REDUCED FIT
plot(pred4, "overall", xlab="Temperature", ylab="RR",
ylim=c(0.8,1.6), main="Overall effects")
lines(redall, ci="lines",col=4,lty=2)
legend("top",c("Original","Reduced"),col=c(2,4),lty=1:2,ins=0.1)
# RECONSTRUCT THE FIT IN TERMS OF ONE-DIMENSIONAL BASIS
b4 <- onebasis(0:30,knots=attributes(cb4)$arglag$knots,int=TRUE,cen=FALSE)
pred4b <- crosspred(b4,coef=coef(redvar),vcov=vcov(redvar),model.link="log",by=1)
# TEST: IDENTICAL FIT BETWEEN ORIGINAL, REDUCED AND RE-CONSTRUCTED
plot(pred4, "slices", var=33, ylab="RR", ylim=c(0.9,1.2),
main="Predictor-specific effects at 33C")
lines(redvar, ci="lines", col=4, lty=2)
points(pred4b, col=1, pch=19, cex=0.6)
legend("top",c("Original","Reduced","Reconstructed"),col=c(2,4,1),lty=c(1:2,NA),
pch=c(NA,NA,19),pt.cex=0.6,ins=0.1)
Run the code above in your browser using DataLab