## Loading required package: nlme
## This is mgcv 1.8-15. For overview type 'help("mgcv-package")'.
Model 1
Mod.Train1 <- gam(flux ~ s(year,bs="cc",k=23)
+s(year,bs="cc",k=23,by=species,m=1)
+species+s(lat,bs="cr")
+s(long,bs="cr")
+indicator, data=Train.data)
summary(Mod.Train1)
##
## Family: gaussian
## Link function: identity
##
## Formula:
## flux ~ s(year, bs = "cc", k = 23) + s(year, bs = "cc", k = 23,
## by = species, m = 1) + species + s(lat, bs = "cr") + s(long,
## bs = "cr") + indicator
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.24141 0.09784 -12.688 < 2e-16 ***
## speciesCSH 0.05682 0.13395 0.424 0.67141
## speciesDBF -0.59941 0.10663 -5.621 1.92e-08 ***
## speciesEBF -1.57711 0.23677 -6.661 2.79e-11 ***
## speciesENF -0.54870 0.10676 -5.140 2.78e-07 ***
## speciesGRA 0.52446 0.10491 4.999 5.81e-07 ***
## speciesMF 0.36561 0.12332 2.965 0.00303 **
## speciesOSH -0.03325 0.13159 -0.253 0.80051
## speciesSAV 0.47162 0.19736 2.390 0.01687 *
## speciesWET 0.21335 0.14591 1.462 0.14369
## speciesWSA 0.24694 0.12745 1.938 0.05269 .
## indicatorW 1.44492 0.03140 46.022 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(year) 5.022e+00 21.000 1.080 2.72e-05 ***
## s(year):speciesCRO 1.777e+01 19.000 4.795 2.91e-12 ***
## s(year):speciesCSH 1.723e+00 17.000 0.393 0.01025 *
## s(year):speciesDBF 8.354e+00 21.000 1.249 0.00016 ***
## s(year):speciesEBF 2.201e-03 10.000 0.000 1.00000
## s(year):speciesENF 1.233e-07 20.000 0.000 0.28593
## s(year):speciesGRA 6.823e+00 18.000 0.795 0.02105 *
## s(year):speciesMF 2.211e+00 20.000 1.255 1.68e-07 ***
## s(year):speciesOSH 1.489e+00 14.000 0.642 0.00107 **
## s(year):speciesSAV 8.127e-09 13.000 0.000 1.00000
## s(year):speciesWET 4.814e+00 20.000 1.285 9.33e-06 ***
## s(year):speciesWSA 2.259e-08 17.000 0.000 0.48048
## s(lat) 7.998e+00 8.585 137.736 < 2e-16 ***
## s(long) 8.712e+00 8.958 121.644 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.205 Deviance explained = 20.8%
## GCV = 4.339 Scale est. = 4.3223 n = 19891
plot(residuals(Mod.Train1),pch=20,cex=0.3)
abline(h=0,col=2)
Prediction of Model 1
plot(Test.data[,2],pch=20,cex=0.1,ylab="CO2 Flux")
lines(predict(Mod.Train1,Test.data[,-c(1,2)]),pch=20,col=2)
legend("bottomleft",legend = c("True","Pred"),col=1:2,lty=2:1)
abline(v=cumsum(table(as.character(Test.data$site))),col="lightblue",lty="dashed")
Model 2
Mod.Train2 <- gam(flux ~ s(year,bs="cc",k=23)
+s(year,bs="cc",k=23,by=species,m=1)
+species+s(lat,bs="cr")
+s(long,bs="cr")
+indicator
+s(airtemp,bs="cr",k=23),data=Train.data)
summary(Mod.Train2)
##
## Family: gaussian
## Link function: identity
##
## Formula:
## flux ~ s(year, bs = "cc", k = 23) + s(year, bs = "cc", k = 23,
## by = species, m = 1) + species + s(lat, bs = "cr") + s(long,
## bs = "cr") + indicator + s(airtemp, bs = "cr", k = 23)
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.48439 0.06109 -7.929 2.33e-15 ***
## speciesCSH -0.18915 0.10389 -1.821 0.068661 .
## speciesDBF -0.74868 0.07154 -10.466 < 2e-16 ***
## speciesEBF -1.99212 0.33102 -6.018 1.80e-09 ***
## speciesENF -0.98470 0.07197 -13.681 < 2e-16 ***
## speciesGRA 0.53037 0.06873 7.717 1.25e-14 ***
## speciesMF 0.02664 0.09346 0.285 0.775643
## speciesOSH -0.11432 0.10551 -1.084 0.278587
## speciesSAV 0.70270 0.17445 4.028 5.65e-05 ***
## speciesWET 0.09814 0.11941 0.822 0.411169
## speciesWSA 0.35556 0.10091 3.524 0.000427 ***
## indicatorW -0.25854 0.04646 -5.565 2.65e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(year) 5.329e+00 21.000 1.085 1.68e-06 ***
## s(year):speciesCRO 1.668e+01 19.000 4.904 7.12e-14 ***
## s(year):speciesCSH 2.065e-07 17.000 0.000 0.789965
## s(year):speciesDBF 9.764e+00 21.000 2.546 3.72e-11 ***
## s(year):speciesEBF 3.372e+00 9.000 0.677 0.106624
## s(year):speciesENF 3.448e+00 20.000 0.743 0.000178 ***
## s(year):speciesGRA 5.897e+00 18.000 1.007 0.000690 ***
## s(year):speciesMF 3.549e-02 19.000 0.003 0.114119
## s(year):speciesOSH 3.310e+00 14.000 3.864 7.86e-15 ***
## s(year):speciesSAV 5.062e-07 8.000 0.000 1.000000
## s(year):speciesWET 5.978e+00 20.000 0.616 0.025940 *
## s(year):speciesWSA 1.923e+00 14.000 0.416 0.023443 *
## s(lat) 8.393e+00 8.830 120.564 < 2e-16 ***
## s(long) 8.844e+00 8.987 153.075 < 2e-16 ***
## s(airtemp) 1.250e+01 14.777 232.728 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.325 Deviance explained = 32.8%
## GCV = 3.712 Scale est. = 3.6936 n = 19624
plot(residuals(Mod.Train2),pch=20,cex=0.3)
abline(h=0,col=2)
Prediction of Model 2
plot(Test.data[,2],pch=20,cex=0.1,ylab="CO2 Flux")
lines(predict(Mod.Train2,Test.data[,-c(1,2)]),pch=20,col=2)
legend("bottomleft",legend = c("True","Pred"),col=1:2,lty=2:1)
abline(v=cumsum(table(as.character(Test.data$site))),col="lightblue",lty="dashed")
