100-seed-weight as dependent variable
Load Libraries
library(ezknitr)
library(knitr)
library(metafor)
library(devtools)
Clear environment and set seed
remove(list=ls())
set.seed(18837)
load(file="data/output_data/data_cleaned.R")
remove(list=c(#"seedwt.data.ROM",
"yield.data.ROM",
"rust.data.ROM"
))
Number of bootstrap simulations:
nsims <- 5000
Empty data frame to hold results
results.seedwt <- as.data.frame(matrix(NA,ncol=17,nrow=nsims))
colnames(results.seedwt) <- c("OVERALL","tau2",
"FLUT","MIXED","TEBU",
"Strobilurin","Triaz_Strob","Triazole",
"R3",
"Application Intrcpt","Application Slope",
"1 Application","2 Applications",
"high","medium","low",
"2006")
Running the bootstraps
for(ii in 1:nsims){
studyIDS <- data.frame(ReferenceNumb = sample(seedwt.data.ROM$ReferenceNumb,
length(unique(seedwt.data.ROM$ReferenceNumb)),
replace=T))
newdata <- merge(studyIDS, seedwt.data.ROM)
# Overall analysis
meta <- rma.uni(yi = yi,
vi = 1/rust.m2i,
data = newdata,
method = "REML")
results.seedwt$OVERALL[ii] <- meta$b
results.seedwt$tau2[ii] <- meta$tau2
# Active Ingredients
ai <- rma.uni(yi = yi,
vi = 1/rust.m2i,
data = newdata,
method = "REML",
mods = ~category_ai-1)
results.seedwt$FLUT[ii] <- ai$b[rownames(ai$b)=="category_aiFLUT"]
results.seedwt$MIXED[ii] <- ai$b[rownames(ai$b)=="category_aiMIXED"]
results.seedwt$TEBU[ii] <- ai$b[rownames(ai$b)=="category_aiTEBU"]
# Classes
class <- rma.uni(yi = yi,
vi = 1/rust.m2i,
data = newdata,
method = "REML",
mods = ~category_class-1)
results.seedwt$Strobilurin[ii] <-
class$b[rownames(class$b)=="category_classstrobilurin"]
results.seedwt$Triaz_Strob[ii] <-
class$b[rownames(class$b)=="category_classtriaz + strob"]
results.seedwt$Triazole[ii] <-
class$b[rownames(class$b)=="category_classtriazole"]
# Growth stage
if(nrow(newdata[!is.na(newdata$category_rstage),])>0){
rstage <- rma.uni(yi = yi,
vi = 1/rust.m2i,
data = newdata[! is.na(newdata$category_rstage),],
method = "REML")
results.seedwt$R3[ii] <- rstage$b
}
# Disease pressure
if(length(unique(newdata$category_pressure))>=2){
pressure <- rma.uni(yi = yi,
vi = (n1i+n2i)/(n1i*n2i),
data = newdata,
method = "REML",
mods = ~category_pressure-1)
if("low" %in% unique(newdata$category_pressure)){
results.seedwt$low[ii] <-
pressure$b[rownames(pressure$b)=="category_pressurelow"]
}
if("medium" %in% unique(newdata$category_pressure)){
results.seedwt$medium[ii] <-
pressure$b[rownames(pressure$b)=="category_pressuremedium"]
}
results.seedwt$high[ii] <- pressure$b[rownames(pressure$b)=="category_pressurehigh"]
}
# Number of applications
applications <- rma.uni(yi = yi,
vi = 1/rust.m2i,
data = newdata,
method = "REML",
mods = ~number_applications)
results.seedwt$`Application Intrcpt`[ii] <-
applications$b[rownames(applications$b)=="intrcpt"]
results.seedwt$`Application Slope`[ii] <-
applications$b[rownames(applications$b)=="number_applications"]
# Number of applications as category
apps.categ <- rma.uni(yi = yi,
vi = 1/rust.m2i,
data = newdata,
method = "REML",
mods = ~as.character(number_applications)-1)
results.seedwt$`1 Application`[ii] <-
apps.categ$b[rownames(apps.categ$b)=="as.character(number_applications)1"]
results.seedwt$`2 Applications`[ii] <-
apps.categ$b[rownames(apps.categ$b)=="as.character(number_applications)2"]
# Study year
# Condition on year1 = 2005
year <- rma.uni(yi = yi,
vi = (n1i+n2i)/(n1i*n2i),
data = newdata[newdata$category_year==2006,],
method = "REML")
results.seedwt$`2006`[ii] <- year$b
}
save(results.seedwt, file="data/output_data/results_seedwt.R")
Spun with ezspin(“programs/data_analysis_bootstrap_seedwt.R”, out_dir=“output”, fig_dir=“figures”, keep_md=FALSE)
Session Info:
devtools::session_info()
## Session info --------------------------------------------------------------
## setting value
## version R version 3.3.2 (2016-10-31)
## system x86_64, darwin13.4.0
## ui RStudio (0.99.902)
## language (EN)
## collate en_US.UTF-8
## tz America/Chicago
## date 2017-01-25
## Packages ------------------------------------------------------------------
## package * version date source
## devtools * 1.12.0 2016-06-24 CRAN (R 3.3.0)
## digest 0.6.10 2016-08-02 CRAN (R 3.3.0)
## doBy * 4.5-15 2016-03-31 CRAN (R 3.3.0)
## evaluate 0.10 2016-10-11 CRAN (R 3.3.0)
## ezknitr * 0.6 2016-09-16 CRAN (R 3.3.0)
## knitr * 1.15.1 2016-11-22 CRAN (R 3.3.2)
## lattice 0.20-34 2016-09-06 CRAN (R 3.3.2)
## magrittr 1.5 2014-11-22 CRAN (R 3.3.0)
## markdown 0.7.7 2015-04-22 CRAN (R 3.3.0)
## MASS 7.3-45 2016-04-21 CRAN (R 3.3.2)
## Matrix * 1.2-7.1 2016-09-01 CRAN (R 3.3.2)
## memoise 1.0.0 2016-01-29 CRAN (R 3.3.0)
## metafor * 1.9-9 2016-09-25 CRAN (R 3.3.0)
## R.methodsS3 1.7.1 2016-02-16 CRAN (R 3.3.0)
## R.oo 1.21.0 2016-11-01 CRAN (R 3.3.0)
## R.utils 2.5.0 2016-11-07 CRAN (R 3.3.0)
## rstudioapi 0.6 2016-06-27 CRAN (R 3.3.0)
## stringi 1.1.2 2016-10-01 CRAN (R 3.3.0)
## stringr 1.1.0 2016-08-19 CRAN (R 3.3.0)
## withr 1.0.2 2016-06-20 CRAN (R 3.3.0)