library(readr)
library(metafor)
library(tidyverse)
library(here)
rust_yld <- read_csv("data/dat-yld.csv")
##
## -- Column specification --------------------------------------------------------
## cols(
## study = col_double(),
## year = col_double(),
## location = col_character(),
## state = col_character(),
## n_spray = col_double(),
## brand_name = col_character(),
## mean_sev = col_double(),
## mean_yld = col_double(),
## v_sev = col_double(),
## v_yld = col_double(),
## yld_check = col_double(),
## v_yld_check = col_double(),
## sev_check = col_double(),
## v_sev_check = col_double(),
## vi_yld = col_double(),
## n2 = col_double(),
## n3 = col_double(),
## design = col_double()
## )
library(metafor)
mv_yld <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name,
random = list(~brand_name | study),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld
)
summary(mv_yld)
##
## Multivariate Meta-Analysis Model (k = 1314; method: ML)
##
## logLik Deviance AIC BIC AICc
## -9215.4624 4054.2612 18538.9248 18818.6897 18543.6429
##
## Variance Components:
##
## outer factor: study (nlvls = 175)
## inner factor: brand_name (nlvls = 9)
##
## estim sqrt k.lvl fixed level
## tau^2.1 636726.0927 797.9512 175 no AACHECK
## tau^2.2 645330.4099 803.3246 145 no AZOX + BENZ
## tau^2.3 779999.1995 883.1756 116 no BIXF + TFLX + PROT
## tau^2.4 676953.2688 822.7717 118 no PICO + BENZ
## tau^2.5 716365.1912 846.3836 175 no PICO + CYPR
## tau^2.6 717892.6609 847.2855 154 no PICO + TEBU
## tau^2.7 745271.0498 863.2908 116 no PYRA + EPOX + FLUX
## tau^2.8 726161.0416 852.1508 146 no TFLX + CYPR
## tau^2.9 703865.3637 838.9668 169 no TFLX + PROT
##
## rho.AACH rho.AZ+B rho.B+T+P rho.PI+B rho.PI+C rho.PI+T
## AACHECK 1 0.8246 0.8342 0.8355 0.9335 0.9119
## AZOX + BENZ 0.8246 1 0.9116 0.9743 0.9132 0.9103
## BIXF + TFLX + PROT 0.8342 0.9116 1 0.9205 0.9393 0.9322
## PICO + BENZ 0.8355 0.9743 0.9205 1 0.9281 0.9263
## PICO + CYPR 0.9335 0.9132 0.9393 0.9281 1 0.9860
## PICO + TEBU 0.9119 0.9103 0.9322 0.9263 0.9860 1
## PYRA + EPOX + FLUX 0.8470 0.9540 0.9744 0.9590 0.9357 0.9312
## TFLX + CYPR 0.9357 0.9216 0.9427 0.9268 0.9958 0.9834
## TFLX + PROT 0.8834 0.9209 0.9703 0.9306 0.9767 0.9670
## rho.P+E+F rho.TF+C rho.TF+P AACH AZ+B B+T+P PI+B
## AACHECK 0.8470 0.9357 0.8834 - no no no
## AZOX + BENZ 0.9540 0.9216 0.9209 145 - no no
## BIXF + TFLX + PROT 0.9744 0.9427 0.9703 116 87 - no
## PICO + BENZ 0.9590 0.9268 0.9306 118 88 116 -
## PICO + CYPR 0.9357 0.9958 0.9767 175 145 116 118
## PICO + TEBU 0.9312 0.9834 0.9670 154 124 95 97
## PYRA + EPOX + FLUX 1 0.9455 0.9475 116 87 114 116
## TFLX + CYPR 0.9455 1 0.9755 146 116 116 118
## TFLX + PROT 0.9475 0.9755 1 169 145 110 112
## PI+C PI+T P+E+F TF+C TF+P
## AACHECK no no no no no
## AZOX + BENZ no no no no no
## BIXF + TFLX + PROT no no no no no
## PICO + BENZ no no no no no
## PICO + CYPR - no no no no
## PICO + TEBU 154 - no no no
## PYRA + EPOX + FLUX 116 95 - no no
## TFLX + CYPR 146 125 116 - no
## TFLX + PROT 169 148 111 140 -
##
## Test for Residual Heterogeneity:
## QE(df = 1305) = 173239.6764, p-val < .0001
##
## Test of Moderators (coefficients 2:9):
## QM(df = 8) = 906.6343, p-val < .0001
##
## Model Results:
##
## estimate se zval pval ci.lb
## intrcpt 2460.9764 60.8511 40.4426 <.0001 2341.7103
## brand_nameAZOX + BENZ 910.6261 38.2758 23.7912 <.0001 835.6070
## brand_nameBIXF + TFLX + PROT 1080.0092 40.9167 26.3953 <.0001 999.8139
## brand_namePICO + BENZ 1010.9269 38.5883 26.1978 <.0001 935.2953
## brand_namePICO + CYPR 600.3912 25.2691 23.7599 <.0001 550.8646
## brand_namePICO + TEBU 682.1111 28.9720 23.5438 <.0001 625.3270
## brand_namePYRA + EPOX + FLUX 981.5161 38.9349 25.2092 <.0001 905.2050
## brand_nameTFLX + CYPR 646.4323 25.3750 25.4752 <.0001 596.6983
## brand_nameTFLX + PROT 891.0381 32.0831 27.7728 <.0001 828.1563
## ci.ub
## intrcpt 2580.2424 ***
## brand_nameAZOX + BENZ 985.6452 ***
## brand_nameBIXF + TFLX + PROT 1160.2044 ***
## brand_namePICO + BENZ 1086.5585 ***
## brand_namePICO + CYPR 649.9178 ***
## brand_namePICO + TEBU 738.8952 ***
## brand_namePYRA + EPOX + FLUX 1057.8271 ***
## brand_nameTFLX + CYPR 696.1663 ***
## brand_nameTFLX + PROT 953.9199 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
yield_res<- data.frame(cbind(mv_yld$b,
mv_yld$ci.lb,
mv_yld$ci.ub)) %>%
mutate(fungicide = c("check", "AZOX + BENZ", "BIXF + TFLX + PROT", "PICO + BENZ", "PICO + CYPR","PICO + TEBU", "PYRA + EPOX + FLUX", "TFLX + CYPR", "TFLX + PROT")) %>%
filter(fungicide != "check")
names (yield_res) = c("yld", "yld_lower", "yld_upper", "fungicide")
yield_res
#write_csv(yield_res, "data/yld_kg.csv")
openxlsx::write.xlsx(yield_res, here("data","yld_kg.xlsx"), colNames = TRUE)
We can set linear contrasts between treatments of interest and get the P-valued using the anova
function.
anova(mv_yld, L = rbind(
c(0, 1, -1, 0, 0, 0, 0, 0, 0),
c(0, 1, 0, -1, 0, 0, 0, 0, 0),
c(0, 1, 0, 0, -1, 0, 0, 0, 0),
c(0, 1, 0, 0, 0, -1, 0, 0, 0),
c(0, 1, 0, 0, 0, 0, -1, 0, 0),
c(0, 1, 0, 0, 0, 0, 0, -1, 0),
c(0, 1, 0, 0, 0, 0, 0, 0, -1),
c(0, 0, 1, -1, 0, 0, 0, 0, 0),
c(0, 0, 1, 0, -1, 0, 0, 0, 0),
c(0, 0, 1, 0, 0, -1, 0, 0, 0),
c(0, 0, 1, 0, 0, 0, -1, 0, 0),
c(0, 0, 1, 0, 0, 0, 0, -1, 0),
c(0, 0, 1, 0, 0, 0, 0, 0, -1),
c(0, 0, 0, 1, -1, 0, 0, 0, 0),
c(0, 0, 0, 1, 0, -1, 0, 0, 0),
c(0, 0, 0, 1, 0, 0, -1, 0, 0),
c(0, 0, 0, 1, 0, 0, 0, -1, 0),
c(0, 0, 0, 1, 0, 0, 0, 0, -1),
c(0, 0, 0, 0, 1, -1, 0, 0, 0),
c(0, 0, 0, 0, 1, 0, -1, 0, 0),
c(0, 0, 0, 0, 1, 0, 0, -1, 0),
c(0, 0, 0, 0, 1, 0, 0, 0, -1),
c(0, 0, 0, 0, 0, 1, -1, 0, 0),
c(0, 0, 0, 0, 0, 1, 0, -1, 0),
c(0, 0, 0, 0, 0, 1, 0, 0, -1),
c(0, 0, 0, 0, 0, 0, 1, -1, 0),
c(0, 0, 0, 0, 0, 0, 1, 0, -1),
c(0, 0, 0, 0, 0, 0, 0, 1, -1)
))
##
## Hypotheses:
## 1: brand_nameAZOX + BENZ - brand_nameBIXF + TFLX + PROT = 0
## 2: brand_nameAZOX + BENZ - brand_namePICO + BENZ = 0
## 3: brand_nameAZOX + BENZ - brand_namePICO + CYPR = 0
## 4: brand_nameAZOX + BENZ - brand_namePICO + TEBU = 0
## 5: brand_nameAZOX + BENZ - brand_namePYRA + EPOX + FLUX = 0
## 6: brand_nameAZOX + BENZ - brand_nameTFLX + CYPR = 0
## 7: brand_nameAZOX + BENZ - brand_nameTFLX + PROT = 0
## 8: brand_nameBIXF + TFLX + PROT - brand_namePICO + BENZ = 0
## 9: brand_nameBIXF + TFLX + PROT - brand_namePICO + CYPR = 0
## 10: brand_nameBIXF + TFLX + PROT - brand_namePICO + TEBU = 0
## 11: brand_nameBIXF + TFLX + PROT - brand_namePYRA + EPOX + FLUX = 0
## 12: brand_nameBIXF + TFLX + PROT - brand_nameTFLX + CYPR = 0
## 13: brand_nameBIXF + TFLX + PROT - brand_nameTFLX + PROT = 0
## 14: brand_namePICO + BENZ - brand_namePICO + CYPR = 0
## 15: brand_namePICO + BENZ - brand_namePICO + TEBU = 0
## 16: brand_namePICO + BENZ - brand_namePYRA + EPOX + FLUX = 0
## 17: brand_namePICO + BENZ - brand_nameTFLX + CYPR = 0
## 18: brand_namePICO + BENZ - brand_nameTFLX + PROT = 0
## 19: brand_namePICO + CYPR - brand_namePICO + TEBU = 0
## 20: brand_namePICO + CYPR - brand_namePYRA + EPOX + FLUX = 0
## 21: brand_namePICO + CYPR - brand_nameTFLX + CYPR = 0
## 22: brand_namePICO + CYPR - brand_nameTFLX + PROT = 0
## 23: brand_namePICO + TEBU - brand_namePYRA + EPOX + FLUX = 0
## 24: brand_namePICO + TEBU - brand_nameTFLX + CYPR = 0
## 25: brand_namePICO + TEBU - brand_nameTFLX + PROT = 0
## 26: brand_namePYRA + EPOX + FLUX - brand_nameTFLX + CYPR = 0
## 27: brand_namePYRA + EPOX + FLUX - brand_nameTFLX + PROT = 0
## 28: brand_nameTFLX + CYPR - brand_nameTFLX + PROT = 0
##
## Results:
## estimate se zval pval
## 1: -169.3830 33.3077 -5.0854 <.0001
## 2: -100.3007 22.8112 -4.3970 <.0001
## 3: 310.2349 29.1941 10.6266 <.0001
## 4: 228.5151 30.0648 7.6008 <.0001
## 5: -70.8899 27.2375 -2.6027 0.0093
## 6: 264.1939 28.4203 9.2959 <.0001
## 7: 19.5880 28.1047 0.6970 0.4858
## 8: 69.0823 32.6511 2.1158 0.0344
## 9: 479.6179 28.4152 16.8789 <.0001
## 10: 397.8981 30.0046 13.2612 <.0001
## 11: 98.4931 21.0235 4.6849 <.0001
## 12: 433.5769 28.1232 15.4171 <.0001
## 13: 188.9710 23.3893 8.0794 <.0001
## 14: 410.5357 28.4939 14.4079 <.0001
## 15: 328.8158 29.3148 11.2167 <.0001
## 16: 29.4108 25.5630 1.1505 0.2499
## 17: 364.4946 29.2291 12.4703 <.0001
## 18: 119.8888 28.2063 4.2504 <.0001
## 19: -81.7199 15.4180 -5.3003 <.0001
## 20: -381.1248 28.5673 -13.3413 <.0001
## 21: -46.0411 11.8314 -3.8914 <.0001
## 22: -290.6469 17.2112 -16.8871 <.0001
## 23: -299.4050 29.7364 -10.0686 <.0001
## 24: 35.6788 16.7421 2.1311 0.0331
## 25: -208.9271 20.1956 -10.3452 <.0001
## 26: 335.0838 27.1044 12.3627 <.0001
## 27: 90.4779 27.0003 3.3510 0.0008
## 28: -244.6059 18.0530 -13.5493 <.0001
mv_yld_year <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.numeric(year),
random = list(~brand_name | study),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld %>% mutate(year= year - 2015))
mv_yld_year
##
## Multivariate Meta-Analysis Model (k = 1314; method: ML)
##
## Variance Components:
##
## outer factor: study (nlvls = 175)
## inner factor: brand_name (nlvls = 9)
##
## estim sqrt k.lvl fixed level
## tau^2.1 620332.6378 787.6120 175 no AACHECK
## tau^2.2 644840.3574 803.0195 145 no AZOX + BENZ
## tau^2.3 766684.4458 875.6052 116 no BIXF + TFLX + PROT
## tau^2.4 675455.2431 821.8608 118 no PICO + BENZ
## tau^2.5 701904.5140 837.7974 175 no PICO + CYPR
## tau^2.6 700466.0081 836.9385 154 no PICO + TEBU
## tau^2.7 741471.5668 861.0874 116 no PYRA + EPOX + FLUX
## tau^2.8 710636.6290 842.9927 146 no TFLX + CYPR
## tau^2.9 693032.8659 832.4860 169 no TFLX + PROT
##
## rho.AACH rho.AZ+B rho.B+T+P rho.PI+B rho.PI+C rho.PI+T
## AACHECK 1 0.8504 0.8305 0.8509 0.9325 0.9139
## AZOX + BENZ 0.8504 1 0.9240 0.9770 0.9318 0.9383
## BIXF + TFLX + PROT 0.8305 0.9240 1 0.9318 0.9391 0.9304
## PICO + BENZ 0.8509 0.9770 0.9318 1 0.9395 0.9477
## PICO + CYPR 0.9325 0.9318 0.9391 0.9395 1 0.9878
## PICO + TEBU 0.9139 0.9383 0.9304 0.9477 0.9878 1
## PYRA + EPOX + FLUX 0.8491 0.9592 0.9759 0.9635 0.9381 0.9365
## TFLX + CYPR 0.9345 0.9413 0.9423 0.9389 0.9960 0.9843
## TFLX + PROT 0.8818 0.9362 0.9701 0.9415 0.9763 0.9691
## rho.P+E+F rho.TF+C rho.TF+P AACH AZ+B B+T+P PI+B
## AACHECK 0.8491 0.9345 0.8818 - no no no
## AZOX + BENZ 0.9592 0.9413 0.9362 145 - no no
## BIXF + TFLX + PROT 0.9759 0.9423 0.9701 116 87 - no
## PICO + BENZ 0.9635 0.9389 0.9415 118 88 116 -
## PICO + CYPR 0.9381 0.9960 0.9763 175 145 116 118
## PICO + TEBU 0.9365 0.9843 0.9691 154 124 95 97
## PYRA + EPOX + FLUX 1 0.9481 0.9490 116 87 114 116
## TFLX + CYPR 0.9481 1 0.9752 146 116 116 118
## TFLX + PROT 0.9490 0.9752 1 169 145 110 112
## PI+C PI+T P+E+F TF+C TF+P
## AACHECK no no no no no
## AZOX + BENZ no no no no no
## BIXF + TFLX + PROT no no no no no
## PICO + BENZ no no no no no
## PICO + CYPR - no no no no
## PICO + TEBU 154 - no no no
## PYRA + EPOX + FLUX 116 95 - no no
## TFLX + CYPR 146 125 116 - no
## TFLX + PROT 169 148 111 140 -
##
## Test for Residual Heterogeneity:
## QE(df = 1296) = 172423.1537, p-val < .0001
##
## Test of Moderators (coefficients 2:18):
## QM(df = 17) = 1016.1004, p-val < .0001
##
## Model Results:
##
## estimate se zval
## intrcpt 2270.5064 107.4955 21.1219
## brand_nameAZOX + BENZ 1161.6923 62.4000 18.6169
## brand_nameBIXF + TFLX + PROT 1025.1549 94.4270 10.8566
## brand_namePICO + BENZ 1336.6621 82.2547 16.2503
## brand_namePICO + CYPR 619.3209 45.1394 13.7202
## brand_namePICO + TEBU 627.3148 50.7278 12.3663
## brand_namePYRA + EPOX + FLUX 1055.0337 91.0230 11.5909
## brand_nameTFLX + CYPR 655.3253 47.8614 13.6922
## brand_nameTFLX + PROT 927.5425 57.2384 16.2049
## as.numeric(year) 80.0760 37.4522 2.1381
## brand_nameAZOX + BENZ:as.numeric(year) -109.8643 22.3741 -4.9103
## brand_nameBIXF + TFLX + PROT:as.numeric(year) 8.8058 29.9421 0.2941
## brand_namePICO + BENZ:as.numeric(year) -119.9367 26.4036 -4.5424
## brand_namePICO + CYPR:as.numeric(year) -7.9894 15.8947 -0.5026
## brand_namePICO + TEBU:as.numeric(year) 30.4794 18.5042 1.6472
## brand_namePYRA + EPOX + FLUX:as.numeric(year) -37.8911 28.6793 -1.3212
## brand_nameTFLX + CYPR:as.numeric(year) -4.8123 16.5700 -0.2904
## brand_nameTFLX + PROT:as.numeric(year) -15.6598 20.1060 -0.7789
## pval ci.lb ci.ub
## intrcpt <.0001 2059.8191 2481.1937
## brand_nameAZOX + BENZ <.0001 1039.3905 1283.9941
## brand_nameBIXF + TFLX + PROT <.0001 840.0814 1210.2284
## brand_namePICO + BENZ <.0001 1175.4459 1497.8783
## brand_namePICO + CYPR <.0001 530.8492 707.7925
## brand_namePICO + TEBU <.0001 527.8901 726.7395
## brand_namePYRA + EPOX + FLUX <.0001 876.6319 1233.4354
## brand_nameTFLX + CYPR <.0001 561.5187 749.1319
## brand_nameTFLX + PROT <.0001 815.3573 1039.7278
## as.numeric(year) 0.0325 6.6711 153.4809
## brand_nameAZOX + BENZ:as.numeric(year) <.0001 -153.7168 -66.0119
## brand_nameBIXF + TFLX + PROT:as.numeric(year) 0.7687 -49.8797 67.4912
## brand_namePICO + BENZ:as.numeric(year) <.0001 -171.6868 -68.1867
## brand_namePICO + CYPR:as.numeric(year) 0.6152 -39.1424 23.1636
## brand_namePICO + TEBU:as.numeric(year) 0.0995 -5.7881 66.7470
## brand_namePYRA + EPOX + FLUX:as.numeric(year) 0.1864 -94.1015 18.3193
## brand_nameTFLX + CYPR:as.numeric(year) 0.7715 -37.2889 27.6643
## brand_nameTFLX + PROT:as.numeric(year) 0.4361 -55.0668 23.7473
##
## intrcpt ***
## brand_nameAZOX + BENZ ***
## brand_nameBIXF + TFLX + PROT ***
## brand_namePICO + BENZ ***
## brand_namePICO + CYPR ***
## brand_namePICO + TEBU ***
## brand_namePYRA + EPOX + FLUX ***
## brand_nameTFLX + CYPR ***
## brand_nameTFLX + PROT ***
## as.numeric(year) *
## brand_nameAZOX + BENZ:as.numeric(year) ***
## brand_nameBIXF + TFLX + PROT:as.numeric(year)
## brand_namePICO + BENZ:as.numeric(year) ***
## brand_namePICO + CYPR:as.numeric(year)
## brand_namePICO + TEBU:as.numeric(year) .
## brand_namePYRA + EPOX + FLUX:as.numeric(year)
## brand_nameTFLX + CYPR:as.numeric(year)
## brand_nameTFLX + PROT:as.numeric(year)
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(mv_yld_year, btt = 11:18)
##
## Test of Moderators (coefficients 11:18):
## QM(df = 8) = 85.4079, p-val < .0001
yld_gain <- rust_yld %>%
mutate(gain = mean_yld - yld_check) %>%
filter(brand_name!= "AACHECK")
yld_gain <- yld_gain %>%
mutate(region = case_when(
state == "MT" ~ "North",
state == "BA" ~ "North",
state == "DF" ~ "North",
state == "MS" ~ "North",
state == "GO" ~ "North",
state == "TO" ~ "North",
state == "PR" ~ "South",
state == "MG" ~ "South",
state == "SP" ~ "South",
state == "RS" ~ "South"))
reg1 = data.frame(mv_yld_year$beta, mv_yld_year$ci.lb, mv_yld_year$ci.ub) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado1, into = c("lixo", "lado3"),sep = "brand_name") %>%
select(-lixo) %>%
filter(lado3 != "NA") %>%
mutate(mod = c(rep("intercept", 8), rep("slope", 8))) %>%
select(-lado2)
names(reg1) = c("fungicide", "mean", "ci.lb", "ci.ub", "mod")
mean = reg1 %>%
group_by(fungicide) %>%
select(1:2,5) %>%
spread(mod, mean)
names(mean) = c("fungicide", "intercept_mean", "slope_mean")
upper = reg1 %>%
group_by(fungicide) %>%
select(1,3,5) %>%
spread(mod, ci.lb)
names(upper) = c("fungicide", "intercept_upper", "slope_upper")
lower = reg1 %>%
group_by(fungicide) %>%
select(1,4:5) %>%
spread(mod, ci.ub)
names(lower) = c("fungicide", "intercept_lower", "slope_lower")
data_model = left_join(mean, lower, by= c("fungicide")) %>%
left_join(upper, by = c("fungicide"))
year = seq(0,6, by = 0.1)
fungicide = NULL
year_col = NULL
for(i in 1:length(data_model$fungicide)){
data_cache = yld_gain %>%
filter(brand_name == data_model$fungicide[i])
years = unique(data_cache$year)-2015
year = sort(years)
year = seq(first(year), last(year), by = 0.1)
year_col = c(year_col,year)
fungicide = c(fungicide, rep(data_model$fungicide[i], length(year)))
}
predicted = data.frame(year_col, fungicide) %>%
mutate(year = year_col) %>%
right_join(data_model, by = "fungicide") %>%
mutate(mean_gain = intercept_mean + slope_mean*year,
CIL = intercept_lower + slope_lower*year,
CIU = intercept_upper + slope_upper*year,
year = year+2015) %>%
mutate(brand_name = fungicide) %>%
filter(year <2020.2) %>%
dplyr::select(-fungicide)
predicted
openxlsx::write.xlsx(predicted, here("data","predicted_yield.xlsx"), colNames = TRUE)
library(janitor)
# create the log of the sev variable
rust_yld1 <- rust_yld %>%
filter(brand_name %in% c("AACHECK", "AZOX + BENZ")) %>%
group_by(study) %>%
mutate(n2 = n()) %>%
filter(n2 != 1)
rust_yld1 = rust_yld1 %>%
mutate(region = case_when(
state == "MT" ~ "North",
state == "BA" ~ "North",
state == "DF" ~ "North",
state == "MS" ~ "North",
state == "GO" ~ "North",
state == "TO" ~ "North",
state == "PR" ~ "South",
state == "MG" ~ "North",
state == "SP" ~ "South",
state == "RS" ~ "South"))
rust_yld1 %>%
tabyl(brand_name, region)
rust_yld1 %>%
tabyl(brand_name, state)
library(metafor)
yld_gain <- rust_yld %>%
mutate(gain = mean_yld - yld_check) %>%
filter(brand_name!= "AACHECK")
rust_yld_nor <- rust_yld1 %>%
filter(region == "North")
mv_yld_year <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.numeric(year),
random = list(~brand_name | study),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld_nor %>% mutate(year= year - 2015))
mv_yld_year
##
## Multivariate Meta-Analysis Model (k = 208; method: ML)
##
## Variance Components:
##
## outer factor: study (nlvls = 104)
## inner factor: brand_name (nlvls = 2)
##
## estim sqrt k.lvl fixed level
## tau^2.1 404762.2394 636.2093 104 no AACHECK
## tau^2.2 407851.7983 638.6328 104 no AZOX + BENZ
##
## rho.AACH rho.AZ+B AACH AZ+B
## AACHECK 1 0.7951 - no
## AZOX + BENZ 0.7951 1 104 -
##
## Test for Residual Heterogeneity:
## QE(df = 204) = 13868.9184, p-val < .0001
##
## Test of Moderators (coefficients 2:4):
## QM(df = 3) = 463.5621, p-val < .0001
##
## Model Results:
##
## estimate se zval pval
## intrcpt 2093.6094 106.3863 19.6793 <.0001
## brand_nameAZOX + BENZ 1162.1697 71.7404 16.1997 <.0001
## as.numeric(year) 145.3687 39.3894 3.6906 0.0002
## brand_nameAZOX + BENZ:as.numeric(year) -126.7211 26.6474 -4.7555 <.0001
## ci.lb ci.ub
## intrcpt 1885.0962 2302.1226 ***
## brand_nameAZOX + BENZ 1021.5611 1302.7783 ***
## as.numeric(year) 68.1669 222.5706 ***
## brand_nameAZOX + BENZ:as.numeric(year) -178.9490 -74.4933 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
reg1 = data.frame(mv_yld_year$beta, mv_yld_year$ci.lb, mv_yld_year$ci.ub, mv_yld_year$se) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado1, into = c("lixo", "lado3"),sep = "brand_name") %>%
select(-lixo) %>%
filter(lado3 != "NA") %>%
mutate(mod = c(rep("intercept", 1), rep("slope", 1))) %>%
mutate(region = rep("North", length(mod))) %>%
select(-lado2)
names(reg1) = c("fungicide", "mean", "ci.lb", "ci.ub", "SE", "mod", "region")
elatus_norte = reg1
library(metafor)
yld_gain <- rust_yld %>%
mutate(gain = mean_yld - yld_check) %>%
filter(brand_name!= "AACHECK")
rust_yld_south <- rust_yld1 %>%
filter(region == "South")
mv_yld_year <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.numeric(year),
random = list(~brand_name | study),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld_south %>% mutate(year= year - 2015))
mv_yld_year
##
## Multivariate Meta-Analysis Model (k = 82; method: ML)
##
## Variance Components:
##
## outer factor: study (nlvls = 41)
## inner factor: brand_name (nlvls = 2)
##
## estim sqrt k.lvl fixed level
## tau^2.1 521876.7119 722.4103 41 no AACHECK
## tau^2.2 636887.6651 798.0524 41 no AZOX + BENZ
##
## rho.AACH rho.AZ+B AACH AZ+B
## AACHECK 1 0.8587 - no
## AZOX + BENZ 0.8587 1 41 -
##
## Test for Residual Heterogeneity:
## QE(df = 78) = 9033.2665, p-val < .0001
##
## Test of Moderators (coefficients 2:4):
## QM(df = 3) = 211.6400, p-val < .0001
##
## Model Results:
##
## estimate se zval pval
## intrcpt 2443.3541 180.4830 13.5379 <.0001
## brand_nameAZOX + BENZ 1217.7822 108.3954 11.2346 <.0001
## as.numeric(year) 174.8334 81.1147 2.1554 0.0311
## brand_nameAZOX + BENZ:as.numeric(year) -149.9365 49.6502 -3.0199 0.0025
## ci.lb ci.ub
## intrcpt 2089.6139 2797.0942 ***
## brand_nameAZOX + BENZ 1005.3311 1430.2332 ***
## as.numeric(year) 15.8516 333.8152 *
## brand_nameAZOX + BENZ:as.numeric(year) -247.2492 -52.6238 **
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
reg1 = data.frame(mv_yld_year$beta, mv_yld_year$ci.lb, mv_yld_year$ci.ub, mv_yld_year$se) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado1, into = c("lixo", "lado3"),sep = "brand_name") %>%
select(-lixo) %>%
filter(lado3 != "NA") %>%
mutate(mod = c(rep("intercept", 1), rep("slope", 1))) %>%
mutate(region = rep("South", length(mod))) %>%
select(-lado2)
names(reg1) = c("fungicide", "mean", "ci.lb", "ci.ub", "SE", "mod", "region")
elatus_sul = reg1
library(janitor)
# create the log of the sev variable
rust_yld1 <- rust_yld %>%
filter(brand_name %in% c("AACHECK", "PICO + BENZ")) %>%
group_by(study) %>%
mutate(n2 = n()) %>%
filter(n2 != 1)
rust_yld1 = rust_yld1 %>%
mutate(region = case_when(
state == "MT" ~ "North",
state == "BA" ~ "North",
state == "DF" ~ "North",
state == "MS" ~ "North",
state == "GO" ~ "North",
state == "TO" ~ "North",
state == "PR" ~ "South",
state == "MG" ~ "North",
state == "SP" ~ "South",
state == "RS" ~ "South"))
rust_yld1 %>%
tabyl(brand_name, region)
rust_yld1 %>%
tabyl(brand_name, state)
library(metafor)
yld_gain <- rust_yld %>%
mutate(gain = mean_yld - yld_check) %>%
filter(brand_name!= "AACHECK")
rust_yld_nor <- rust_yld1 %>%
filter(region == "North")
mv_yld_year <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.numeric(year),
random = list(~brand_name | study),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld_nor %>% mutate(year= year - 2015))
mv_yld_year
##
## Multivariate Meta-Analysis Model (k = 176; method: ML)
##
## Variance Components:
##
## outer factor: study (nlvls = 88)
## inner factor: brand_name (nlvls = 2)
##
## estim sqrt k.lvl fixed level
## tau^2.1 694041.3466 833.0914 88 no AACHECK
## tau^2.2 798894.3137 893.8089 88 no PICO + BENZ
##
## rho.AACH rho.PI+B AACH PI+B
## AACHECK 1 0.8468 - no
## PICO + BENZ 0.8468 1 88 -
##
## Test for Residual Heterogeneity:
## QE(df = 172) = 22658.9911, p-val < .0001
##
## Test of Moderators (coefficients 2:4):
## QM(df = 3) = 336.0714, p-val < .0001
##
## Model Results:
##
## estimate se zval pval
## intrcpt 1779.3432 290.0553 6.1345 <.0001
## brand_namePICO + BENZ 1486.4364 173.8693 8.5492 <.0001
## as.numeric(year) 179.3583 82.3555 2.1779 0.0294
## brand_namePICO + BENZ:as.numeric(year) -154.9267 49.5106 -3.1292 0.0018
## ci.lb ci.ub
## intrcpt 1210.8453 2347.8411 ***
## brand_namePICO + BENZ 1145.6589 1827.2139 ***
## as.numeric(year) 17.9444 340.7721 *
## brand_namePICO + BENZ:as.numeric(year) -251.9657 -57.8877 **
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
reg1 = data.frame(mv_yld_year$beta, mv_yld_year$ci.lb, mv_yld_year$ci.ub, mv_yld_year$se) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado1, into = c("lixo", "lado3"),sep = "brand_name") %>%
select(-lixo) %>%
filter(lado3 != "NA") %>%
mutate(mod = c(rep("intercept", 1), rep("slope", 1))) %>%
mutate(region = rep("North", length(mod))) %>%
select(-lado2)
names(reg1) = c("fungicide", "mean", "ci.lb", "ci.ub", "SE", "mod", "region")
vess_norte = reg1
library(metafor)
yld_gain <- rust_yld %>%
mutate(gain = mean_yld - yld_check) %>%
filter(brand_name!= "AACHECK")
rust_yld_south <- rust_yld1 %>%
filter(region == "South")
mv_yld_year <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.numeric(year),
random = list(~brand_name | study),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld_south %>% mutate(year= year - 2015))
mv_yld_year
##
## Multivariate Meta-Analysis Model (k = 60; method: ML)
##
## Variance Components:
##
## outer factor: study (nlvls = 30)
## inner factor: brand_name (nlvls = 2)
##
## estim sqrt k.lvl fixed level
## tau^2.1 518346.7061 719.9630 30 no AACHECK
## tau^2.2 535819.5912 731.9970 30 no PICO + BENZ
##
## rho.AACH rho.PI+B AACH PI+B
## AACHECK 1 0.8879 - no
## PICO + BENZ 0.8879 1 30 -
##
## Test for Residual Heterogeneity:
## QE(df = 56) = 4877.1263, p-val < .0001
##
## Test of Moderators (coefficients 2:4):
## QM(df = 3) = 160.2478, p-val < .0001
##
## Model Results:
##
## estimate se zval pval
## intrcpt 3528.7061 441.8376 7.9864 <.0001
## brand_namePICO + BENZ 699.2209 230.0585 3.0393 0.0024
## as.numeric(year) -167.8161 134.7084 -1.2458 0.2128
## brand_namePICO + BENZ:as.numeric(year) 52.9786 70.8458 0.7478 0.4546
## ci.lb ci.ub
## intrcpt 2662.7202 4394.6919 ***
## brand_namePICO + BENZ 248.3145 1150.1274 **
## as.numeric(year) -431.8396 96.2074
## brand_namePICO + BENZ:as.numeric(year) -85.8766 191.8338
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
reg1 = data.frame(mv_yld_year$beta, mv_yld_year$ci.lb, mv_yld_year$ci.ub, mv_yld_year$se) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado1, into = c("lixo", "lado3"),sep = "brand_name") %>%
select(-lixo) %>%
filter(lado3 != "NA") %>%
mutate(mod = c(rep("intercept", 1), rep("slope", 1))) %>%
mutate(region = rep("South", length(mod))) %>%
select(-lado2)
names(reg1) = c("fungicide", "mean", "ci.lb", "ci.ub", "SE", "mod", "region")
vess_sul = reg1
data_two_fung = rbind(elatus_norte, elatus_sul, vess_norte, vess_sul)
openxlsx::write.xlsx(data_two_fung, here("data","Elatus_Vessarya_yld_simulation.xlsx"), colNames = TRUE)
rust_yld1 = rust_yld %>%
mutate(region = case_when(
state == "MT" ~ "North",
state == "BA" ~ "North",
state == "DF" ~ "North",
state == "MS" ~ "North",
state == "GO" ~ "North",
state == "TO" ~ "North",
state == "PR" ~ "South",
state == "MG" ~ "North",
state == "SP" ~ "South",
state == "RS" ~ "South"))
rust_yld1 %>%
tabyl(brand_name, region)
mv_yld_region <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.factor(region),
random = list(~factor(brand_name) | factor(study)),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld1)
mv_yld_region
##
## Multivariate Meta-Analysis Model (k = 1314; method: ML)
##
## Variance Components:
##
## outer factor: factor(study) (nlvls = 175)
## inner factor: factor(brand_name) (nlvls = 9)
##
## estim sqrt k.lvl fixed level
## tau^2.1 606007.6711 778.4649 175 no AACHECK
## tau^2.2 601937.4814 775.8463 145 no AZOX + BENZ
## tau^2.3 721679.6495 849.5173 116 no BIXF + TFLX + PROT
## tau^2.4 655626.7851 809.7078 118 no PICO + BENZ
## tau^2.5 653264.0342 808.2475 175 no PICO + CYPR
## tau^2.6 661692.6217 813.4449 154 no PICO + TEBU
## tau^2.7 720031.9783 848.5470 116 no PYRA + EPOX + FLUX
## tau^2.8 656849.4796 810.4625 146 no TFLX + CYPR
## tau^2.9 628275.7729 792.6385 169 no TFLX + PROT
##
## rho.AACH rho.AZ+B rho.B+T+P rho.PI+B rho.PI+C rho.PI+T
## AACHECK 1 0.8142 0.8273 0.8341 0.9320 0.9082
## AZOX + BENZ 0.8142 1 0.9057 0.9715 0.9069 0.9038
## BIXF + TFLX + PROT 0.8273 0.9057 1 0.9198 0.9353 0.9290
## PICO + BENZ 0.8341 0.9715 0.9198 1 0.9324 0.9317
## PICO + CYPR 0.9320 0.9069 0.9353 0.9324 1 0.9849
## PICO + TEBU 0.9082 0.9038 0.9290 0.9317 0.9849 1
## PYRA + EPOX + FLUX 0.8453 0.9502 0.9749 0.9580 0.9384 0.9355
## TFLX + CYPR 0.9350 0.9164 0.9391 0.9321 0.9952 0.9821
## TFLX + PROT 0.8801 0.9169 0.9685 0.9349 0.9747 0.9651
## rho.P+E+F rho.TF+C rho.TF+P AACH AZ+B B+T+P PI+B
## AACHECK 0.8453 0.9350 0.8801 - no no no
## AZOX + BENZ 0.9502 0.9164 0.9169 145 - no no
## BIXF + TFLX + PROT 0.9749 0.9391 0.9685 116 87 - no
## PICO + BENZ 0.9580 0.9321 0.9349 118 88 116 -
## PICO + CYPR 0.9384 0.9952 0.9747 175 145 116 118
## PICO + TEBU 0.9355 0.9821 0.9651 154 124 95 97
## PYRA + EPOX + FLUX 1 0.9499 0.9504 116 87 114 116
## TFLX + CYPR 0.9499 1 0.9728 146 116 116 118
## TFLX + PROT 0.9504 0.9728 1 169 145 110 112
## PI+C PI+T P+E+F TF+C TF+P
## AACHECK no no no no no
## AZOX + BENZ no no no no no
## BIXF + TFLX + PROT no no no no no
## PICO + BENZ no no no no no
## PICO + CYPR - no no no no
## PICO + TEBU 154 - no no no
## PYRA + EPOX + FLUX 116 95 - no no
## TFLX + CYPR 146 125 116 - no
## TFLX + PROT 169 148 111 140 -
##
## Test for Residual Heterogeneity:
## QE(df = 1296) = 155126.5694, p-val < .0001
##
## Test of Moderators (coefficients 2:18):
## QM(df = 17) = 977.5113, p-val < .0001
##
## Model Results:
##
## estimate se
## intrcpt 2351.0840 69.9976
## brand_nameAZOX + BENZ 894.3456 45.1105
## brand_nameBIXF + TFLX + PROT 1032.7464 47.2686
## brand_namePICO + BENZ 1009.8853 44.6194
## brand_namePICO + CYPR 554.8111 29.0685
## brand_namePICO + TEBU 645.6712 33.8301
## brand_namePYRA + EPOX + FLUX 967.0207 44.9301
## brand_nameTFLX + CYPR 593.7882 28.9020
## brand_nameTFLX + PROT 832.3037 36.8105
## as.factor(region)South 394.8497 132.2577
## brand_nameAZOX + BENZ:as.factor(region)South 58.9984 85.1539
## brand_nameBIXF + TFLX + PROT:as.factor(region)South 150.1983 90.9738
## brand_namePICO + BENZ:as.factor(region)South -36.6158 85.3171
## brand_namePICO + CYPR:as.factor(region)South 159.7342 54.6664
## brand_namePICO + TEBU:as.factor(region)South 129.9113 63.1823
## brand_namePYRA + EPOX + FLUX:as.factor(region)South 10.4437 86.4463
## brand_nameTFLX + CYPR:as.factor(region)South 182.2832 54.5530
## brand_nameTFLX + PROT:as.factor(region)South 202.9433 69.3671
## zval pval ci.lb
## intrcpt 33.5881 <.0001 2213.8912
## brand_nameAZOX + BENZ 19.8257 <.0001 805.9307
## brand_nameBIXF + TFLX + PROT 21.8485 <.0001 940.1016
## brand_namePICO + BENZ 22.6333 <.0001 922.4328
## brand_namePICO + CYPR 19.0863 <.0001 497.8378
## brand_namePICO + TEBU 19.0857 <.0001 579.3654
## brand_namePYRA + EPOX + FLUX 21.5228 <.0001 878.9594
## brand_nameTFLX + CYPR 20.5449 <.0001 537.1414
## brand_nameTFLX + PROT 22.6105 <.0001 760.1565
## as.factor(region)South 2.9855 0.0028 135.6293
## brand_nameAZOX + BENZ:as.factor(region)South 0.6928 0.4884 -107.9002
## brand_nameBIXF + TFLX + PROT:as.factor(region)South 1.6510 0.0987 -28.1072
## brand_namePICO + BENZ:as.factor(region)South -0.4292 0.6678 -203.8342
## brand_namePICO + CYPR:as.factor(region)South 2.9220 0.0035 52.5900
## brand_namePICO + TEBU:as.factor(region)South 2.0561 0.0398 6.0762
## brand_namePYRA + EPOX + FLUX:as.factor(region)South 0.1208 0.9038 -158.9880
## brand_nameTFLX + CYPR:as.factor(region)South 3.3414 0.0008 75.3613
## brand_nameTFLX + PROT:as.factor(region)South 2.9256 0.0034 66.9863
## ci.ub
## intrcpt 2488.2767 ***
## brand_nameAZOX + BENZ 982.7605 ***
## brand_nameBIXF + TFLX + PROT 1125.3911 ***
## brand_namePICO + BENZ 1097.3378 ***
## brand_namePICO + CYPR 611.7843 ***
## brand_namePICO + TEBU 711.9769 ***
## brand_namePYRA + EPOX + FLUX 1055.0819 ***
## brand_nameTFLX + CYPR 650.4350 ***
## brand_nameTFLX + PROT 904.4509 ***
## as.factor(region)South 654.0701 **
## brand_nameAZOX + BENZ:as.factor(region)South 225.8970
## brand_nameBIXF + TFLX + PROT:as.factor(region)South 328.5038 .
## brand_namePICO + BENZ:as.factor(region)South 130.6027
## brand_namePICO + CYPR:as.factor(region)South 266.8784 **
## brand_namePICO + TEBU:as.factor(region)South 253.7463 *
## brand_namePYRA + EPOX + FLUX:as.factor(region)South 179.8754
## brand_nameTFLX + CYPR:as.factor(region)South 289.2051 ***
## brand_nameTFLX + PROT:as.factor(region)South 338.9004 **
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(mv_yld_region, btt = 11:18)
##
## Test of Moderators (coefficients 11:18):
## QM(df = 8) = 34.0610, p-val < .0001
reg1 = data.frame(mv_yld_region$beta, mv_yld_region$ci.lb, mv_yld_region$ci.ub, mv_yld_region$se) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado2, into = c("lixo","lado3"),sep = "region") %>%
select(-lixo) %>%
separate(lado1, into = c("lixo","lado1"),sep = "brand_name" ) %>%
select(-lixo) %>%
filter(lado1 != "NA") %>%
mutate(n = seq(1:16))
names(reg1) = c("fungicide", "region", "mean_yld", "ci.lb_yld", "ci.ub_yld", "SE_yld","n")
reg2 = reg1 %>%
filter(n < 9) %>%
mutate(region = rep("North", length(fungicide)))
reg3 = reg1 %>%
filter(n > 8) %>%
mutate(region = rep("South", length(fungicide)))
reg4 = rbind(reg2,reg3) %>%
select(-"n")
openxlsx::write.xlsx(reg4, here("data","yld_SE_region.xlsx"), colNames = TRUE)
eff_log = read_excel("data/eff_log_region.xlsx")
yld_reg <- read_excel("data/yld_SE_region.xlsx")
data_simulation = full_join(eff_log, yld_reg, by = c("fungicide", "region"))
openxlsx::write.xlsx(data_simulation, here("data","data_simulation.xlsx"), colNames = TRUE)
reg1 = data.frame(mv_yld_region$beta, mv_yld_region$ci.lb, mv_yld_region$ci.ub) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado2, into = c("lixo","lado3"),sep = "region") %>%
select(-lixo) %>%
separate(lado1, into = c("lixo","lado1"),sep = "brand_name" ) %>%
select(-lixo) %>%
filter(lado1 != "NA") %>%
mutate(n = seq(1:16))
names(reg1) = c("fungicide", "class", "mean", "ci.lb", "ci.ub", "n")
reg2 = reg1 %>%
filter(n < 9) %>%
mutate(class = rep("North", length(fungicide)))
reg3 = reg1 %>%
filter(n > 8) %>%
mutate(class = rep("South", length(fungicide)))
reg4 = rbind(reg2,reg3)
mean = reg4%>%
group_by(fungicide) %>%
select(1:3) %>%
spread(class, mean) %>%
mutate(mean = North + South) %>%
select(1,4)
lower = reg4%>%
group_by(fungicide) %>%
select(1,2,4) %>%
spread(class, ci.lb) %>%
mutate(lower = North + South) %>%
select(1,4)
upper = reg4%>%
group_by(fungicide) %>%
select(1,2,5) %>%
spread(class, ci.ub) %>%
mutate(upper = North + South) %>%
select(1,4)
reg5 = left_join(mean, lower, by= c("fungicide")) %>%
left_join(upper, by = c("fungicide")) %>%
mutate(class = rep("South", length("fungicide"))) %>%
select("fungicide", "class", "mean", "lower", "upper")
North = reg4 %>%
filter(class == "North") %>%
select(1:5)
names(North) = c("fungicide", "class", "mean", "lower", "upper")
reg6 = full_join(North,reg5)
## Joining, by = c("fungicide", "class", "mean", "lower", "upper")
reg6
openxlsx::write.xlsx(reg6, here("data","yld_region.xlsx"), colNames = TRUE)
library(tidyverse)
rust_yld <- rust_yld %>%
mutate(
sev_check_class = case_when(
sev_check < 70 ~ "Low",
sev_check >= 70 ~ "High"))
rust_yld %>%
tabyl(brand_name, sev_check_class)
mv_yld_sev_class <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * as.factor(sev_check_class) ,
random = list(~factor(brand_name) | factor(study)),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld)
mv_yld_sev_class
##
## Multivariate Meta-Analysis Model (k = 1314; method: ML)
##
## Variance Components:
##
## outer factor: factor(study) (nlvls = 175)
## inner factor: factor(brand_name) (nlvls = 9)
##
## estim sqrt k.lvl fixed level
## tau^2.1 588997.6993 767.4619 175 no AACHECK
## tau^2.2 634330.7301 796.4488 145 no AZOX + BENZ
## tau^2.3 762949.9431 873.4701 116 no BIXF + TFLX + PROT
## tau^2.4 668811.0794 817.8087 118 no PICO + BENZ
## tau^2.5 691060.7120 831.3006 175 no PICO + CYPR
## tau^2.6 691070.0371 831.3062 154 no PICO + TEBU
## tau^2.7 735711.1697 857.7361 116 no PYRA + EPOX + FLUX
## tau^2.8 704707.8950 839.4688 146 no TFLX + CYPR
## tau^2.9 683393.4981 826.6762 169 no TFLX + PROT
##
## rho.AACH rho.AZ+B rho.B+T+P rho.PI+B rho.PI+C rho.PI+T
## AACHECK 1 0.8255 0.8343 0.8413 0.9336 0.9100
## AZOX + BENZ 0.8255 1 0.9120 0.9728 0.9121 0.9097
## BIXF + TFLX + PROT 0.8343 0.9120 1 0.9193 0.9381 0.9314
## PICO + BENZ 0.8413 0.9728 0.9193 1 0.9289 0.9278
## PICO + CYPR 0.9336 0.9121 0.9381 0.9289 1 0.9855
## PICO + TEBU 0.9100 0.9097 0.9314 0.9278 0.9855 1
## PYRA + EPOX + FLUX 0.8525 0.9568 0.9746 0.9585 0.9365 0.9324
## TFLX + CYPR 0.9371 0.9210 0.9413 0.9268 0.9958 0.9831
## TFLX + PROT 0.8829 0.9206 0.9696 0.9306 0.9760 0.9663
## rho.P+E+F rho.TF+C rho.TF+P AACH AZ+B B+T+P PI+B
## AACHECK 0.8525 0.9371 0.8829 - no no no
## AZOX + BENZ 0.9568 0.9210 0.9206 145 - no no
## BIXF + TFLX + PROT 0.9746 0.9413 0.9696 116 87 - no
## PICO + BENZ 0.9585 0.9268 0.9306 118 88 116 -
## PICO + CYPR 0.9365 0.9958 0.9760 175 145 116 118
## PICO + TEBU 0.9324 0.9831 0.9663 154 124 95 97
## PYRA + EPOX + FLUX 1 0.9455 0.9479 116 87 114 116
## TFLX + CYPR 0.9455 1 0.9745 146 116 116 118
## TFLX + PROT 0.9479 0.9745 1 169 145 110 112
## PI+C PI+T P+E+F TF+C TF+P
## AACHECK no no no no no
## AZOX + BENZ no no no no no
## BIXF + TFLX + PROT no no no no no
## PICO + BENZ no no no no no
## PICO + CYPR - no no no no
## PICO + TEBU 154 - no no no
## PYRA + EPOX + FLUX 116 95 - no no
## TFLX + CYPR 146 125 116 - no
## TFLX + PROT 169 148 111 140 -
##
## Test for Residual Heterogeneity:
## QE(df = 1296) = 167547.7250, p-val < .0001
##
## Test of Moderators (coefficients 2:18):
## QM(df = 17) = 983.1316, p-val < .0001
##
## Model Results:
##
## estimate se
## intrcpt 2272.7953 77.4679
## brand_nameAZOX + BENZ 1005.3144 49.6013
## brand_nameBIXF + TFLX + PROT 1154.1764 53.9432
## brand_namePICO + BENZ 1119.6025 50.1055
## brand_namePICO + CYPR 651.0735 32.9249
## brand_namePICO + TEBU 728.4738 37.7547
## brand_namePYRA + EPOX + FLUX 1088.0312 50.7783
## brand_nameTFLX + CYPR 707.2705 33.1093
## brand_nameTFLX + PROT 956.3937 41.7684
## as.factor(sev_check_class)Low 438.0983 118.3458
## brand_nameAZOX + BENZ:as.factor(sev_check_class)Low -220.7779 75.6310
## brand_nameBIXF + TFLX + PROT:as.factor(sev_check_class)Low -171.4304 81.4070
## brand_namePICO + BENZ:as.factor(sev_check_class)Low -250.2076 75.5580
## brand_namePICO + CYPR:as.factor(sev_check_class)Low -117.3651 50.2694
## brand_namePICO + TEBU:as.factor(sev_check_class)Low -107.1294 58.1481
## brand_namePYRA + EPOX + FLUX:as.factor(sev_check_class)Low -237.2981 76.3145
## brand_nameTFLX + CYPR:as.factor(sev_check_class)Low -137.6340 50.2045
## brand_nameTFLX + PROT:as.factor(sev_check_class)Low -151.6628 63.7592
## zval pval
## intrcpt 29.3385 <.0001
## brand_nameAZOX + BENZ 20.2679 <.0001
## brand_nameBIXF + TFLX + PROT 21.3962 <.0001
## brand_namePICO + BENZ 22.3449 <.0001
## brand_namePICO + CYPR 19.7745 <.0001
## brand_namePICO + TEBU 19.2949 <.0001
## brand_namePYRA + EPOX + FLUX 21.4271 <.0001
## brand_nameTFLX + CYPR 21.3617 <.0001
## brand_nameTFLX + PROT 22.8976 <.0001
## as.factor(sev_check_class)Low 3.7018 0.0002
## brand_nameAZOX + BENZ:as.factor(sev_check_class)Low -2.9191 0.0035
## brand_nameBIXF + TFLX + PROT:as.factor(sev_check_class)Low -2.1058 0.0352
## brand_namePICO + BENZ:as.factor(sev_check_class)Low -3.3115 0.0009
## brand_namePICO + CYPR:as.factor(sev_check_class)Low -2.3347 0.0196
## brand_namePICO + TEBU:as.factor(sev_check_class)Low -1.8424 0.0654
## brand_namePYRA + EPOX + FLUX:as.factor(sev_check_class)Low -3.1095 0.0019
## brand_nameTFLX + CYPR:as.factor(sev_check_class)Low -2.7415 0.0061
## brand_nameTFLX + PROT:as.factor(sev_check_class)Low -2.3787 0.0174
## ci.lb
## intrcpt 2120.9609
## brand_nameAZOX + BENZ 908.0976
## brand_nameBIXF + TFLX + PROT 1048.4498
## brand_namePICO + BENZ 1021.3974
## brand_namePICO + CYPR 586.5418
## brand_namePICO + TEBU 654.4759
## brand_namePYRA + EPOX + FLUX 988.5077
## brand_nameTFLX + CYPR 642.3774
## brand_nameTFLX + PROT 874.5292
## as.factor(sev_check_class)Low 206.1448
## brand_nameAZOX + BENZ:as.factor(sev_check_class)Low -369.0118
## brand_nameBIXF + TFLX + PROT:as.factor(sev_check_class)Low -330.9852
## brand_namePICO + BENZ:as.factor(sev_check_class)Low -398.2986
## brand_namePICO + CYPR:as.factor(sev_check_class)Low -215.8914
## brand_namePICO + TEBU:as.factor(sev_check_class)Low -221.0976
## brand_namePYRA + EPOX + FLUX:as.factor(sev_check_class)Low -386.8717
## brand_nameTFLX + CYPR:as.factor(sev_check_class)Low -236.0329
## brand_nameTFLX + PROT:as.factor(sev_check_class)Low -276.6287
## ci.ub
## intrcpt 2424.6296 ***
## brand_nameAZOX + BENZ 1102.5311 ***
## brand_nameBIXF + TFLX + PROT 1259.9031 ***
## brand_namePICO + BENZ 1217.8075 ***
## brand_namePICO + CYPR 715.6051 ***
## brand_namePICO + TEBU 802.4717 ***
## brand_namePYRA + EPOX + FLUX 1187.5548 ***
## brand_nameTFLX + CYPR 772.1636 ***
## brand_nameTFLX + PROT 1038.2582 ***
## as.factor(sev_check_class)Low 670.0518 ***
## brand_nameAZOX + BENZ:as.factor(sev_check_class)Low -72.5439 **
## brand_nameBIXF + TFLX + PROT:as.factor(sev_check_class)Low -11.8756 *
## brand_namePICO + BENZ:as.factor(sev_check_class)Low -102.1166 ***
## brand_namePICO + CYPR:as.factor(sev_check_class)Low -18.8388 *
## brand_namePICO + TEBU:as.factor(sev_check_class)Low 6.8388 .
## brand_namePYRA + EPOX + FLUX:as.factor(sev_check_class)Low -87.7245 **
## brand_nameTFLX + CYPR:as.factor(sev_check_class)Low -39.2350 **
## brand_nameTFLX + PROT:as.factor(sev_check_class)Low -26.6970 *
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(mv_yld_sev_class, btt = 11:18)
##
## Test of Moderators (coefficients 11:18):
## QM(df = 8) = 15.1418, p-val = 0.0564
reg1 = data.frame(mv_yld_sev_class$beta, mv_yld_sev_class$ci.lb, mv_yld_sev_class$ci.ub) %>%
rownames_to_column("trat") %>%
separate(trat, into = c("lado1", "lado2"), sep = ":") %>%
separate(lado2, into = c("lixo","lado3"),sep = "sev_check_class") %>%
select(-lixo) %>%
separate(lado1, into = c("lixo","lado1"),sep = "brand_name" ) %>%
select(-lixo) %>%
filter(lado1 != "NA") %>%
mutate(n = seq(1:16))
names(reg1) = c("fungicide", "class", "mean", "ci.lb", "ci.ub", "n")
reg2 = reg1 %>%
filter(n < 9) %>%
mutate(class = rep("High", length(fungicide)))
reg3 = reg1 %>%
filter(n > 8) %>%
mutate(class = rep("Low", length(fungicide)))
reg4 = rbind(reg2,reg3)
mean = reg4%>%
group_by(fungicide) %>%
select(1:3) %>%
spread(class, mean) %>%
mutate(mean = High + Low) %>%
select(1,4)
lower = reg4%>%
group_by(fungicide) %>%
select(1,2,4) %>%
spread(class, ci.lb) %>%
mutate(lower = High + Low) %>%
select(1,4)
upper = reg4%>%
group_by(fungicide) %>%
select(1,2,5) %>%
spread(class, ci.ub) %>%
mutate(upper = High + Low) %>%
select(1,4)
reg5 = left_join(mean, lower, by= c("fungicide")) %>%
left_join(upper, by = c("fungicide")) %>%
mutate(class = rep("Low", length("fungicide"))) %>%
select("fungicide", "class", "mean", "lower", "upper")
High = reg4 %>%
filter(class == "High") %>%
select(1:5)
names(High) = c("fungicide", "class", "mean", "lower", "upper")
reg6 = full_join(High,reg5)
## Joining, by = c("fungicide", "class", "mean", "lower", "upper")
reg6
openxlsx::write.xlsx(reg6, here("data","yld_baseline.xlsx"), colNames = TRUE)
We used a factorial-type ANOVA model to determine the significance of the treatment x design interaction, evaluated based on the Wald test statistic.
mv_yld_design <- rma.mv(mean_yld, vi_yld,
mods = ~brand_name * design,
random = list(~ 1 | study / design / brand_name),
struct = "UN",
method = "ML",
control = list(optimizer = "nlm"),
data = rust_yld)
mv_yld_design
##
## Multivariate Meta-Analysis Model (k = 1314; method: ML)
##
## Variance Components:
##
## estim sqrt nlvls fixed factor
## sigma^2.1 326626.5655 571.5125 175 no study
## sigma^2.2 326626.5651 571.5125 175 no study/design
## sigma^2.3 44199.7132 210.2373 1314 no study/design/brand_name
##
## Test for Residual Heterogeneity:
## QE(df = 1296) = 169160.9325, p-val < .0001
##
## Test of Moderators (coefficients 2:18):
## QM(df = 17) = 2377.4326, p-val < .0001
##
## Model Results:
##
## estimate se zval pval
## intrcpt 2593.4693 106.3961 24.3756 <.0001
## brand_nameAZOX + BENZ 757.4275 42.5418 17.8043 <.0001
## brand_nameBIXF + TFLX + PROT 987.2565 49.6496 19.8845 <.0001
## brand_namePICO + BENZ 1028.2639 47.7568 21.5313 <.0001
## brand_namePICO + CYPR 587.2696 41.6389 14.1039 <.0001
## brand_namePICO + TEBU 678.9830 42.9294 15.8163 <.0001
## brand_namePYRA + EPOX + FLUX 911.4023 48.1021 18.9473 <.0001
## brand_nameTFLX + CYPR 631.5014 42.8403 14.7408 <.0001
## brand_nameTFLX + PROT 855.6228 41.6484 20.5439 <.0001
## design -29.2953 19.6717 -1.4892 0.1364
## brand_nameAZOX + BENZ:design 40.2080 7.7787 5.1690 <.0001
## brand_nameBIXF + TFLX + PROT:design 26.6692 16.7424 1.5929 0.1112
## brand_namePICO + BENZ:design -42.0625 15.0385 -2.7970 0.0052
## brand_namePICO + CYPR:design 0.1022 7.7117 0.0133 0.9894
## brand_namePICO + TEBU:design -0.7833 7.7578 -0.1010 0.9196
## brand_namePYRA + EPOX + FLUX:design 3.3399 15.3704 0.2173 0.8280
## brand_nameTFLX + CYPR:design 2.8844 9.1329 0.3158 0.7521
## brand_nameTFLX + PROT:design 7.3169 7.7320 0.9463 0.3440
## ci.lb ci.ub
## intrcpt 2384.9368 2802.0019 ***
## brand_nameAZOX + BENZ 674.0471 840.8079 ***
## brand_nameBIXF + TFLX + PROT 889.9450 1084.5679 ***
## brand_namePICO + BENZ 934.6624 1121.8655 ***
## brand_namePICO + CYPR 505.6588 668.8804 ***
## brand_namePICO + TEBU 594.8428 763.1231 ***
## brand_namePYRA + EPOX + FLUX 817.1240 1005.6806 ***
## brand_nameTFLX + CYPR 547.5360 715.4668 ***
## brand_nameTFLX + PROT 773.9934 937.2523 ***
## design -67.8511 9.2606
## brand_nameAZOX + BENZ:design 24.9619 55.4540 ***
## brand_nameBIXF + TFLX + PROT:design -6.1454 59.4838
## brand_namePICO + BENZ:design -71.5374 -12.5875 **
## brand_namePICO + CYPR:design -15.0124 15.2168
## brand_namePICO + TEBU:design -15.9884 14.4217
## brand_namePYRA + EPOX + FLUX:design -26.7854 33.4652
## brand_nameTFLX + CYPR:design -15.0157 20.7845
## brand_nameTFLX + PROT:design -7.8376 22.4714
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(mv_yld_design, btt = 11:18)
##
## Test of Moderators (coefficients 11:18):
## QM(df = 8) = 55.7813, p-val < .0001