Make a table for mvmeta output for the meta-regression section in the manuscript. And figures if applicable.
# lpath <- "/home/ad.abt.local/layc/R/library" # assign library for user w/o admin priveleges on ACE3
# .libPaths(lpath)
library(here) # convenience for filepaths dynamic
## here() starts at C:/Users/LayC/Projects/cira3-mortality
here <- here::here # to avoid overwrite
if(grepl("C://Users/", here())){
lpath = paste("C://Users/", Sys.info( )["login"], "/R.lib", sep = "")
.libPaths(lpath)}
#library(chron)
library(tictoc) # for system time
library(readxl) # for reading in excel
library(readr) # for csv
library(survival) # from original vodonos code. Loaded for compatibility
library(splines) # for splines
library(MASS) # for mvmeta, likely
#library(nlme)
library(mgcv) # from original vodonos code. Loaded for compatibility
## Loading required package: nlme
## This is mgcv 1.8-31. For overview type 'help("mgcv-package")'.
library(dlnm) # from original vodonos code. Loaded for compatibility
## This is dlnm 2.3.9. For details: help(dlnm) and vignette('dlnmOverview').
library(tsModel) # from original vodonos code. Loaded for compatibility
## Time Series Modeling for Air Pollution and Health (0.6)
library(mvmeta) # from original vodonos code. Loaded for compatibility
## This is mvmeta 1.0.3. For an overview type: help('mvmeta-package').
library(mixmeta) # UPDATED MVMETA package with better random specification syntax
## This is mixmeta 1.0.7. For an overview type: help('mixmeta-package').
library(data.table) # fread for faster read (read.table takes >15 min)
library(tidyverse)
## -- Attaching packages ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- tidyverse 1.3.0 --
## v ggplot2 3.2.1 v dplyr 0.8.4
## v tibble 2.1.3 v stringr 1.4.0
## v tidyr 1.0.2 v forcats 0.4.0
## v purrr 0.3.3
## -- Conflicts -------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- tidyverse_conflicts() --
## x dplyr::between() masks data.table::between()
## x dplyr::collapse() masks nlme::collapse()
## x dplyr::filter() masks stats::filter()
## x dplyr::first() masks data.table::first()
## x dplyr::lag() masks stats::lag()
## x dplyr::last() masks data.table::last()
## x dplyr::select() masks MASS::select()
## x purrr::transpose() masks data.table::transpose()
library(scales)
##
## Attaching package: 'scales'
## The following object is masked from 'package:purrr':
##
## discard
## The following object is masked from 'package:readr':
##
## col_factor
select <- dplyr::select # to re-assign select to specify dplyr as is more commonly used than raster and MASS functions in script below.
# Temporary fix due to changes in nest/unnest syntax (deprecated arguments)
nest <- nest_legacy
unnest <- unnest_legacy
#library(RcppRoll) # rolling averages from vodonos? Likey not needed
library(purrr) # for iterative processing.
library(gridExtra) # for grid.arrange
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
source(here::here("/code/functions/01_cluster_load_cp.r"))# load cluster functions
rm("cluster_load_cp"); rm("cluster_load_mv") # remove extra functions used for meta analysis
source(here::here("code/functions/03_lag_calc.R"))
source(here::here("code/functions/06_rr_calc.R"))
source(here::here("code/functions/metareg/000_ftab_mvmeta.R"))
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/"), recursive = TRUE)}
ymat_get <- function(tp = 2, dframe = cluster_2_parms, ymat_name = "ymat_1"){
# function to extract ymat from dframe
cs = dframe %>% filter(time_period == tp ) %>%
pluck("parms" ,1, ymat_name ) #cs[!rowSums( is.na(cs)),]
}
Slist_get <- function(tp =2, dframe = cluster_2_parms, Slist_name = "Slist_1"){
# functio to get slist from dframe
cs = dframe %>% filter(time_period == tp ) %>%
pluck("parms" , 1, Slist_name )
}
cluster_fits <- seq(1,9) %>%
map_dfr(cluster_load,
nsim = 1000,
out_folder = "output/vodonos_out/2019-06-13_mvcl_4_time/",
time_periods = c(1,2,3,4)
)
## [1] 1
## [1] 2
## [1] 3
## [1] 4
## [1] 5
## [1] 6
## [1] 7
## [1] 8
## [1] 9
Look at some of the mvmeta output from the cluster specific meta-smoothing models
cluster_fits %>% pluck("mv_out", 1, "ymat_1" ) # ymat of all cities
## spl1 spl2 spl3 spl4 spl5
## allentown 0.223362177 0.198929909 0.302113223 0.37267640 0.64065875
## annandale 1.592067947 1.102159027 1.398238321 1.41230645 1.82157873
## atlantic city 0.077264947 -0.016694798 -0.002552391 0.07002456 0.15926094
## baltimore -0.066703992 -0.039379716 -0.007150994 0.07048845 0.21609647
## barnstable 0.199860389 0.089040430 0.117558485 0.29077362 1.46001543
## bergen -0.026331468 0.048608556 0.113890584 0.19651228 0.33974409
## boston 0.029881558 0.056607030 0.048018196 0.22324269 0.64192642
## carlisle -0.200703738 -0.053618663 0.004072488 -0.05877374 0.10789410
## dover 0.794429794 0.797074626 0.841675695 0.67448954 0.85258381
## elizabeth 0.191301233 0.210001140 0.279130764 0.34673129 0.46687332
## gettysburg -0.068329262 0.027458093 -0.051612631 0.16795089 0.50103005
## harrisburg 0.365601390 0.292608332 0.379424365 0.28845810 0.40461245
## hartford 0.088259436 0.075639814 0.125614587 0.08145279 0.38426543
## jersy city 0.005533399 0.033438390 0.056586700 0.15170943 0.52639732
## lancast -0.058229108 0.005920053 0.037262813 0.13721602 0.24517864
## marlboro -0.197370303 -0.061205581 -0.104535652 0.13385565 0.48303328
## middles 0.381306724 0.306319053 0.345796044 0.33511468 0.55882694
## monmouth -0.076003285 0.005408390 -0.041782831 0.07157760 0.10381079
## nassau -0.165873679 -0.122767045 -0.122366542 -0.01200808 0.01801486
## newark 0.031663413 0.050541883 0.076475016 0.16638079 0.46967362
## newburgh 0.023559619 0.026587694 0.056810913 0.15977627 0.79748784
## newhaven -0.015451022 -0.022268768 -0.034749404 0.11310331 0.48635866
## newlond -0.010379816 0.001305222 -0.016449195 0.20137046 0.30852881
## newyork 0.112892307 0.141771793 0.002081847 0.24007191 0.84848270
## philly 0.080621509 0.083594127 0.095244735 0.20285817 0.36605050
## plymouth -0.027413608 0.012622576 0.054512148 0.13879260 0.39423790
## reading 0.044844118 0.103148246 0.119939733 0.35204257 0.94036644
## richmond -0.056445197 -0.037093835 -0.008144682 0.05923435 0.10311170
## rockville 0.295507678 0.333940778 0.376909419 0.47534103 0.27696733
## essex 0.055368475 0.034925302 0.061694096 0.14939194 0.78228920
## springfieldma -0.076787507 -0.021978272 0.022010201 0.13568199 0.58115590
## stamford -0.040104706 0.040993705 0.018454321 0.04744804 0.42693811
## trenton -0.081119457 -0.037277705 -0.022051622 0.06720402 0.23277161
## wdc 0.111375764 0.157426972 0.188700594 0.29083444 0.31139213
## wilmington -0.080215580 0.006804586 0.085853112 0.08159035 0.26840127
## york NA NA NA NA NA
cluster_fits %>% pluck("mv_out", 1, "Slist_1", "allentown") # specific vcov
## btmean0b1 btmean0b2 btmean0b3 btmean0b4 btmean0b5
## btmean0b1 0.01627959 0.01156877 0.01318139 0.01247519 0.01307927
## btmean0b2 0.01156877 0.01132552 0.01125925 0.01130633 0.01120745
## btmean0b3 0.01318139 0.01125925 0.01468001 0.01303859 0.01457783
## btmean0b4 0.01247519 0.01130633 0.01303859 0.01588292 0.01249334
## btmean0b5 0.01307927 0.01120745 0.01457783 0.01249334 0.04486401
Pulls in data used for mvmeta fit. Refits models fitted by Vodonos for QA and to identify missing city-specific data. Visual check against vodonos output in meta-regression_Cluster1.md. In cluster 1, winter and summer mean temp have r of 0.96. Avoided including both in models. Population variables are also tightly linked (0.91), but were left in for consistency with Vodonos results.
# Year 1: had data corrections, may change if data obtained from Harvard.
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster1/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster1/"), recursive = T)}
source(here("code/scripts/cluster_1_meta_reg_data_correct.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_con totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ----------- ----------- ----------- ----------- ----------- ------------
## summer_mean 1.0000000 -0.1604216 0.9557826 0.0873320 -0.0200521 -0.0839655 0.5021941 0.1770252 -0.0200521 -0.0839655
## summer_std -0.1604216 1.0000000 -0.3379536 0.4091625 -0.0633439 -0.0544270 -0.0625313 -0.0778238 -0.0633439 -0.0544270
## winter_mean 0.9557826 -0.3379536 1.0000000 -0.0875072 -0.0284958 -0.0888064 0.4374341 0.1207286 -0.0284958 -0.0888064
## winter_std 0.0873320 0.4091625 -0.0875072 1.0000000 -0.0184793 -0.1527144 -0.1937093 -0.3628200 -0.0184793 -0.1527144
## totalpop -0.0200521 -0.0633439 -0.0284958 -0.0184793 1.0000000 0.9114431 0.1262878 0.0712323 1.0000000 0.9114431
## popover65 -0.0839655 -0.0544270 -0.0888064 -0.1527144 0.9114431 1.0000000 0.1490393 0.2495949 0.9114431 1.0000000
## air_cond 0.5021941 -0.0625313 0.4374341 -0.1937093 0.1262878 0.1490393 1.0000000 0.7371860 0.1262878 0.1490393
## time_con 0.1770252 -0.0778238 0.1207286 -0.3628200 0.0712323 0.2495949 0.7371860 1.0000000 0.0712323 0.2495949
## totpop1 -0.0200521 -0.0633439 -0.0284958 -0.0184793 1.0000000 0.9114431 0.1262878 0.0712323 1.0000000 0.9114431
## popover65_1 -0.0839655 -0.0544270 -0.0888064 -0.1527144 0.9114431 1.0000000 0.1490393 0.2495949 0.9114431 1.0000000
meta_reg_tab_cluster1 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster1/cluster_1_meta_reg.csv"))
meta_reg_tab_cluster1
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 mod1_lag0 1 lag 0
## 4 mod2_lag0 1 lag 0
## 5 mod_ac_int_lag0 1 lag 0
## 6 modac_lag0 1 lag 0
## 7 int_ma15 1 MA15
## 8 time_ma15 1 MA15
## 9 mod1_ma15 1 MA15
## 10 mod2_ma15 1 MA15
## 11 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + summer_mean + time1 + time2 + time4 + totpop1 + popover65_1
## 4 (Intercept) + time1 + time2 + time4 + totpop1 + popover65_1
## 5 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 6 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean
## 7 (Intercept)
## 8 (Intercept) + time1 + time2 + time4
## 9 (Intercept) + time1 + time2 + time4 + winter_mean + totpop1 + popover65_1
## 10 (Intercept) + time1 + time2 + time4 + totpop1 + popover65_1
## 11 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 1196.9373 705 0.000000e+00 41.09967 -1629.921 -1538.616
## 2 (Intercept) 963.2552 690 2.370126e-11 28.36789 -1682.015 -1522.231
## 3 (Intercept) 851.1906 675 4.298623e-06 20.69931 -1728.296 -1500.032
## 4 (Intercept) 949.4542 680 3.185829e-11 28.37990 -1670.272 -1464.835
## 5 (Intercept) 794.8092 675 9.560312e-04 15.07396 -1770.727 -1542.464
## 6 (Intercept) 849.7423 680 9.130792e-06 19.97574 -1738.166 -1532.729
## 7 1108.1146 846 2.802937e-09 23.65411 -2096.715 -1968.530
## 8 (Intercept) 1008.5797 828 1.545127e-05 17.90436 -2114.704 -1901.063
## 9 (Intercept) 903.9894 810 1.169714e-02 10.39718 -2138.207 -1839.109
## 10 (Intercept) 949.8648 816 7.772441e-04 14.09303 -2113.013 -1842.401
## 11 (Intercept) 912.3799 810 6.944024e-03 11.22120 -2129.805 -1830.707
Note that in cluster 2, total and 65+ populations had pearson r of 0.98; removed the total population from models.
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster2/"))){dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster2/"))}
source(here("code/scripts/cluster_2_meta_reg.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_con totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ---------- ----------- ---------- ----------- ---------- ------------
## summer_mean 1.0000000 0.3114383 0.4975557 0.3773406 0.1164384 0.1242576 0.5140942 0.1031048 0.1164384 0.1242576
## summer_std 0.3114383 1.0000000 -0.4136432 0.6228724 0.0271125 0.0026857 0.1209993 -0.0993570 0.0271125 0.0026857
## winter_mean 0.4975557 -0.4136432 1.0000000 -0.5157627 0.0482538 0.0659956 0.2700248 0.2087518 0.0482538 0.0659956
## winter_std 0.3773406 0.6228724 -0.5157627 1.0000000 0.0136148 -0.0097478 0.0714751 -0.2538976 0.0136148 -0.0097478
## totalpop 0.1164384 0.0271125 0.0482538 0.0136148 1.0000000 0.9882329 0.0893624 0.0197545 1.0000000 0.9882329
## popover65 0.1242576 0.0026857 0.0659956 -0.0097478 0.9882329 1.0000000 0.1263842 0.0793820 0.9882329 1.0000000
## air_cond 0.5140942 0.1209993 0.2700248 0.0714751 0.0893624 0.1263842 1.0000000 0.7729171 0.0893624 0.1263842
## time_con 0.1031048 -0.0993570 0.2087518 -0.2538976 0.0197545 0.0793820 0.7729171 1.0000000 0.0197545 0.0793820
## totpop1 0.1164384 0.0271125 0.0482538 0.0136148 1.0000000 0.9882329 0.0893624 0.0197545 1.0000000 0.9882329
## popover65_1 0.1242576 0.0026857 0.0659956 -0.0097478 0.9882329 1.0000000 0.1263842 0.0793820 0.9882329 1.0000000
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## ref Q df_q Q_p I2 AIC BIC
## 1 1222.8547 828 0.000000e+00 32.28959 -1417.068 -1350.934
## 2 (Intercept) 1054.2999 816 2.942610e-08 22.60267 -1474.335 -1351.515
## 4 (Intercept) 909.3318 804 5.598028e-03 11.58343 -1550.106 -1370.600
## model_name cluster lag
## 1 int_ma15 1 MA15
## 2 time_ma15 1 MA15
## 3 mod1_ma15 1 MA15
## 5 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + popover65_1
## 5 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 1108.1146 846 2.802937e-09 23.65411 -2096.715 -1968.530
## 2 (Intercept) 1008.5797 828 1.545127e-05 17.90436 -2114.704 -1901.063
## 3 (Intercept) 913.3065 810 6.543176e-03 11.31126 -2130.044 -1830.946
## 5 (Intercept) 912.3799 810 6.944024e-03 11.22120 -2129.805 -1830.707
meta_reg_tab_cluster2 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster2/cluster_2_metareg.csv"))
meta_reg_tab_cluster2
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 modac_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## 5 mod1_lag0 1 lag 0
## 6 mod2_lag0 1 lag 0
## 7 int_ma15 1 MA15
## 8 time_ma15 1 MA15
## 9 mod1_ma15 1 MA15
## 10 mod2_ma15 1 MA15
## 11 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 5 (Intercept) + summer_mean + time1 + time2 + time4 + popover65_1
## 6 (Intercept) + time1 + time2 + time4 + popover65_1
## 7 (Intercept)
## 8 (Intercept) + time1 + time2 + time4
## 9 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + popover65_1
## 10 (Intercept) + time1 + time2 + time4 + popover65_1
## 11 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 1222.8547 828 0.000000e+00 32.28959 -1417.068 -1350.934
## 2 (Intercept) 1054.2999 816 2.942610e-08 22.60267 -1474.335 -1351.515
## 3 (Intercept) 946.5510 808 5.125053e-04 14.63746 -1526.392 -1365.782
## 4 (Intercept) 909.3318 804 5.598028e-03 11.58343 -1550.106 -1370.600
## 5 (Intercept) 952.6739 808 3.126156e-04 15.18609 -1526.245 -1365.635
## 6 (Intercept) 1040.2301 812 8.851209e-08 21.94035 -1473.302 -1331.587
## 7 1108.1146 846 2.802937e-09 23.65411 -2096.715 -1968.530
## 8 (Intercept) 1008.5797 828 1.545127e-05 17.90436 -2114.704 -1901.063
## 9 (Intercept) 913.3065 810 6.543176e-03 11.31126 -2130.044 -1830.946
## 10 (Intercept) 968.5495 822 2.955231e-04 15.13082 -2110.686 -1868.559
## 11 (Intercept) 912.3799 810 6.944024e-03 11.22120 -2129.805 -1830.707
This cluster had a frame-shift issue similar to Cluster 1. St. Louis was removed from the meta-regression, but the data had Wichita removed instead. Removed both cities.
# Year 1 missing the last city in the time period [1:26,]
# Vodonos lag0 - mod 1 summer_mean+as.factor(time)+totpop1+popover65_1
# Vodonos lag04 mod1 summer_mean+winter_mean+as.factor(time)+totpop1+popover65_1
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster3/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster3/"), recursive = T)}
source(here("code/scripts/cluster_3_meta_reg.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_rev totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ----------- ----------- ----------- ----------- ----------- ------------
## summer_mean 1.0000000 0.3332959 0.8261337 0.1633614 -0.3777460 -0.3941669 0.6116407 0.1508590 -0.3777460 -0.3941669
## summer_std 0.3332959 1.0000000 -0.1411438 0.4852300 -0.2931388 -0.3464019 0.1852764 -0.1525042 -0.2931388 -0.3464019
## winter_mean 0.8261337 -0.1411438 1.0000000 -0.2544158 -0.2619416 -0.2377128 0.5307216 0.2852607 -0.2619416 -0.2377128
## winter_std 0.1633614 0.4852300 -0.2544158 1.0000000 -0.2754959 -0.3960677 -0.1888755 -0.4961840 -0.2754959 -0.3960677
## totalpop -0.3777460 -0.2931388 -0.2619416 -0.2754959 1.0000000 0.9504208 -0.2083284 -0.0036072 1.0000000 0.9504208
## popover65 -0.3941669 -0.3464019 -0.2377128 -0.3960677 0.9504208 1.0000000 -0.1380769 0.1325274 0.9504208 1.0000000
## air_cond 0.6116407 0.1852764 0.5307216 -0.1888755 -0.2083284 -0.1380769 1.0000000 0.7708572 -0.2083284 -0.1380769
## time_rev 0.1508590 -0.1525042 0.2852607 -0.4961840 -0.0036072 0.1325274 0.7708572 1.0000000 -0.0036072 0.1325274
## totpop1 -0.3777460 -0.2931388 -0.2619416 -0.2754959 1.0000000 0.9504208 -0.2083284 -0.0036072 1.0000000 0.9504208
## popover65_1 -0.3941669 -0.3464019 -0.2377128 -0.3960677 0.9504208 1.0000000 -0.1380769 0.1325274 0.9504208 1.0000000
meta_reg_tab_cluster_3 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster3/cluster_3_metareg.csv"))
meta_reg_tab_cluster_3
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 mod1_lag0 1 lag 0
## 4 mod2_lag0 1 lag 0
## 5 mod_ac_int_lag0 1 lag 0
## 6 modac_lag0 1 lag 0
## 7 int_ma15 1 MA15
## 8 time_ma15 1 MA15
## 9 mod1_ma15 1 MA15
## 10 mod2_ma15 1 MA15
## 11 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + summer_mean + time1 + time2 + time4 + popover65_1
## 4 (Intercept) + time1 + time2 + time4 + popover65_1
## 5 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 6 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean
## 7 (Intercept)
## 8 (Intercept) + time1 + time2 + time4
## 9 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + popover65_1
## 10 (Intercept) + time1 + time2 + time4 + popover65_1
## 11 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 565.4460 372 3.589863e-10 34.21122 -799.3973 -744.3831
## 2 (Intercept) 474.6819 360 4.546195e-05 24.15975 -815.6369 -713.4676
## 3 (Intercept) 410.6064 352 1.693118e-02 14.27314 -823.3892 -689.7832
## 4 (Intercept) 417.4161 356 1.369693e-02 14.71340 -826.4175 -708.5298
## 5 (Intercept) 414.9380 348 7.854455e-03 16.13205 -815.0832 -665.7588
## 6 (Intercept) 422.5535 352 5.794518e-03 16.69694 -818.0682 -684.4622
## 7 391.7059 372 2.313483e-01 5.03078 -846.9324 -791.9182
## 8 (Intercept) 357.5769 360 5.261694e-01 0.00000 -850.9162 -748.7469
## 9 (Intercept) 344.8606 348 5.374712e-01 0.00000 -838.3713 -689.0469
## 10 (Intercept) 351.7207 356 5.540928e-01 0.00000 -848.2381 -730.3504
## 11 (Intercept) 342.1529 348 5.783125e-01 0.00000 -840.9580 -691.6336
#vodonos lag0 mod1 summer_mean + winter_mean + as.factor(time)
#vodonos lag05 mod1 summer_mean+winter_mean+ as.factor(time)+totpop1+popover65_1+totpop1+air_cond
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster4/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster4/"), recursive = T)}
source(here("code/scripts/cluster_4_meta_reg.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_con totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ---------- ---------- ----------- ----------- ---------- ------------
## summer_mean 1.0000000 -0.3035588 0.8367802 -0.0342598 0.5439749 0.4941585 0.2955859 0.1241268 0.5439749 0.4941585
## summer_std -0.3035588 1.0000000 -0.7404784 0.7548768 0.1813262 0.2362585 0.0326036 -0.1698272 0.1813262 0.2362585
## winter_mean 0.8367802 -0.7404784 1.0000000 -0.4691200 0.2335434 0.1589356 0.1140392 0.1209720 0.2335434 0.1589356
## winter_std -0.0342598 0.7548768 -0.4691200 1.0000000 0.2438651 0.2501911 -0.1521836 -0.4196718 0.2438651 0.2501911
## totalpop 0.5439749 0.1813262 0.2335434 0.2438651 1.0000000 0.9696580 0.3350578 0.1500953 1.0000000 0.9696580
## popover65 0.4941585 0.2362585 0.1589356 0.2501911 0.9696580 1.0000000 0.4330929 0.2764142 0.9696580 1.0000000
## air_cond 0.2955859 0.0326036 0.1140392 -0.1521836 0.3350578 0.4330929 1.0000000 0.8187204 0.3350578 0.4330929
## time_con 0.1241268 -0.1698272 0.1209720 -0.4196718 0.1500953 0.2764142 0.8187204 1.0000000 0.1500953 0.2764142
## totpop1 0.5439749 0.1813262 0.2335434 0.2438651 1.0000000 0.9696580 0.3350578 0.1500953 1.0000000 0.9696580
## popover65_1 0.4941585 0.2362585 0.1589356 0.2501911 0.9696580 1.0000000 0.4330929 0.2764142 0.9696580 1.0000000
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 modac_lag0 1 lag 0
## mod ref
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4 (Intercept)
## 3 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean (Intercept)
## Q df_q Q_p I2 AIC BIC
## 1 577.1849 420 5.122327e-07 27.23302 -836.0124 -779.3161
## 2 500.7304 408 1.138837e-03 18.51902 -864.5751 -759.2820
## 3 477.3856 400 4.656503e-03 16.21030 -862.4814 -724.7905
## model_name cluster lag
## 1 int_ma15 1 MA15
## 2 time_ma15 1 MA15
## 3 mod1_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond + popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 512.5305 420 0.001319493 18.05366 -843.3958 -786.6995
## 2 (Intercept) 473.6110 408 0.013651237 13.85335 -855.0696 -749.7765
## 3 (Intercept) 445.5598 392 0.031753332 12.02078 -848.1147 -678.0259
meta_reg_tab_cluster_4 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster4/cluster_4_metareg.csv"))
meta_reg_tab_cluster_4
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 modac_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## 5 mod1_lag0 1 lag 0
## 6 mod2_lag0 1 lag 0
## 7 mod3_lag0 1 lag 0
## 8 int_ma15 1 MA15
## 9 time_ma15 1 MA15
## 10 mod1_ma15 1 MA15
## 11 mod2_ma15 1 MA15
## 12 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 5 (Intercept) + air_cond + summer_mean + time1 + time2 + time4 + totpop1 + popover65_1
## 6 (Intercept) + time1 + time2 + time4 + air_cond + totpop1 + popover65_1
## 7 (Intercept) + time1 + time2 + time4 + totpop1 + popover65_1
## 8 (Intercept)
## 9 (Intercept) + time1 + time2 + time4
## 10 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond + popover65_1
## 11 (Intercept) + time1 + time2 + time4 + air_cond + popover65_1
## 12 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 577.1849 420 5.122327e-07 27.23302 -836.0124 -779.3161
## 2 (Intercept) 500.7304 408 1.138837e-03 18.51902 -864.5751 -759.2820
## 3 (Intercept) 477.3856 400 4.656503e-03 16.21030 -862.4814 -724.7905
## 4 (Intercept) 473.9355 396 4.277516e-03 16.44433 -856.9782 -703.0883
## 5 (Intercept) 469.2884 392 4.387395e-03 16.46927 -852.3012 -682.2124
## 6 (Intercept) 477.1715 396 3.155248e-03 17.01096 -855.9226 -702.0327
## 7 (Intercept) 488.2034 400 1.651568e-03 18.06694 -858.8667 -721.1758
## 8 512.5305 420 1.319493e-03 18.05366 -843.3958 -786.6995
## 9 (Intercept) 473.6110 408 1.365124e-02 13.85335 -855.0696 -749.7765
## 10 (Intercept) 445.5598 392 3.175333e-02 12.02078 -848.1147 -678.0259
## 11 (Intercept) 460.2860 400 1.988300e-02 13.09751 -850.9217 -713.2308
## 12 (Intercept) 453.9418 396 2.331430e-02 12.76415 -848.2626 -694.3728
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster6/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster6/"), recursive = T)}
source(here("code/scripts/cluster_6_meta_reg.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_con totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ---------- ----------- ----------- ----------- ---------- ------------
## summer_mean 1.0000000 -0.2441263 0.6311149 -0.1060338 0.4714084 0.5206866 0.3145885 0.3328939 0.4714084 0.5206866
## summer_std -0.2441263 1.0000000 -0.8641236 0.8133605 0.0763300 0.0021514 0.0128317 -0.1562215 0.0763300 0.0021514
## winter_mean 0.6311149 -0.8641236 1.0000000 -0.6526725 0.1462948 0.2140265 0.0307875 0.1549982 0.1462948 0.2140265
## winter_std -0.1060338 0.8133605 -0.6526725 1.0000000 0.0434398 -0.0950771 -0.2803974 -0.4932011 0.0434398 -0.0950771
## totalpop 0.4714084 0.0763300 0.1462948 0.0434398 1.0000000 0.9623808 0.2219028 0.1542825 1.0000000 0.9623808
## popover65 0.5206866 0.0021514 0.2140265 -0.0950771 0.9623808 1.0000000 0.3116030 0.3067505 0.9623808 1.0000000
## air_cond 0.3145885 0.0128317 0.0307875 -0.2803974 0.2219028 0.3116030 1.0000000 0.8423861 0.2219028 0.3116030
## time_con 0.3328939 -0.1562215 0.1549982 -0.4932011 0.1542825 0.3067505 0.8423861 1.0000000 0.1542825 0.3067505
## totpop1 0.4714084 0.0763300 0.1462948 0.0434398 1.0000000 0.9623808 0.2219028 0.1542825 1.0000000 0.9623808
## popover65_1 0.5206866 0.0021514 0.2140265 -0.0950771 0.9623808 1.0000000 0.3116030 0.3067505 0.9623808 1.0000000
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 5 mod1_lag0 1 lag 0
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 5 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond
## ref Q df_q Q_p I2 AIC BIC
## 1 332.4330 295 0.06583941 11.260301 -553.8002 -479.7246
## 2 (Intercept) 300.9885 280 0.18574372 6.973187 -547.7063 -418.0739
## 5 (Intercept) 278.3906 265 0.27392769 4.810012 -539.5465 -354.3574
## model_name cluster lag
## 1 int_ma15 1 MA15
## 2 time_ma15 1 MA15
## 3 moda_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond + popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 283.7322 236 0.01810864 16.822977 -432.8219 -384.0930
## 2 (Intercept) 256.5406 224 0.06682955 12.684397 -424.5546 -334.0580
## 3 (Intercept) 227.7070 208 0.16609123 8.654551 -417.4109 -271.2241
meta_reg_tab_cluster_6 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster6/cluster_6_metareg.csv"))
meta_reg_tab_cluster_6
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 modac_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## 5 mod1_lag0 1 lag 0
## 6 mod2_lag0 1 lag 0
## 7 mod3_lag0 1 lag 0
## 8 int_ma15 1 MA15
## 9 time_ma15 1 MA15
## 10 moda_ma15 1 MA15
## 11 mod1_ma15 1 MA15
## 12 mod2_ma15 1 MA15
## 13 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 5 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond
## 6 (Intercept) + time1 + time2 + time4 + air_cond + totpop1 + popover65_1
## 7 (Intercept) + time1 + time2 + time4 + totpop1 + popover65_1
## 8 (Intercept)
## 9 (Intercept) + time1 + time2 + time4
## 10 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond + popover65_1
## 11 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond
## 12 (Intercept) + time1 + time2 + time4 + air_cond
## 13 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 332.4330 295 0.06583941 11.260301 -553.8002 -479.7246
## 2 (Intercept) 300.9885 280 0.18574372 6.973187 -547.7063 -418.0739
## 3 (Intercept) 291.5469 270 0.17562688 7.390533 -536.3533 -369.6831
## 4 (Intercept) 285.3658 265 0.18627170 7.136738 -532.7325 -347.5434
## 5 (Intercept) 278.3906 265 0.27392769 4.810012 -539.5465 -354.3574
## 6 (Intercept) 288.0457 265 0.15813082 8.000703 -528.0336 -342.8445
## 7 (Intercept) 295.3076 270 0.13890279 8.569909 -532.5744 -365.9041
## 8 283.7322 236 0.01810864 16.822977 -432.8219 -384.0930
## 9 (Intercept) 256.5406 224 0.06682955 12.684397 -424.5546 -334.0580
## 10 (Intercept) 227.7070 208 0.16609123 8.654551 -417.4109 -271.2241
## 11 (Intercept) 234.9262 212 0.13399160 9.758912 -419.9512 -287.6869
## 12 (Intercept) 244.6650 220 0.12179011 10.081136 -428.3016 -323.8824
## 13 (Intercept) 234.6059 212 0.13715607 9.635708 -418.2535 -285.9892
# vodonos mod1 lag0 = summer_mean+winter_mean+ as.factor(time)+air_cond
# vodonos mod2 lag0 = as.factor(time)+totpop1+popover65_1+air_cond
# vodonos mod1 lag05 = summer_mean+winter_mean+ as.factor(time)+air_cond
# vodonos lag 0 mod1 summer_mean+ as.factor(time)+air_cond
# vodonos lag05 mod1 summer_mean+winter_mean+ as.factor(time)+air_cond
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster7/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster7/"), recursive = T)}
source(here("code/scripts/cluster_7_meta_reg.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_con totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ----------- ----------- ----------- ----------- ----------- ------------
## summer_mean 1.0000000 0.2281936 0.3987825 0.1759789 0.1774208 0.1109428 0.1683549 0.3447052 0.1774208 0.1109428
## summer_std 0.2281936 1.0000000 -0.7658609 0.9370874 -0.3084551 -0.4413713 -0.2600212 -0.0402018 -0.3084551 -0.4413713
## winter_mean 0.3987825 -0.7658609 1.0000000 -0.7493665 0.3371723 0.4173951 0.2532085 0.1197161 0.3371723 0.4173951
## winter_std 0.1759789 0.9370874 -0.7493665 1.0000000 -0.3801667 -0.5208813 -0.4884183 -0.2407328 -0.3801667 -0.5208813
## totalpop 0.1774208 -0.3084551 0.3371723 -0.3801667 1.0000000 0.8906829 0.4449854 0.4609526 1.0000000 0.8906829
## popover65 0.1109428 -0.4413713 0.4173951 -0.5208813 0.8906829 1.0000000 0.4998701 0.4344477 0.8906829 1.0000000
## air_cond 0.1683549 -0.2600212 0.2532085 -0.4884183 0.4449854 0.4998701 1.0000000 0.8427727 0.4449854 0.4998701
## time_con 0.3447052 -0.0402018 0.1197161 -0.2407328 0.4609526 0.4344477 0.8427727 1.0000000 0.4609526 0.4344477
## totpop1 0.1774208 -0.3084551 0.3371723 -0.3801667 1.0000000 0.8906829 0.4449854 0.4609526 1.0000000 0.8906829
## popover65_1 0.1109428 -0.4413713 0.4173951 -0.5208813 0.8906829 1.0000000 0.4998701 0.4344477 0.8906829 1.0000000
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 5 mod1_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 5 (Intercept) + summer_mean + time1 + time2 + time4 + air_cond
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## ref Q df_q Q_p I2 AIC BIC
## 1 275.2446 236 0.04047038 14.258079 -565.9020 -517.1731
## 2 (Intercept) 242.8344 224 0.18476715 7.756062 -568.6033 -478.1067
## 5 (Intercept) 223.9754 216 0.34050737 3.560847 -570.4603 -452.1186
## 4 (Intercept) 219.7155 212 0.34364581 3.511590 -566.4482 -434.1839
## model_name cluster lag
## 1 int_ma15 1 MA15
## 2 time_ma15 1 MA15
## 4 mod1_ma15 1 MA15
## mod ref
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4 (Intercept)
## 4 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean (Intercept)
## Q df_q Q_p I2 AIC BIC
## 1 298.8246 236 0.003478418 21.02391 -538.2779 -489.5489
## 2 266.4864 224 0.027187196 15.94319 -535.6081 -445.1115
## 4 242.0625 216 0.107756183 10.76686 -539.9236 -421.5819
meta_reg_tab_cluster_7 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster7/cluster_7_metareg.csv"))
meta_reg_tab_cluster_7
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 mod_add_all_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## 5 mod1_lag0 1 lag 0
## 6 mod2_lag0 1 lag 0
## 7 mod3_lag0 1 lag 0
## 8 int_ma15 1 MA15
## 9 time_ma15 1 MA15
## 10 moda_ma15 1 MA15
## 11 mod1_ma15 1 MA15
## 12 mod2_ma15 1 MA15
## 13 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + air_cond + summer_mean + time1 + time2 + time4 + totpop1 + popover65_1
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 5 (Intercept) + summer_mean + time1 + time2 + time4 + air_cond
## 6 (Intercept) + time1 + time2 + time4 + air_cond
## 7 (Intercept) + time1 + time2 + time4 + totpop1 + popover65_1
## 8 (Intercept)
## 9 (Intercept) + time1 + time2 + time4
## 10 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + air_cond + totpop1 + popover65_1
## 11 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean
## 12 (Intercept) + time1 + time2 + time4 + air_cond
## 13 (Intercept) + time1 + time2 + time4 + air_cond + popover65_1 + air_cond:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 275.2446 236 0.040470384 14.258079 -565.9020 -517.1731
## 2 (Intercept) 242.8344 224 0.184767148 7.756062 -568.6033 -478.1067
## 3 (Intercept) 217.9350 208 0.304282858 4.558691 -560.4659 -414.2791
## 4 (Intercept) 219.7155 212 0.343645813 3.511590 -566.4482 -434.1839
## 5 (Intercept) 223.9754 216 0.340507369 3.560847 -570.4603 -452.1186
## 6 (Intercept) 237.3567 220 0.200952685 7.312500 -565.6956 -461.2764
## 7 (Intercept) 235.4792 216 0.173029211 8.272139 -559.9172 -441.5755
## 8 298.8246 236 0.003478418 21.023914 -538.2779 -489.5489
## 9 (Intercept) 266.4864 224 0.027187196 15.943193 -535.6081 -445.1115
## 10 (Intercept) 236.2389 204 0.060382169 13.646728 -520.8527 -360.7433
## 11 (Intercept) 242.0625 216 0.107756183 10.766862 -539.9236 -421.5819
## 12 (Intercept) 265.3413 220 0.019700280 17.087915 -528.4953 -424.0761
## 13 (Intercept) 250.4450 212 0.036189344 15.350672 -524.8139 -392.5496
# vodonos lag0 mod1 summer_mean+winter_mean+ as.factor(time)
# vodonos lag5 mod1summer_mean+winter_mean+as.factor(time)
if(!dir.exists(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster9/"))){ dir.create(here("output/y4-final-manuscript-output/2020-03-03-Manu-Meta-reg/Cluster9/"), recursive = T)}
source(here("code/scripts/cluster_9_meta_reg.r"))
## Parsed with column specification:
## cols(
## city = col_character(),
## summer_mean = col_double(),
## summer_std = col_double(),
## winter_mean = col_double(),
## winter_std = col_double(),
## time = col_double(),
## totalpop = col_double(),
## popover65 = col_double(),
## fips_sm = col_double(),
## air_cond = col_double()
## )
##
##
## summer_mean summer_std winter_mean winter_std totalpop popover65 air_cond time_con totpop1 popover65_1
## ------------ ------------ ----------- ------------ ----------- ----------- ----------- ----------- ----------- ----------- ------------
## summer_mean 1.0000000 0.1627083 0.5568355 0.3320857 0.2044237 0.0733911 0.2656498 0.1961129 0.2044237 0.0733911
## summer_std 0.1627083 1.0000000 -0.4086842 0.5305974 0.0577095 0.0030144 -0.0858308 0.1328253 0.0577095 0.0030144
## winter_mean 0.5568355 -0.4086842 1.0000000 -0.4888799 0.2722466 0.2014908 0.3210037 0.1699737 0.2722466 0.2014908
## winter_std 0.3320857 0.5305974 -0.4888799 1.0000000 -0.1127992 -0.1947994 -0.1107561 -0.1613095 -0.1127992 -0.1947994
## totalpop 0.2044237 0.0577095 0.2722466 -0.1127992 1.0000000 0.9342209 0.2383239 0.2978862 1.0000000 0.9342209
## popover65 0.0733911 0.0030144 0.2014908 -0.1947994 0.9342209 1.0000000 0.2363786 0.4136327 0.9342209 1.0000000
## air_cond 0.2656498 -0.0858308 0.3210037 -0.1107561 0.2383239 0.2363786 1.0000000 0.5835344 0.2383239 0.2363786
## time_con 0.1961129 0.1328253 0.1699737 -0.1613095 0.2978862 0.4136327 0.5835344 1.0000000 0.2978862 0.4136327
## totpop1 0.2044237 0.0577095 0.2722466 -0.1127992 1.0000000 0.9342209 0.2383239 0.2978862 1.0000000 0.9342209
## popover65_1 0.0733911 0.0030144 0.2014908 -0.1947994 0.9342209 1.0000000 0.2363786 0.4136327 0.9342209 1.0000000
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 mod_add_all_lag0 1 lag 0
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + air_cond + summer_mean + time1 + time2 + time4 + totpop1 + popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 198.6080 204 0.5933798 0.000000 -213.4489 -166.7234
## 2 (Intercept) 181.8878 192 0.6883580 0.000000 -202.5256 -115.7496
## 3 (Intercept) 217.9350 208 0.3042829 4.558691 -560.4659 -414.2791
## model_name cluster lag mod ref
## 1 int_ma15 1 MA15 (Intercept)
## 2 time_ma15 1 MA15 (Intercept) + time1 + time2 + time4 (Intercept)
## Q df_q Q_p I2 AIC BIC
## 1 171.8391 204 0.9506366 0 -212.5059 -165.7803
## 2 160.7221 192 0.9513013 0 -199.0393 -112.2633
meta_reg_tab_cluster_9 %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/Cluster9/cluster_9_metareg.csv"))
meta_reg_tab_cluster_9
## model_name cluster lag
## 1 int_lag0 1 lag 0
## 2 time_lag0 1 lag 0
## 3 mod_add_all_lag0 1 lag 0
## 4 modac_int_lag0 1 lag 0
## 5 modac_lag0 1 lag 0
## 6 mod1_lag0 1 lag 0
## 7 mod2_lag0 1 lag 0
## 8 mod3_lag0 1 lag 0
## 9 int_ma15 1 MA15
## 10 time_ma15 1 MA15
## 11 moda_ma15 1 MA15
## 12 mod1_ma15 1 MA15
## 13 mod2_ma15 1 MA15
## 14 mod3_ma15 1 MA15
## mod
## 1 (Intercept)
## 2 (Intercept) + time1 + time2 + time4
## 3 (Intercept) + air_cond + summer_mean + time1 + time2 + time4 + totpop1 + popover65_1
## 4 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + air_cond:summer_mean
## 5 (Intercept) + time1 + time2 + time4 + air_cond + summer_mean + winter_mean + totpop1 + popover65_1
## 6 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean
## 7 (Intercept) + time1 + time2 + time4 + air_cond
## 8 (Intercept) + time1 + time2 + time4 + totpop1 + popover65_1
## 9 (Intercept)
## 10 (Intercept) + time1 + time2 + time4
## 11 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean + totpop1 + popover65_1
## 12 (Intercept) + time1 + time2 + time4 + summer_mean + winter_mean
## 13 (Intercept) + time1 + time2 + time4 + air_cond
## 14 (Intercept) + time1 + time2 + time4 + winter_mean + popover65_1 + winter_mean:popover65_1
## ref Q df_q Q_p I2 AIC BIC
## 1 198.6080 204 0.5933798 0.000000 -213.4489 -166.72338
## 2 (Intercept) 181.8878 192 0.6883580 0.000000 -202.5256 -115.74957
## 3 (Intercept) 217.9350 208 0.3042829 4.558691 -560.4659 -414.27909
## 4 (Intercept) 285.3658 265 0.1862717 7.136738 -532.7325 -347.54338
## 5 (Intercept) 146.1928 172 0.9238431 0.000000 -196.2975 -42.77075
## 6 (Intercept) 157.4825 184 0.9221583 0.000000 -209.6957 -96.21943
## 7 (Intercept) 178.2531 188 0.6833469 0.000000 -197.7292 -97.60310
## 8 (Intercept) 167.5943 184 0.8015045 0.000000 -198.8747 -85.39838
## 9 171.8391 204 0.9506366 0.000000 -212.5059 -165.78034
## 10 (Intercept) 160.7221 192 0.9513013 0.000000 -199.0393 -112.26327
## 11 (Intercept) 147.5325 176 0.9420005 0.000000 -180.2066 -40.02997
## 12 (Intercept) 151.1997 184 0.9631286 0.000000 -192.5394 -79.06307
## 13 (Intercept) 158.4390 188 0.9426909 0.000000 -193.3001 -93.17392
## 14 (Intercept) 152.9035 180 0.9293812 0.000000 -182.8355 -56.00910
meta_reg_tab <- bind_rows(meta_reg_tab_cluster1, meta_reg_tab_cluster2, meta_reg_tab_cluster_3, meta_reg_tab_cluster_4,
meta_reg_tab_cluster_6, meta_reg_tab_cluster_7, meta_reg_tab_cluster_9)
meta_reg_tab %>% write_csv(here("output/y4-final-manuscript-output/2020-03-03-Manu-meta-reg/2020-03-03-Meta-Compiled.csv"))
meta_reg_tab %>% filter(cluster ==3 & lag == "lag 0")
## [1] model_name cluster lag mod ref Q
## [7] df_q Q_p I2 AIC BIC
## <0 rows> (or 0-length row.names)
print(cluster_load)
## function (c_cluster = 6, nsim = 1000, out_folder = "output/vodonos_out/2019-06-13_mvcl_4_time",
## time_periods = c(3, 4))
## {
## require(here)
## require(dplyr)
## require(tidyr)
## require(purrr)
## require(MASS)
## print(c_cluster)
## cluster_f <- tibble(cluster = c_cluster, time_period = time_periods,
## cp_fname = here::here(out_folder, gsub("X", c_cluster,
## paste0("mv_clX/mv_clX/cp.clusterX.years", time_period,
## ".Rdata"))), cp_out = map(cp_fname, .f = load_parms),
## mv_fname = here::here(out_folder, gsub("X", c_cluster,
## paste0("mv_clX/mv_clX/mv_clX_y", time_period, ".Rdata"))),
## mv_out = map(mv_fname, .f = load_parms), daily_arg = map(cp_out,
## ~pluck_all(., col = "cp_0", ns = nsim)), ma15_arg = map(cp_out,
## ~pluck_all(., col = "cp_05", ns = nsim)))
## }
## <bytecode: 0x00000000127c5778>
print(lag_calc)
## function (tmean_2lag = tmean_daily$doy_tmean0[[1]])
## {
## lagged_tmean <- tmean_2lag %>% mutate(tmeanL1 = ifelse(is.na(lag(tmean0)),
## tmean0, lag(tmean0)), tmeanL2 = ifelse(is.na(lag(tmean0,
## 2)), tmean0, lag(tmean0, 2)), tmeanL3 = ifelse(is.na(lag(tmean0,
## 3)), tmean0, lag(tmean0, 3)), tmeanL4 = ifelse(is.na(lag(tmean0,
## 4)), tmean0, lag(tmean0, 4)), tmean15 = (tmean0 + tmeanL1 +
## tmeanL2 + tmeanL3 + tmeanL4)/5) %>% select(-starts_with("tmeanL"))
## return(lagged_tmean)
## }
print(calc_rr_se)
## function (tmean = seq(-25, 42, 1), coef_mat = pluck(cl4_y3, "cp_out",
## 1, "cp_0", "coefficients"), vcov_mat = pluck(cl4_y3, "cp_out",
## 1, "cp_0", "vcov"), center_val = 15.6, degree_val = 2, knot_vec = pluck(cl4_y3,
## "cp_out", 1, "knots_1"), bound_vec = c(-25, 42))
## {
## require(dlnm)
## require(splines)
## bs_ce <- bs(center_val, degree = degree_val, knots = knot_vec,
## Bound = bound_vec, intercept = FALSE)
## bs_var <- suppressWarnings(bs(tmean, degree = degree_val,
## knots = knot_vec, Bound = bound_vec, intercept = FALSE))
## bs_vc <- scale(bs_var, center = bs_ce, scale = FALSE)
## rr_log <- as.numeric(bs_vc %*% coef_mat)
## rr_log_se <- sqrt(pmax(0, rowSums((bs_vc %*% vcov_mat) *
## bs_vc)))
## RR_an <- 1 - exp(-rr_log)
## RR_an_low <- 1 - exp(-(rr_log - 1.96 * rr_log_se))
## RR_an_high <- 1 - exp(-(rr_log + 1.96 * rr_log_se))
## tibble(rr_log = rr_log, rr_log_se = rr_log_se, RR = RR_an,
## RR_lci = RR_an_low, RR_uci = RR_an_high)
## }
print(ftab)
## function (m = time_lag0, mref = int_lag0, tests = c("LR", "wald"))
## {
## q <- qtest(m)
## het <- data.frame(Q = q$Q[1], df = q$df[1], Q_p = q$pvalue[1],
## I2 = (q$Q[1] - q$df[1])/q$Q[1] * 100)
## if (het$I2 < 0) {
## het$I2 = 0
## }
## if (m$method != "reml") {
## het <- cbind(het, data.frame(AIC = AIC(m), BIC = BIC(m)))
## }
## if ("LR" %in% tests & !is.null(mref) & m$method != "reml") {
## lrstat <- -2 * (logLik(mref) - logLik(m))[1]
## df <- attr(logLik(m), "df") - attr(logLik(mref), "df")
## pvalue <- 1 - pchisq(lrstat, df)
## lr <- data.frame(LR = lrstat, LR_df = df, LR_p = pvalue,
## stringsAsFactors = F)
## }
## if ("wald" %in% tests & !is.null(mref) & m$method != "reml") {
## coef <- coef(m)[-grep("Int", names(coef(m)))]
## vcov <- vcov(m)[-grep("Int", names(coef(m))), -grep("Int",
## names(coef(m)))]
## waldstat <- coef %*% solve(vcov) %*% coef
## df <- length(coef)
## pvalue <- 1 - pchisq(waldstat, df)
## wald <- data.frame(wald = waldstat, wald_df = df, wald_p = pvalue,
## stringsAsFactors = F)
## }
## if (class(m) == "mvmeta") {
## labs = data.frame(mod = paste(rownames(m$coefficients),
## collapse = " + "), ref = paste(rownames(mref$coefficients),
## collapse = " + "))
## }
## else if (class(m) == "mixmeta") {
## labs = data.frame(mod = paste(m$call$formula[3], collapse = ""),
## ref = paste(mref$call$formula[3], collapse = ""))
## }
## else {
## stop("model class not supported")
## }
## if (m$method != "reml") {
## if (!is.null(mref)) {
## return(cbind(labs, het, lr, wald))
## }
## else {
## return(cbind(labs, het, data.frame(LR = NA, LR_df = NA,
## LR_p = NA, wald = NA, wald_df = NA, wald_p = NA,
## stringsAsFactors = F)))
## }
## }
## else {
## if (!is.null(mref)) {
## return(cbind(labs, het))
## }
## else {
## return(cbind(labs, het))
## }
## }
## }
## <bytecode: 0x000000001db777d8>
devtools::session_info()
## - Session info ---------------------------------------------------------------
## setting value
## version R version 3.6.2 (2019-12-12)
## os Windows 10 x64
## system x86_64, mingw32
## ui RTerm
## language (EN)
## collate English_United States.1252
## ctype English_United States.1252
## tz America/Denver
## date 2020-03-04
##
## - Packages -------------------------------------------------------------------
## package * version date lib source
## assertthat 0.2.1 2019-03-21 [1] CRAN (R 3.6.0)
## backports 1.1.5 2019-10-02 [1] CRAN (R 3.6.1)
## broom 0.5.4 2020-01-27 [1] CRAN (R 3.6.2)
## callr 3.4.2 2020-02-12 [1] CRAN (R 3.6.2)
## cellranger 1.1.0 2016-07-27 [1] CRAN (R 3.6.0)
## cli 2.0.1 2020-01-08 [1] CRAN (R 3.6.2)
## colorspace 1.4-1 2019-03-18 [1] CRAN (R 3.6.0)
## crayon 1.3.4 2017-09-16 [1] CRAN (R 3.6.0)
## data.table * 1.12.8 2019-12-09 [1] CRAN (R 3.6.2)
## DBI 1.1.0 2019-12-15 [1] CRAN (R 3.6.2)
## dbplyr 1.4.2 2019-06-17 [1] CRAN (R 3.6.2)
## desc 1.2.0 2018-05-01 [1] CRAN (R 3.6.1)
## devtools 2.2.2 2020-02-17 [1] CRAN (R 3.6.2)
## digest 0.6.24 2020-02-12 [1] CRAN (R 3.6.2)
## dlnm * 2.3.9 2019-03-11 [1] CRAN (R 3.6.0)
## dplyr * 0.8.4 2020-01-31 [1] CRAN (R 3.6.2)
## ellipsis 0.3.0 2019-09-20 [1] CRAN (R 3.6.2)
## evaluate 0.14 2019-05-28 [1] CRAN (R 3.6.2)
## fansi 0.4.1 2020-01-08 [1] CRAN (R 3.6.2)
## forcats * 0.4.0 2019-02-17 [1] CRAN (R 3.6.0)
## fs 1.3.1 2019-05-06 [1] CRAN (R 3.6.0)
## generics 0.0.2 2018-11-29 [1] CRAN (R 3.6.0)
## ggplot2 * 3.2.1 2019-08-10 [1] CRAN (R 3.6.2)
## glue 1.3.1 2019-03-12 [1] CRAN (R 3.6.0)
## gridExtra * 2.3 2017-09-09 [1] CRAN (R 3.6.0)
## gtable 0.3.0 2019-03-25 [1] CRAN (R 3.6.0)
## haven 2.2.0 2019-11-08 [1] CRAN (R 3.6.2)
## here * 0.1 2017-05-28 [1] CRAN (R 3.6.0)
## highr 0.8 2019-03-20 [1] CRAN (R 3.6.0)
## hms 0.5.3 2020-01-08 [1] CRAN (R 3.6.2)
## htmltools 0.4.0 2019-10-04 [1] CRAN (R 3.6.2)
## httr 1.4.1 2019-08-05 [1] CRAN (R 3.6.2)
## jsonlite 1.6.1 2020-02-02 [1] CRAN (R 3.6.2)
## knitr 1.28 2020-02-06 [1] CRAN (R 3.6.2)
## lattice 0.20-38 2018-11-04 [2] CRAN (R 3.6.2)
## lazyeval 0.2.2 2019-03-15 [1] CRAN (R 3.6.0)
## lifecycle 0.1.0 2019-08-01 [1] CRAN (R 3.6.2)
## lubridate 1.7.4 2018-04-11 [1] CRAN (R 3.6.0)
## magrittr 1.5 2014-11-22 [1] CRAN (R 3.6.0)
## MASS * 7.3-51.5 2019-12-20 [1] CRAN (R 3.6.2)
## Matrix 1.2-18 2019-11-27 [2] CRAN (R 3.6.2)
## memoise 1.1.0 2017-04-21 [1] CRAN (R 3.6.1)
## mgcv * 1.8-31 2019-11-09 [1] CRAN (R 3.6.2)
## mixmeta * 1.0.7 2019-12-09 [1] CRAN (R 3.6.2)
## modelr 0.1.5 2019-08-08 [1] CRAN (R 3.6.2)
## munsell 0.5.0 2018-06-12 [1] CRAN (R 3.6.0)
## mvmeta * 1.0.3 2019-12-10 [1] CRAN (R 3.6.2)
## nlme * 3.1-144 2020-02-06 [1] CRAN (R 3.6.2)
## pillar 1.4.3 2019-12-20 [1] CRAN (R 3.6.2)
## pkgbuild 1.0.6 2019-10-09 [1] CRAN (R 3.6.2)
## pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 3.6.2)
## pkgload 1.0.2 2018-10-29 [1] CRAN (R 3.6.1)
## prettyunits 1.1.1 2020-01-24 [1] CRAN (R 3.6.2)
## processx 3.4.2 2020-02-09 [1] CRAN (R 3.6.2)
## ps 1.3.2 2020-02-13 [1] CRAN (R 3.6.2)
## purrr * 0.3.3 2019-10-18 [1] CRAN (R 3.6.2)
## R6 2.4.1 2019-11-12 [1] CRAN (R 3.6.2)
## Rcpp 1.0.3 2019-11-08 [1] CRAN (R 3.6.2)
## readr * 1.3.1 2018-12-21 [1] CRAN (R 3.6.0)
## readxl * 1.3.1 2019-03-13 [1] CRAN (R 3.6.0)
## remotes 2.1.1 2020-02-15 [1] CRAN (R 3.6.2)
## reprex 0.3.0 2019-05-16 [1] CRAN (R 3.6.0)
## rlang 0.4.4 2020-01-28 [1] CRAN (R 3.6.2)
## rmarkdown 2.1 2020-01-20 [1] CRAN (R 3.6.2)
## rprojroot 1.3-2 2018-01-03 [1] CRAN (R 3.6.0)
## rstudioapi 0.11 2020-02-07 [1] CRAN (R 3.6.2)
## rvest 0.3.5 2019-11-08 [1] CRAN (R 3.6.2)
## scales * 1.1.0 2019-11-18 [1] CRAN (R 3.6.2)
## sessioninfo 1.1.1 2018-11-05 [1] CRAN (R 3.6.1)
## stringi 1.4.3 2019-03-12 [1] CRAN (R 3.6.0)
## stringr * 1.4.0 2019-02-10 [1] CRAN (R 3.6.0)
## survival * 3.1-8 2019-12-03 [1] CRAN (R 3.6.2)
## testthat 2.3.1 2019-12-01 [1] CRAN (R 3.6.2)
## tibble * 2.1.3 2019-06-06 [1] CRAN (R 3.6.2)
## tictoc * 1.0 2014-06-17 [1] CRAN (R 3.6.0)
## tidyr * 1.0.2 2020-01-24 [1] CRAN (R 3.6.2)
## tidyselect 1.0.0 2020-01-27 [1] CRAN (R 3.6.2)
## tidyverse * 1.3.0 2019-11-21 [1] CRAN (R 3.6.2)
## tsModel * 0.6 2013-06-24 [1] CRAN (R 3.6.0)
## usethis 1.5.1 2019-07-04 [1] CRAN (R 3.6.1)
## vctrs 0.2.3 2020-02-20 [1] CRAN (R 3.6.2)
## withr 2.1.2 2018-03-15 [1] CRAN (R 3.6.0)
## xfun 0.12 2020-01-13 [1] CRAN (R 3.6.2)
## xml2 1.2.2 2019-08-09 [1] CRAN (R 3.6.2)
## yaml 2.2.1 2020-02-01 [1] CRAN (R 3.6.2)
##
## [1] C:/Users/LayC/R.Lib
## [2] C:/Program Files/R/R-3.6.2/library
## [3] C:/Users/LayC/Documents/R/win-library/3.6