The metan (multi-environment trials analysis) package provides useful functions for analyzing multi-environment trial data using parametric and nonparametric methods, including, but not limited to:

  • Within-environment analysis of variance;
  • Estimation using AMMI considering different numbers of interaction principal component axes;
  • AMMI-based stability indexes;
  • GGE biplot analysis;
  • Prediction in mixed-effect models;
  • BLUP-based stability indexes;
  • Variance components and genetic parameters in mixed-effect models;
  • Cross-validation procedures for AMMI-family and BLUP models;
  • Graphics tools for generating biplots;
  • Parametric and nonparametric stability statistics

For more details see the complete vignette.

Installation

The latest development version can be download from GitHub by running

Getting started

Here, we will use the example dataset data_ge that contains data on two variables assessed in 10 genotypes growing in 14 environments. For more details see ?data_ge. The tables of this vignette were produced with the package KableExtra

AMMI model

Fitting the model

The AMMI model is fitted with the function performs_ammi(). For more details, see the complete vignette.

Predicting the response variable

The S3 method predict() is implemented for objects of class performs_ammi and may be used to estimate the response of each genotype in each environment considering different number of Interaction Principal Component Axis (IPCA). For example, we will use four IPCA (number of significant IPCAs) to estimate the variable GY using the model object. Note that $GY was used because using the predict() function to a model of class waas return a list, with one tbl_df for each variable.

predict(model, naxis = 4)$GY %>% 
  head() %>% 
  print_table()
ENV GEN Y resOLS Ypred ResAMMI YpredAMMI AMMI0
E1 G1 2.366 \-0.084 2.450 0.07115484 2.521273 2.450
E1 G10 1.974 \-0.344 2.318 \-0.35391141 1.963751 2.318
E1 G2 2.902 0.311 2.591 0.29035016 2.880939 2.591
E1 G3 2.889 0.087 2.802 \-0.04518795 2.756598 2.802
E1 G4 2.589 0.100 2.488 0.04942370 2.537781 2.488
E1 G5 2.188 \-0.196 2.384 \-0.07091881 2.312867 2.384

Biplots

ggplot2-based graphics are easily obtained in metan package. For example, the well-known AMMI2 biplot may be obtained as follows. Please, note that since waas() function allows analyzing multiple variables at the same time, e.g., resp = c(v1, v2, ...), the output model is a list, in this case with one element, GY.

GGE model

The GGE model is fitted with the function gge(). For more details, see the complete vignette.

model <- gge(data_ge, ENV, GEN, GY)
model2 <- gge(data_ge, ENV, GEN, GY, svp = "symmetrical")
c <- plot(model)
d <- plot(model2,
          type = 2,
          col.gen = "black",
          col.env = "gray70",
          axis.expand = 1.5)
arrange_ggplot(c, d, labels = letters[3:4])

BLUP model

Linear-mixed effect models to predict the response variable in METs are fitted using the function waasb(). Here we will obtain the predicted means for genotypes in the variables GY and HM. For more details, see the complete vignette.

model2 <- waasb(data_ge, ENV, GEN, REP,
                resp = c(GY, HM),
                verbose = FALSE)

get_model_data(model2, what = "genpar") %>% 
print_table()
Parameters GY HM
GEI variance 0.057 2.189
GEI (%) 31.259 39.665
Genotypic variance 0.028 0.490
Gen (%) 15.447 8.876
Residual variance 0.097 2.840
Res (%) 53.294 51.458
Phenotypic variance 0.181 5.519
Heritability 0.154 0.089
GEIr2 0.313 0.397
Heribatility of means 0.815 0.686
Accuracy 0.903 0.828
rge 0.370 0.435
CVg 6.260 1.456
CVr 11.628 3.504
CV ratio 0.538 0.415

get_model_data(model2, what = "blupg") %>% 
print_table()
gen GY HM
G1 2.617 47.396
G10 2.509 48.375
G2 2.731 47.108
G3 2.903 47.756
G4 2.648 48.050
G5 2.563 48.918
G6 2.560 48.530
G7 2.729 48.004
G8 2.943 48.784
G9 2.541 47.962

Plotting the BLUPs for genotypes

To produce a plot with the predicted means, use the function plot_blup().

e <- plot_blup(model2$GY)
f <- plot_blup(model2$GY,
               prob = 0.1,
               col.shape  =  c("gray20", "gray80")) + coord_flip()
arrange_ggplot(e, f, labels = letters[5:6])

BLUPS for genotype-vs-environment interaction

The object BLUPgge contains the blups for the genotype-vs-environment interaction. In the following example, the values for GY are shown.

ENV GEN BLUPge BLUPg BLUPg+ge Predicted LL UL
E1 G1 \-0.062 \-0.058 \-0.120 2.401 2.298 2.505
E1 G10 \-0.243 \-0.166 \-0.409 2.112 2.009 2.216
E1 G2 0.207 0.057 0.264 2.784 2.681 2.888
E1 G3 0.088 0.229 0.318 2.838 2.735 2.942
E1 G4 0.060 \-0.026 0.034 2.554 2.451 2.658

When more than one variable is fitted, the predicted means for genotype-vs-environment combination may be obtained for all variables in the model using get_model_data().

get_model_data(model2, what = "blupge") %>% 
  head() %>% 
  print_table()
ENV GEN GY HM
E1 G1 2.401 46.587
E1 G10 2.112 47.152
E1 G2 2.784 45.664
E1 G3 2.838 46.249
E1 G4 2.554 48.017
E1 G5 2.268 49.390

Computing parametric and non-parametric stability indexes

stats <- ge_stats(data_ge, ENV, GEN, REP, GY)
get_model_data(stats, "stats") %>% 
  print_table()
var gen Y Var Shukla Wi\_g Wi\_f Wi\_u Ecoval bij Sij R2 ASV SIPC EV ZA WAAS HMGV RPGV HMRPGV Pi\_a Pi\_f Pi\_u Gai S1 S2 S3 S6 N1 N2 N3 N4
GY G1 2.604 10.897 0.028 84.350 89.172 81.062 1.219 1.063 \-0.001 0.966 0.346 0.463 0.015 0.100 0.151 2.316 0.969 0.967 0.169 0.228 0.125 2.449 0.000 8.571 16.581 6.064 2.571 0.367 0.429 0.000
GY G10 2.471 14.237 0.244 59.241 64.605 54.394 7.955 1.122 0.177 0.823 1.225 2.069 0.210 0.437 0.652 2.112 0.913 0.896 0.344 0.475 0.245 2.247 0.022 16.879 32.276 9.655 3.714 0.531 0.577 0.003
GY G2 2.744 11.345 0.086 82.797 95.294 75.616 3.032 1.054 0.050 0.913 0.249 1.544 0.179 0.216 0.283 2.466 1.025 1.020 0.126 0.149 0.108 2.594 0.066 9.604 17.000 4.773 2.429 0.540 0.634 0.014
GY G3 2.955 10.131 0.012 103.522 99.678 106.863 0.725 1.031 \-0.013 0.977 0.113 0.552 0.021 0.080 0.106 2.681 1.105 1.104 0.041 0.072 0.018 2.824 0.022 5.297 1.830 1.525 2.000 0.667 0.862 0.008
GY G4 2.642 8.934 0.064 85.924 79.549 91.944 2.343 0.937 0.030 0.917 0.594 1.036 0.052 0.219 0.326 2.388 0.990 0.988 0.173 0.289 0.085 2.515 0.088 6.269 11.200 4.400 1.929 0.321 0.402 0.015
GY G5 2.537 7.822 0.048 82.717 82.217 82.350 1.844 0.887 0.009 0.937 0.430 0.997 0.043 0.186 0.270 2.302 0.954 0.952 0.240 0.382 0.133 2.422 0.000 8.571 16.804 6.471 2.714 0.339 0.384 0.000
GY G6 2.534 7.339 0.047 82.995 83.668 81.753 1.806 0.861 0.003 0.942 0.265 1.140 0.091 0.172 0.233 2.304 0.954 0.952 0.238 0.377 0.134 2.423 0.022 8.725 13.982 5.228 2.429 0.347 0.411 0.003
GY G7 2.741 7.331 0.122 83.876 77.552 93.387 4.162 0.819 0.058 0.852 0.663 1.787 0.191 0.303 0.428 2.520 1.037 1.031 0.149 0.318 0.021 2.641 0.011 8.308 15.306 4.777 2.286 0.508 0.564 0.002
GY G8 3.004 10.817 0.071 98.799 90.521 107.307 2.569 1.034 0.038 0.922 0.574 1.178 0.067 0.225 0.327 2.740 1.126 1.122 0.041 0.088 0.006 2.877 0.033 6.571 5.732 1.951 2.000 1.000 1.116 0.015
GY G9 2.510 14.745 0.167 68.849 68.891 70.308 5.563 1.192 0.094 0.897 0.983 1.495 0.131 0.336 0.507 2.177 0.926 0.917 0.291 0.390 0.217 2.303 0.066 9.632 18.276 6.414 2.500 0.357 0.436 0.010
get_model_data(stats, "ranks") %>% 
  print_table()
var gen Y\_R Var\_R Shukla\_R Wi\_g\_R Wi\_f\_R Wi\_u\_R Ecoval\_R Sij\_R R2\_R ASV\_R SIPC\_R EV\_R ZA\_R WAAS\_R HMGV\_R RPGV\_R HMRPGV\_R Pi\_a\_R Pi\_f\_R Pi\_u\_R Gai\_R S1\_R S2\_R S3\_R S6\_R N1\_R N2\_R N3\_R N4\_R
GY G1 6 7 2 4 4 7 2 1 2 4 1 1 2 2 6 6 6 5 4 6 6 1.5 5.5 6 7 8.0 5 4 1.5
GY G10 10 9 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 5.0 10.0 10 10 10.0 7 7 4.5
GY G2 3 8 7 7 2 8 7 7 7 2 8 8 5 5 4 4 4 3 3 5 4 8.5 8.0 8 4 5.5 8 8 8.0
GY G3 2 5 1 1 1 2 1 4 1 1 2 2 1 1 2 2 2 1 1 2 2 5.0 1.0 1 1 2.5 9 9 6.0
GY G4 5 4 5 3 7 4 5 5 6 7 4 4 6 6 5 5 5 6 5 4 5 10.0 2.0 3 3 1.0 1 2 9.0
GY G5 7 3 4 8 6 5 4 3 4 5 3 3 4 4 8 7 8 8 8 7 8 1.5 5.5 7 9 9.0 2 1 1.5
GY G6 8 2 3 6 5 6 3 2 3 3 5 6 3 3 7 8 7 7 7 8 7 5.0 7.0 4 6 5.5 3 3 4.5
GY G7 4 1 8 5 8 3 8 8 9 8 9 9 8 8 3 3 3 4 6 3 3 3.0 4.0 5 5 4.0 6 6 3.0
GY G8 1 6 6 2 3 1 6 6 5 6 6 5 7 7 1 1 1 2 2 1 1 7.0 3.0 2 2 2.5 10 10 10.0
GY G9 9 10 9 9 9 9 9 9 8 9 7 7 9 9 9 9 9 9 9 9 9 8.5 9.0 9 8 7.0 4 5 7.0

Getting help

  • If you encounter a clear bug, please file a minimal reproducible example on github

  • Suggestions and criticisms to improve the quality and usability of the package are welcome!