fit <-
bayes_ammi(
.data = Maiz
, .y = y
, .gen = entry
, .env = site
, .rep = rep
, .nIter = 6000
, .k = 2
, .burnin = 1000
, .thin = 1
, .chains = 2
, .priors = bam_priors(
mu_mean = 0
, mu_var = 1e+8
, tau_mean = 0
, tau_var = 1e+8
, delta_mean = 0
, delta_var = 1e+8
, lambda_mean = 0
, lambda_var = 1e+8
, nu = 1
, s2 = NULL
)
, .level = 0.95
, .seed = 4
, .align = TRUE
, .model = "AMMI"
, .het = FALSE
, .waic = TRUE
)\[ \newcommand{\vy}{\mathbf{y}} \newcommand{\vx}{\mathbf{x}} \newcommand{\vb}{\mathbf{b}} \newcommand{\ve}{\mathbf{e}} \newcommand{\vone}{\mathbf{1}} \newcommand{\vbeta}{\boldsymbol{\beta}} \newcommand{\vmu}{\boldsymbol{\mu}} \newcommand{\mX}{\mathbf{X}} \newcommand{\mZ}{\mathbf{Z}} \newcommand{\mV}{\mathbf{V}} \newcommand{\mI}{\mathbf{I}} \newcommand{\mJ}{\mathbf{J}} \newcommand{\mSigma}{\boldsymbol{\Sigma}} \newcommand{\mY}{\mathbf{Y}} \newcommand{\mK}{\mathbf{K}} \newcommand{\mH}{\mathbf{H}} \newcommand{\va}{\boldsymbol{\alpha}} \newcommand{\vg}{\boldsymbol{\gamma}} \newcommand{\vt}{\boldsymbol{\tau}} \newcommand{\vd}{\boldsymbol{\delta}} \newcommand{\vtheta}{\boldsymbol{\theta}} \newcommand{\vv}{\mathbf{v}} \newcommand{\E}{\mathrm{E}} \newcommand{\Var}{\mathrm{Var}} \newcommand{\Cov}{\mathrm{Cov}} \newcommand{\tr}{\operatorname{tr}} \]
1 Abstract
A least-squares AMMI biplot shows point estimates only. Its Bayesian counterpart gives every genotype and environment a highest posterior density (HPD) region. A region that contains the origin marks a stable genotype or environment, regions that overlap mark a homogeneous group, and the marginal interval on the first term decides the direction of the interaction. This vignette shows how biplots() builds the regions and how to read them.
2 The coordinates
For draw \(s\) the point of genotype \(i\) is \(\left(\lambda_1^{(s)1/2}\alpha^{(s)}_{i1}, \lambda_2^{(s)1/2}\alpha^{(s)}_{i2}\right)\). The draws of a genotype form a cloud whose shape shows the uncertainty of its interaction.
3 HPD regions
A bivariate kernel density \(\hat{f}\) is estimated from the cloud. The \(100(1-p)\%\) HPD region is the set where \(\hat{f}\) exceeds the \(p\)-quantile of \(\hat{f}\) at the draws, so it holds \(100(1-p)\%\) of them and is the smallest such set (Hyndman 1996). It contains the origin exactly when \(\hat{f}(\mathbf{0})\) reaches that quantile.
4 Stable genotypes and environments
bp$biplot_data[, list(byplot, name, x, y, stable, axis1_lower, axis1_upper, group)] byplot name x y stable axis1_lower
<char> <char> <num> <num> <lgcl> <num>
1: Genotype 1 -27.39165676 5.76438212 FALSE -38.27946935
2: Genotype 2 -17.88925076 8.40300672 FALSE -28.32124518
3: Genotype 3 -11.78567779 3.78193269 TRUE -22.73533406
4: Genotype 4 37.11735449 0.09401861 FALSE 26.60001160
5: Genotype 5 20.43844153 -9.10673615 FALSE 9.03528012
6: Genotype 6 19.54334049 -5.78161457 FALSE 8.91289231
7: Genotype 7 -8.63141072 9.96597984 TRUE -19.98012491
8: Genotype 8 -34.47493736 -18.02866280 FALSE -45.52039819
9: Genotype 9 23.07379688 4.90769355 FALSE 12.32697873
10: Environment 1 -15.50657921 -5.55082820 FALSE -26.34872405
11: Environment 2 1.17617282 3.68818862 TRUE -9.69787599
12: Environment 3 -17.92256587 -3.66236650 FALSE -28.83775413
13: Environment 4 10.29797179 -1.87419076 TRUE -0.05076488
14: Environment 5 5.90081820 -1.84457566 TRUE -5.09290073
15: Environment 6 -0.37997269 4.32224532 TRUE -10.77625345
16: Environment 7 -5.64668555 -2.18620563 TRUE -16.65112969
17: Environment 8 -46.94161098 8.79385892 FALSE -56.69780534
18: Environment 9 2.67647433 0.83060165 TRUE -8.35287890
19: Environment 10 -15.28337774 -1.34734626 FALSE -25.39295576
20: Environment 11 27.11174574 1.51012117 FALSE 16.85410381
21: Environment 12 11.17082006 6.69505350 TRUE -0.28560818
22: Environment 13 17.28972011 5.07043354 FALSE 5.48043438
23: Environment 14 -0.54832049 0.49671271 TRUE -11.70112079
24: Environment 15 14.55743946 -5.97775500 FALSE 3.28371897
25: Environment 16 5.66452684 -0.03517786 TRUE -4.62622081
26: Environment 17 -0.01626159 -9.16130318 TRUE -11.30940962
27: Environment 18 15.08382629 7.31282769 FALSE 3.72421775
28: Environment 19 -6.06837963 -1.59997864 TRUE -17.33199141
29: Environment 20 -2.61576189 -5.48031542 TRUE -13.60933486
byplot name x y stable axis1_lower
<char> <char> <num> <num> <lgcl> <num>
axis1_upper group
<num> <char>
1: -16.186961 negative
2: -6.547069 negative
3: -1.082809 negative
4: 46.697197 positive
5: 31.512507 positive
6: 30.093179 positive
7: 4.090398 not significant
8: -22.302477 negative
9: 33.828219 positive
10: -4.743675 negative
11: 12.000181 not significant
12: -6.725737 negative
13: 21.354530 not significant
14: 16.566620 not significant
15: 10.981398 not significant
16: 4.768966 not significant
17: -36.925315 negative
18: 13.265884 not significant
19: -4.598291 negative
20: 37.375407 positive
21: 21.643977 not significant
22: 27.665171 positive
23: 9.279971 not significant
24: 25.392452 positive
25: 16.839005 not significant
26: 10.851369 not significant
27: 25.565660 positive
28: 4.525090 not significant
29: 7.968492 not significant
axis1_upper group
<num> <char>-
stable = TRUE: the 95% region contains the origin, so the interaction of the genotype (environment) is not different from zero on the two terms jointly. -
group: the marginal 95% HPD interval of the first coordinate. “positive” and “negative” groups interact in opposite directions; the paper uses this criterion.
A genotype can be stable in the joint sense and still have a first coordinate that differs from zero. This happens when the second term is so uncertain that the region stretches over the origin along the second axis.
5 Only the unstable ones
bp_u <-
biplots(
model = fit
, burnin = 0
, thin = 0.5
, pb = 0.05
, plot_stable = FALSE
, plot_unstable = TRUE
, ncolors = 5
, .which = "both"
, .grid = 100
)
bp_u$biplot
6 Reading genotype and environment together
A genotype yields more than its main effects predict in environments whose first scores have the same sign as its own, and less where the signs differ. The genotypes in the positive group are therefore adapted to the environments in the positive group.
bam_interpret(
x = bp
, higher_better = TRUE
, top = 5
)Biplot with highest posterior density regions
7 of 9 genotypes and 8 of 20 environments have a 95% HPD region that excludes
the origin (significant interaction).
Findings:
- Genotypes whose region contains the origin (interaction not different from
zero, stable): 3 and 7.
- Groups on the first term (marginal HPD interval excludes zero): positive 4,
5, 6 and 9; negative 1, 2, 3 and 8.
- Environments on the first term: positive 11, 13, 15 and 18; negative 1, 3,
8 and 10; the others do not differ from zero.
- Overlapping regions mark genotypes (environments) whose interaction
parameters do not differ: a homogeneous group.
Recommendations:
- Genotypes 4, 5, 6 and 9 are adapted to environments 11, 13, 15 and 18 (both
positive on the first term).
- Genotypes 1, 2, 3 and 8 are adapted to environments 1, 3, 8 and 10 (both
negative on the first term).
Cautions:
- The regions are estimated from the draws by a kernel density; with few
draws their outline is rough.
- The second term is usually estimated with much more uncertainty than the
first: long regions along the second axis.
Sources:
- Crossa, J., Perez-Elizalde, S., Jarquin, D., Cotes, J. M., Viele, K., Liu,
G. and Cornelius, P. L. (2011). Bayesian estimation of the additive main
effects and multiplicative interaction model. Crop Science, 51(4),
1458-1469. doi:10.2135/cropsci2010.06.03437 References
Crossa, J., Perez-Elizalde, S., Jarquin, D., Cotes, J. M., Viele, K., Liu, G. and Cornelius, P. L. (2011). Bayesian estimation of the additive main effects and multiplicative interaction model. Crop Science, 51(4), 1458–1469. doi:10.2135/cropsci2010.06.0343
Hyndman, R. J. (1996). Computing and graphing highest density regions. The American Statistician, 50(2), 120–126.