library(readxl)
obs1 <- read_excel("C:/Users/wsand/Dropbox/Centro Oriente/Mayo/obs1 (6).xlsx")
library(DT)
DT::datatable(obs1)
library(dlnm)
library(splines)
library(MASS)
library(ggplot2)
library(quantmod)
obs1$fecha=as.Date(obs1$fecha)
obs1$tipodia=as.factor(obs1$tipodia)
##################MODELO CON DLMN#######################
argvar1_pm25 <- list(fun="lin")
argvar1_pm25
## $fun
## [1] "lin"
templin <- list(fun="lin")
HRlin <- list(fun="lin")
### sin rezago####
maxlag <-0
cb1_pm25 <- crossbasis(obs1$PM2.5,maxlag,argvar1_pm25)
templ <- crossbasis(obs1$temperatura,maxlag,templin)
HRl <- crossbasis(obs1$HR,maxlag,HRlin)
cb1_pm25=crossbasis(x=obs1$PM2.5, lag = 7,arglag =list(fun="ns"))
temp1 <- crossbasis(obs1$temperatura, lag=14, argvar=list("ns"))
# arglag=list(fun="strata",breaks=1))
summary(temp1)
## CROSSBASIS FUNCTIONS
## observations: 1826
## range: 12.21 to 17.43333
## lag period: 0 14
## total df: 1
##
## BASIS FOR VAR:
## fun: ns
## knots:
## intercept: FALSE
## Boundary.knots: 12.21 17.43333
##
## BASIS FOR LAG:
## fun: strata
## df: 1
## ref: 1
## intercept: TRUE
m1 <- glm(respira ~ temp1 + diasem + tipodia + ns(fecha,df=round(6*length(fecha)/365.25)),
data=obs1, family=quasipoisson)
varcen <- round(mean(obs1$temperatura),0)
#pred.o3 <- crosspred(basis.o3, model, at=c(0:65,40.3,50.3))
pred.temp <- crosspred(temp1, m1, by=1,cen = varcen)
pred.temp
## $predvar
## [1] 13 14 15 16 17
##
## $cen
## [1] 15
##
## $lag
## [1] 0 14
##
## $bylag
## [1] 1
##
## $coefficients
## temp1
## 0.03246167
##
## $vcov
## temp1
## temp1 0.0001558197
##
## $matfit
## lag0 lag1 lag2 lag3 lag4
## 13 -0.009965759 -0.009965759 -0.009965759 -0.009965759 -0.009965759
## 14 -0.004982879 -0.004982879 -0.004982879 -0.004982879 -0.004982879
## 15 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
## 16 0.004982879 0.004982879 0.004982879 0.004982879 0.004982879
## 17 0.009965759 0.009965759 0.009965759 0.009965759 0.009965759
## lag5 lag6 lag7 lag8 lag9
## 13 -0.009965759 -0.009965759 -0.009965759 -0.009965759 -0.009965759
## 14 -0.004982879 -0.004982879 -0.004982879 -0.004982879 -0.004982879
## 15 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
## 16 0.004982879 0.004982879 0.004982879 0.004982879 0.004982879
## 17 0.009965759 0.009965759 0.009965759 0.009965759 0.009965759
## lag10 lag11 lag12 lag13 lag14
## 13 -0.009965759 -0.009965759 -0.009965759 -0.009965759 -0.009965759
## 14 -0.004982879 -0.004982879 -0.004982879 -0.004982879 -0.004982879
## 15 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
## 16 0.004982879 0.004982879 0.004982879 0.004982879 0.004982879
## 17 0.009965759 0.009965759 0.009965759 0.009965759 0.009965759
##
## $matse
## lag0 lag1 lag2 lag3 lag4 lag5
## 13 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222
## 14 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111
## 15 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
## 16 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111
## 17 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222
## lag6 lag7 lag8 lag9 lag10 lag11
## 13 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222
## 14 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111
## 15 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
## 16 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111 0.001916111
## 17 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222 0.003832222
## lag12 lag13 lag14
## 13 0.003832222 0.003832222 0.003832222
## 14 0.001916111 0.001916111 0.001916111
## 15 0.000000000 0.000000000 0.000000000
## 16 0.001916111 0.001916111 0.001916111
## 17 0.003832222 0.003832222 0.003832222
##
## $allfit
## 13 14 15 16 17
## -0.14948638 -0.07474319 0.00000000 0.07474319 0.14948638
##
## $allse
## 13 14 15 16 17
## 0.05748333 0.02874167 0.00000000 0.02874167 0.05748333
##
## $matRRfit
## lag0 lag1 lag2 lag3 lag4 lag5 lag6
## 13 0.9900837 0.9900837 0.9900837 0.9900837 0.9900837 0.9900837 0.9900837
## 14 0.9950295 0.9950295 0.9950295 0.9950295 0.9950295 0.9950295 0.9950295
## 15 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
## 16 1.0049953 1.0049953 1.0049953 1.0049953 1.0049953 1.0049953 1.0049953
## 17 1.0100156 1.0100156 1.0100156 1.0100156 1.0100156 1.0100156 1.0100156
## lag7 lag8 lag9 lag10 lag11 lag12 lag13
## 13 0.9900837 0.9900837 0.9900837 0.9900837 0.9900837 0.9900837 0.9900837
## 14 0.9950295 0.9950295 0.9950295 0.9950295 0.9950295 0.9950295 0.9950295
## 15 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
## 16 1.0049953 1.0049953 1.0049953 1.0049953 1.0049953 1.0049953 1.0049953
## 17 1.0100156 1.0100156 1.0100156 1.0100156 1.0100156 1.0100156 1.0100156
## lag14
## 13 0.9900837
## 14 0.9950295
## 15 1.0000000
## 16 1.0049953
## 17 1.0100156
##
## $matRRlow
## lag0 lag1 lag2 lag3 lag4 lag5 lag6
## 13 0.9826751 0.9826751 0.9826751 0.9826751 0.9826751 0.9826751 0.9826751
## 14 0.9912997 0.9912997 0.9912997 0.9912997 0.9912997 0.9912997 0.9912997
## 15 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
## 16 1.0012281 1.0012281 1.0012281 1.0012281 1.0012281 1.0012281 1.0012281
## 17 1.0024578 1.0024578 1.0024578 1.0024578 1.0024578 1.0024578 1.0024578
## lag7 lag8 lag9 lag10 lag11 lag12 lag13
## 13 0.9826751 0.9826751 0.9826751 0.9826751 0.9826751 0.9826751 0.9826751
## 14 0.9912997 0.9912997 0.9912997 0.9912997 0.9912997 0.9912997 0.9912997
## 15 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
## 16 1.0012281 1.0012281 1.0012281 1.0012281 1.0012281 1.0012281 1.0012281
## 17 1.0024578 1.0024578 1.0024578 1.0024578 1.0024578 1.0024578 1.0024578
## lag14
## 13 0.9826751
## 14 0.9912997
## 15 1.0000000
## 16 1.0012281
## 17 1.0024578
##
## $matRRhigh
## lag0 lag1 lag2 lag3 lag4 lag5 lag6
## 13 0.9975483 0.9975483 0.9975483 0.9975483 0.9975483 0.9975483 0.9975483
## 14 0.9987734 0.9987734 0.9987734 0.9987734 0.9987734 0.9987734 0.9987734
## 15 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
## 16 1.0087767 1.0087767 1.0087767 1.0087767 1.0087767 1.0087767 1.0087767
## 17 1.0176304 1.0176304 1.0176304 1.0176304 1.0176304 1.0176304 1.0176304
## lag7 lag8 lag9 lag10 lag11 lag12 lag13
## 13 0.9975483 0.9975483 0.9975483 0.9975483 0.9975483 0.9975483 0.9975483
## 14 0.9987734 0.9987734 0.9987734 0.9987734 0.9987734 0.9987734 0.9987734
## 15 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000
## 16 1.0087767 1.0087767 1.0087767 1.0087767 1.0087767 1.0087767 1.0087767
## 17 1.0176304 1.0176304 1.0176304 1.0176304 1.0176304 1.0176304 1.0176304
## lag14
## 13 0.9975483
## 14 0.9987734
## 15 1.0000000
## 16 1.0087767
## 17 1.0176304
##
## $allRRfit
## 13 14 15 16 17
## 0.8611502 0.9279818 1.0000000 1.0776074 1.1612377
##
## $allRRlow
## 13 14 15 16 17
## 0.7693943 0.8771513 1.0000000 1.0185811 1.0375074
##
## $allRRhigh
## 13 14 15 16 17
## 0.9638485 0.9817579 1.0000000 1.1400542 1.2997236
##
## $ci.level
## [1] 0.95
##
## $model.class
## [1] "glm" "lm"
##
## $model.link
## [1] "log"
##
## attr(,"class")
## [1] "crosspred"
plot(pred.temp, "slices", type="p", pch=19, cex=1.5, var=15, ci="bars",
ylab="RR",main="Lag-specific effects")
plot(pred.temp, "overall", ci="lines", ylim=c(0.95,1.25), lwd=2, col=4,
xlab="Temperatura", ylab="RR", main="Overall effect")
plot(pred.temp, xlab="Temperature",zlab="RR", main="3D graph")
plot(pred.temp, "contour", plot.title=title(xlab="Temperature",
ylab="Lag", main="Contour graph"), key.title=title("RR"))
plot(pred.temp, "slices", var=14, ci="n",ylim=c(0.99,1.005), lwd=1.5)
for(i in 1:2) lines(pred.temp, "slices", var=c(14,16)[i], col=i+2, lwd=1.5)
legend("topright",paste("Temperature =",c(14,16,17)), col=2:4, lwd=1.5)
basis.o3 <- crossbasis(chicagoNMMAPS$o3, vartype="hthr", varknots=40.3,
lagtype="strata", lagknots=c(2,6), maxlag=10)
basis.temp <- crossbasis(chicagoNMMAPS$temp, vartype="bs", vardegree=3,
vardf=6, cenvalue=25, lagdf=5, maxlag=30)