Meta-analysis: Soybean yield ~ Target spot
Oct 2016
1 Introduction
Soybean target spot (TS) caused by Corynespora cassiicola has been considered a disease of limited importance in Brazil since it first report in 1976. However, due the massive adoption of no-till cultivation practices and the sowing of susceptible cultivars, in recent years this disease has spread throughout the Brazilian soybean-growing area with a considerably incidence.
By the moment, there is not a study quantifying the relationship between soybean yield (W, kg/ha) and target spot severity (TSs) in Brazil combining several growing seasons.
Our main objectives in this work were: using meta-analytic techniques
estimate effect size of metrics that characterize the relationship W~TSs, and
to identify variables that can moderate the effects in case of significant heterogeneity among studies.
A total of 33 Uniform Fungicide Trials of TS carried out across five brazilian states (Goiás, Mato Grosso, Mato Grosso do Sul, Paraná, and Tocantins) during four growing seasons (2012-2015, years of harvesting) was available to study the relationship W~TSs.
Each of the 33 selected trials (i-studies) constituted an independent study for this analysis.
Three effect-sizes related to the relationship W~TSs were estimated for each study: Pearson’s r correlation coefficient (r); and both linear regression coefficients, the intercept (β0) and slope (β1). The robust dispersion estimator, interquartile range (IQR), also called the middle fifty, being equal to the difference between the upper and lower quartiles, IQR=Q3−Q1, was estimated for each parameter.
2 Effect size based on correlations
2.1 Dataset
r: Pearson’s r coefficients (estimated with n pairs of plots)
yi: Fisher’s z coefficients
vi: sampling variance (within-study variance for study (i))
wts: weights (1/vi, %)
## study year state cultivar Dis_level Class r pairs yi vi wts
## 1 1 2012 MT TMG803 3High Low -0.181 36 -0.1830 0.0303 2.867072
## 2 2 2012 MT TMG803 2Medium High -0.525 36 -0.5832 0.0303 2.867072
## 3 3 2012 MT TMG803 3High High -0.349 36 -0.3643 0.0303 2.867072
## 4 4 2012 MS 5G830_RR 2Medium Low -0.186 36 -0.1882 0.0303 2.867072
## 5 5 2012 MS BMX_Potencia_RR 2Medium Low -0.187 45 -0.1892 0.0238 3.649001
## 6 6 2012 MT TMG803 3High Low -0.493 36 -0.5400 0.0303 2.867072- Description of Pearson’s r
## [[1]]
## [1] "range: -0.9 to 0.15"
##
## [[2]]
## [1] "median: -0.5"
##
## [[3]]
## [1] "IQR: 0.4"- Description of Fisher’s z
## [[1]]
## [1] "range: -1.3 to 0.15"
##
## [[2]]
## [1] "median: -0.54"
##
## [[3]]
## [1] "IQR: 0.53"The sampling distribution of the Pearson’s r correlation coefficients is typically very skewed. However, standard meta-analytic models assume that the sampling distribution of the outcome measure is (at least approximately) normal and that assumption is much more appropriate for the transformed Fisher’s z coefficient.
2.2 Random effect model
\[ \gamma_i = \mu + b_i + e \] where:
\(\gamma_i\) is the vector of estimated z;
\(\mu\) is the population average of z;
\(b_i\) is the estimated random effect of the i-study (between-study variance)
Random-effects model allows \(\mu_i\) to differ between studies under the assumption of a normal distribution around \(\mu\)
The i-studies included in the meta-analysis are assumed to be a random sample of the relevant distribution of effects, and the combined effect estimates the “mean effect” in this distribution (a weighted mean).
Hypothesis:
\(H_0 : \mu = 0\) ()
\(H_1 : \mu \neq 0\) ()
To test this hypothesis we compute Q, which is basically a weighted sum of squares (we compute the difference of every effect size from the mean effect size, square that difference, assign larger weights to more precise studies, and then sum these weighted values).
If the null hypothesis is true (that all the variation in effects is due to sampling error), the expected value of Q is equal to the number of studies minus 1 (here, 33-1 = 32).
##
## Random-Effects Model (k = 33; tau^2 estimator: ML)
##
## tau^2 (estimated amount of total heterogeneity): 0.087 (SE = 0.028)
## tau (square root of estimated tau^2 value): 0.294
## I^2 (total heterogeneity / total variability): 75.13%
## H^2 (total variability / sampling variability): 4.02
##
## Test for Heterogeneity:
## Q(df = 32) = 134.559, p-val < .001
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## -0.532 0.059 -9.001 <.001 -0.648 -0.417 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##
## estimate ci.lb ci.ub
## tau^2 0.09 0.05 0.18
## tau 0.29 0.22 0.42
## I^2(%) 75.13 62.55 86.08
## H^2 4.02 2.67 7.18
The mean effect size is -0.532 (-0.648 to -0.417) and the observed Q value is 134.559. This is more than we would expect if the null is true (32). Therefore, we reject the null. We have evidence that the true effect size varies from study to study.
\(I^2\) = 75.13%. This tells us that about 75% of the variance that we see in the forest plot reflects difference in the true effect sizes, while the other 25% reflects sampling error. Importantly, \(I^2\) is a proportion – it tells us what proportion of the observed variance is real (if our esimates are correct) but does not tell us how much variance there is.
- Random forest
- Prediction intervals for z
## pred se ci.lb ci.ub cr.lb cr.ub
## -0.532 0.059 -0.648 -0.417 -1.121 0.056- Prediction intervals for r
## pred ci.lb ci.ub cr.lb cr.ub
## -0.487 -0.571 -0.394 -0.808 0.056- Describe the eblups for the z vlaues
## [[1]]
## [1] "range: -1.1 to -0.01"
##
## [[2]]
## [1] "median: -0.54"
##
## [[3]]
## [1] "IQR: 0.39"2.3 Mixed model - inclusion of moderator variables
\[ \gamma_i = \mu + b_i + \delta_k + e \] where:
\(\gamma_i\) is the vector of estimated z;
\(\mu\) is the population average of z;
\(b_i\) is the estimated random effect (perturbation) of the i-study.
\(\delta_k\) is the fixed effect of k-level of the moderator variable.
\(e\) is the residual
## # A tibble: 3 × 3
## Moderator `Test of moderator` R2
## <chr> <dbl> <dbl>
## 1 Season 5.195488e-01 1.568190
## 2 Disease pressure 4.121975e-01 7.086112
## 3 Yield response 3.985073e-08 63.382219- Yield response -> significative moderator variable
##
## Mixed-Effects Model (k = 33; tau^2 estimator: ML)
##
## logLik deviance AIC BIC AICc
## -0.431 57.281 6.862 11.352 7.690
##
## tau^2 (estimated amount of residual heterogeneity): 0.032 (SE = 0.015)
## tau (square root of estimated tau^2 value): 0.178
## I^2 (residual heterogeneity / unaccounted variability): 52.51%
## H^2 (unaccounted variability / sampling variability): 2.11
## R^2 (amount of heterogeneity accounted for): 63.38%
##
## Test for Residual Heterogeneity:
## QE(df = 31) = 69.546, p-val < .001
##
## Test of Moderators (coefficient(s) 2):
## QM(df = 1) = 30.157, p-val < .001
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt -0.761 0.060 -12.719 <.001 -0.878 -0.644 ***
## ClassLow 0.470 0.086 5.492 <.001 0.303 0.638 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## df AIC BIC AICc logLik LRT pval QE tau^2 R^2
## Full 3 6.8622 11.3517 7.6897 -0.4311 69.5455 0.0317
## Reduced 2 26.1847 29.1777 26.5847 -11.0924 21.3225 <.0001 134.5586 0.0866 63.38%
## pred se pi.lb pi.ub
## 1 -0.236 0.128 -0.487 0.015
## 2 -0.670 0.128 -0.921 -0.420
## 3 -0.558 0.128 -0.809 -0.308
## 4 -0.238 0.128 -0.489 0.013
## 5 -0.233 0.120 -0.467 0.002
## 6 -0.418 0.128 -0.669 -0.167
## 7 -0.301 0.128 -0.552 -0.050
## 8 -0.202 0.124 -0.445 0.041
## 9 -0.919 0.124 -1.162 -0.676
## 10 -0.911 0.124 -1.153 -0.668
## 11 -0.753 0.124 -0.996 -0.510
## 12 -0.575 0.124 -0.818 -0.332
## 13 -0.865 0.124 -1.107 -0.622
## 14 -0.541 0.124 -0.784 -0.298
## 15 -0.446 0.132 -0.705 -0.186
## 16 -0.312 0.128 -0.563 -0.061
## 17 -0.714 0.130 -0.969 -0.459
## 18 -0.273 0.128 -0.524 -0.022
## 19 -0.919 0.132 -1.178 -0.660
## 20 -0.772 0.131 -1.029 -0.515
## 21 -0.660 0.135 -0.924 -0.396
## 22 -0.755 0.128 -1.006 -0.505
## 23 -0.718 0.124 -0.961 -0.475
## 24 -0.573 0.124 -0.817 -0.330
## 25 -0.196 0.124 -0.439 0.048
## 26 -0.234 0.124 -0.477 0.009
## 27 -0.784 0.124 -1.027 -0.541
## 28 -0.456 0.124 -0.699 -0.213
## 29 -0.796 0.124 -1.039 -0.553
## 30 -0.054 0.124 -0.297 0.189
## 31 -0.186 0.124 -0.429 0.057
## 32 -0.293 0.124 -0.536 -0.050
## 33 -1.030 0.124 -1.273 -0.788## [[1]]
## [1] "range: -0.6 to -0.28"
##
## [[2]]
## [1] "median: -0.64"
##
## [[3]]
## [1] "IQR: 0.36"- Fisher’s z - Output
## pred se ci.lb ci.ub cr.lb cr.ub
## 1 -0.2907 0.0613 -0.4108 -0.1705 -0.6599 0.0785
## 2 -0.7612 0.0598 -0.8785 -0.6439 -1.1294 -0.3929
## 3 -0.7612 0.0598 -0.8785 -0.6439 -1.1294 -0.3929
## 4 -0.2907 0.0613 -0.4108 -0.1705 -0.6599 0.0785
## 5 -0.2907 0.0613 -0.4108 -0.1705 -0.6599 0.0785
## 6 -0.2907 0.0613 -0.4108 -0.1705 -0.6599 0.0785- Pearson’s r - Output
## est low up
## ClassHigh -0.6417667 -0.7056489 -0.5675355
## ClassLow -0.2827587 -0.3891921 -0.16887503 Effect size of linear regression parameters
3.1 Dataset
## study Class b0_est b1_est Varb0 Covb0b1 Varb1
## 1 1 Low 2898.754 -7.584689 63541.922 -1698.2064 50.193294
## 2 2 High 3927.619 -24.469702 40514.167 -1306.7423 46.393218
## 3 3 High 3251.294 -16.095129 51817.014 -1633.7182 55.069153
## 4 4 Low 3445.393 -4.360786 7543.201 -264.3786 15.526315
## 5 5 Low 4218.295 -9.249662 16681.045 -860.1509 55.137876
## 6 6 Low 3610.224 -7.227121 4751.120 -129.5709 4.781703Description of the i-study intercepts and slopes
## [[1]]
## [1] "range: 1612.1 to 4849.9"
##
## [[2]]
## [1] "median: 3601.63"
##
## [[3]]
## [1] "IQR: 632.33"## [[1]]
## [1] "range: -43.5 to 7"
##
## [[2]]
## [1] "median: -16.07"
##
## [[3]]
## [1] "IQR: 21.4"Plot of individual linear regressions
3.2 Random model
Definiton of the model outcomes, moderator and var-cov
Statistics of the model
## Call: mvmeta(formula = outcomes ~ 1, S = varcov, method = "ml", bscov = "un")
##
## Multivariate random-effects meta-analysis
## Dimension: 2
## Estimation method: ML
##
## Fixed-effects coefficients
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## b0_est 3507.1 118.4 29.6 0.0 3275.1 3739.1 ***
## b1_est -17.1 2.1 -8.0 0.0 -21.2 -12.9 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Between-study random-effects (co)variance components
## Structure: General positive-definite
## Std. Dev Corr
## b0_est 669.5 b0_est
## b1_est 11.0 -0.1
##
## Multivariate Cochran Q-test for heterogeneity:
## Q = 11951.4 (df = 64), p-value = 0.0
## I-square statistic = 99.5%
##
## 33 studies, 66 observations, 2 fixed and 3 random-effects parameters
## logLik AIC BIC
## -392.1 794.3 805.2## Multivariate Cochran Q-test for heterogeneity
##
## Overall test:
## Q = 11951.444 (df = 64), p-value = 0.000
##
## Tests on single outcome parameters:
## Q df p-value
## b0_est 2473.625 32 0.000 ***
## b1_est 537.313 32 0.000 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## b0_est b1_est
## b0_est 448202.7186 -753.4152
## b1_est -753.4152 120.8191Description of the EBLUPS
## [[1]]
## [1] "range: 1631.5 to 4617.21"
##
## [[2]]
## [1] "median: 3601.49"
##
## [[3]]
## [1] "IQR: 623.69"## [[1]]
## [1] "range: -34 to -0.85"
##
## [[2]]
## [1] "median: -16.11"
##
## [[3]]
## [1] "IQR: 18.54"Output table: fitted’s + random’s = EBLUPs
## study b0_fitted b1_fitted b0_random b1_random b0_blup b1_blup
## 1 1 3507.093 -17.06217 -475.36741 5.72561309 3031.726 -11.3365557
## 12 2 3507.093 -17.06217 345.28735 -4.82988503 3852.380 -21.8920538
## 23 3 3507.093 -17.06217 -231.43229 0.20430500 3275.661 -16.8578638
## 28 4 3507.093 -17.06217 -37.48272 11.26598604 3469.610 -5.7961827
## 29 5 3507.093 -17.06217 734.09824 5.98100929 4241.191 -11.0811595
## 30 6 3507.093 -17.06217 111.79587 9.49762306 3618.889 -7.5645457
## 31 7 3507.093 -17.06217 357.98755 4.33421038 3865.081 -12.7279584
## 32 8 3507.093 -17.06217 409.43860 11.78264621 3916.532 -5.2795225
## 33 9 3507.093 -17.06217 212.81831 -5.34437451 3719.911 -22.4065433
## 2 10 3507.093 -17.06217 718.43690 -16.39223358 4225.530 -33.4544023
## 3 11 3507.093 -17.06217 95.09096 0.94826522 3602.184 -16.1139035
## 4 12 3507.093 -17.06217 -1875.54367 -0.01953077 1631.549 -17.0816995
## 5 13 3507.093 -17.06217 -174.94323 -16.25559327 3332.150 -33.3177620
## 6 14 3507.093 -17.06217 -1491.90434 1.79091819 2015.189 -15.2712506
## 7 15 3507.093 -17.06217 380.07254 5.92334730 3887.166 -11.1388215
## 8 16 3507.093 -17.06217 53.59440 10.71016860 3560.687 -6.3520002
## 9 17 3507.093 -17.06217 216.29994 -8.13095750 3723.393 -25.1931263
## 10 18 3507.093 -17.06217 337.30200 9.58276645 3844.395 -7.4794023
## 11 19 3507.093 -17.06217 879.70583 -12.13746833 4386.799 -29.1996371
## 13 20 3507.093 -17.06217 -265.70125 -16.89862558 3241.392 -33.9607943
## 14 21 3507.093 -17.06217 -1282.94604 -9.21226451 2224.147 -26.2744333
## 15 22 3507.093 -17.06217 -530.31966 -9.04505987 2976.773 -26.1072286
## 16 23 3507.093 -17.06217 337.16523 -11.76603735 3844.258 -28.8282061
## 17 24 3507.093 -17.06217 -592.26862 5.08484965 2914.824 -11.9773191
## 18 25 3507.093 -17.06217 820.49235 10.74000615 4327.585 -6.3221626
## 19 26 3507.093 -17.06217 68.08135 13.08604845 3575.174 -3.9761203
## 20 27 3507.093 -17.06217 94.40121 -9.80731869 3601.494 -26.8694874
## 21 28 3507.093 -17.06217 885.60670 -5.72552251 4392.700 -22.7876913
## 22 29 3507.093 -17.06217 1110.12198 -12.22309142 4617.215 -29.2852602
## 24 30 3507.093 -17.06217 -15.91517 16.21492923 3491.178 -0.8472395
## 25 31 3507.093 -17.06217 -334.43687 14.25731416 3172.656 -2.8048546
## 26 32 3507.093 -17.06217 -233.29116 6.48319622 3273.802 -10.5789725
## 27 33 3507.093 -17.06217 -626.24490 -5.82523977 2880.848 -22.8874085Mean of random effects should be ~= 0
## [1] 7.267125e-15## [1] 2.383972e-12- Damage coefficient for overall model
## [1] -0.49- Relative yield Loss at TSs=50
## [1] 24.5Plot of fitted (and 95%CI) and study specific blups lines
## [[1]]
## NULL
##
## [[2]]
## NULL
3.3 Mixed model - moderator variables (only the significative one)
## Call: mvmeta(formula = outcomes ~ metacov, S = varcov, method = "ml")
##
## Multivariate random-effects meta-regression
## Dimension: 2
## Estimation method: ML
##
## Fixed-effects coefficients
## b0_est :
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## (Intercept) 3382.5607 161.9093 20.8917 0.0000 3065.2243 3699.8971 ***
## metacovLow 248.9863 232.0009 1.0732 0.2832 -205.7272 703.6997
## b1_est :
## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## (Intercept) -25.9173 1.7016 -15.2314 0.0000 -29.2523 -22.5822 ***
## metacovLow 18.6306 2.2910 8.1321 0.0000 14.1404 23.1209 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Between-study random-effects (co)variance components
## Structure: General positive-definite
## Std. Dev Corr
## b0_est 657.6944 b0_est
## b1_est 4.4716 -0.5524
##
## Multivariate Cochran Q-test for residual heterogeneity:
## Q = 10526.1480 (df = 62), p-value = 0.0000
## I-square statistic = 99.4%
##
## 33 studies, 66 observations, 4 fixed and 3 random-effects parameters
## logLik AIC BIC
## -369.0428 752.0857 767.4132## Estimate Std. Error z Pr(>|z|) 95%ci.lb 95%ci.ub
## b0_est.(Intercept) 3382.56073 161.909309 20.891700 0.000000e+00 3065.22431 3699.89714
## b0_est.metacovLow 248.98626 232.000932 1.073212 2.831759e-01 -205.72721 703.69973
## b1_est.(Intercept) -25.91726 1.701569 -15.231388 0.000000e+00 -29.25227 -22.58224
## b1_est.metacovLow 18.63065 2.290993 8.132128 4.440892e-16 14.14038 23.12091## b0_est b1_est
## b0_est 1.0000000 -0.5523883
## b1_est -0.5523883 1.0000000## b0_est b1_est
## b0_est 432561.941 -1624.53412
## b1_est -1624.534 19.99495Improve model0 or over-parametrization?
## df AIC
## model_YR 7 752.0857
## model0 5 794.2754## b0_est b1_est b0_est.1 b1_est.2 b0_est.2 b1_est.1
## 1 3631.547 -7.28661 3305.875 -4.279926 3957.219 -10.29329
## 2 3382.561 -25.91726 3065.224 -22.582243 3699.897 -29.252273.3.1 Relative yield loss for model with RY as moderator
## b0_est b1_est b0_low b1_low b0_up b1_up dam_coef L50
## 1 3631.547 -7.28661 3305.875 -10.29329 3957.219 -4.279926 -0.20 -10.0
## 2 3382.561 -25.91726 3065.224 -29.25227 3699.897 -22.582243 -0.77 -38.5Ouput table
## study Class b0_fitted b1_fitted b0_blup b1_blup b0_random b1_random
## 1 1 Low 3631.547 -7.28661 2854.392 -5.876931 -777.15530 1.40967900
## 2 2 High 3382.561 -25.91726 3986.561 -26.795600 604.00076 -0.87834284
## 3 3 High 3382.561 -25.91726 3481.497 -23.951925 98.93621 1.96533212
## 4 4 Low 3631.547 -7.28661 3466.583 -5.589150 -164.96358 1.69745939
## 5 5 Low 3631.547 -7.28661 4212.545 -9.178439 580.99844 -1.89182923
## 6 6 Low 3631.547 -7.28661 3610.293 -7.227183 -21.25375 0.05942703
## 7 7 Low 3631.547 -7.28661 3737.203 -8.792558 105.65623 -1.50594786
## 8 8 Low 3631.547 -7.28661 3941.126 -6.279310 309.57943 1.00729980
## 9 9 High 3382.561 -25.91726 3780.192 -24.320062 397.63116 1.59719467
## 10 10 High 3382.561 -25.91726 4178.042 -31.764019 795.48143 -5.84676162
## 11 11 High 3382.561 -25.91726 3687.577 -20.405663 305.01639 5.51159415
## 12 12 High 3382.561 -25.91726 1653.034 -19.023039 -1729.52623 6.89421800
## 13 13 High 3382.561 -25.91726 3299.851 -29.400885 -82.70967 -3.48362836
## 14 14 High 3382.561 -25.91726 2053.870 -19.490996 -1328.69120 6.42626147
## 15 15 Low 3631.547 -7.28661 3872.037 -9.728736 240.49032 -2.44212623
## 16 16 Low 3631.547 -7.28661 3555.557 -6.158361 -75.98999 1.12824840
## 17 17 High 3382.561 -25.91726 3745.608 -27.490541 363.04709 -1.57328408
## 18 18 Low 3631.547 -7.28661 3840.426 -7.173270 208.87937 0.11334003
## 19 19 High 3382.561 -25.91726 4395.632 -30.178152 1013.07169 -4.26089480
## 20 20 High 3382.561 -25.91726 3048.987 -28.079063 -333.57383 -2.16180638
## 21 21 High 3382.561 -25.91726 2196.940 -23.111366 -1185.62097 2.80589081
## 22 22 High 3382.561 -25.91726 2970.991 -25.688776 -411.56952 0.22848146
## 23 23 High 3382.561 -25.91726 3845.121 -28.912922 462.56024 -2.99566450
## 24 24 Low 3631.547 -7.28661 2883.324 -10.184356 -748.22278 -2.89774609
## 25 25 Low 3631.547 -7.28661 4336.038 -7.791315 704.49095 -0.50470562
## 26 26 Low 3631.547 -7.28661 3579.760 -4.523661 -51.78687 2.76294855
## 27 27 High 3382.561 -25.91726 3607.815 -27.525757 225.25389 -1.60850038
## 28 28 Low 3631.547 -7.28661 4271.740 -13.505490 640.19307 -6.21888011
## 29 29 High 3382.561 -25.91726 4676.851 -31.081196 1294.29020 -5.16393927
## 30 30 Low 3631.547 -7.28661 3529.438 -4.185128 -102.10911 3.10148203
## 31 31 Low 3631.547 -7.28661 3184.282 -3.233679 -447.26485 4.05293121
## 32 32 Low 3631.547 -7.28661 3230.005 -7.158190 -401.54158 0.12841970
## 33 33 High 3382.561 -25.91726 2894.963 -23.373407 -487.59764 2.54384958## [[1]]
## [1] "range: 1653 to 4676.85"
##
## [[2]]
## [1] "median: 3607.81"
##
## [[3]]
## [1] "IQR: 687.76"## [[1]]
## [1] "range: -31.8 to -3.23"
##
## [[2]]
## [1] "median: -19.02"
##
## [[3]]
## [1] "IQR: 19.62"Plot of fitted (and 95%CI) and study specific blups lines
## [[1]]
## NULL
##
## [[2]]
## NULL## [[1]]
## NULL
##
## [[2]]
## NULL
Histograms
