Abstracts
The additive main effect and multiplicative interaction (AMMI) models allows analysts to detect interactions between rows and columns in a twoway table. However, there are many methods proposed in the literature to determine the number of multiplicative components to include in the AMMI model. These methods typically give different results for any particular data set, so the user needs some guidance as to which methods to use. In this paper we compare four commonly used methods using simulated data based on real experiments, and provide some general recommendations.
twoway table; principal components; simulation
Os modelos de efeitos principais aditivos e interação significativa (AMMI) permitem ao analista detectar interações entre linhas e colunas em uma tabela de dupla entrada. Entretanto, existem muitos métodos propostos na literatura para determinar o número de componentes multiplicativos a incluir nos modelos AMMI. Esses métodos fornecem diferentes resultados para um particular conjunto de dados. Assim, o usuário necessita de orientação sobre quais métodos usar. Nesse trabalho investigamos quatro métodos comumente utilizados em um amplo espectro de condições usando dados simulados baseados em dados de ensaios reais. O método de EastmentKrzanowski apresentou melhor desempenho; cada um dos outros métodos exibiram baixo desempenho em alguma das situações estudadas.
tabela de dupla entrada; componentes principais; simulação
SOILS AND PLANT NUTRITION
Choosing components in the additive main effect and multiplicative interaction (AMMI) models
Escolha de componentes nos modelos de efeitos principais aditivos e interação multiplicativa (AMMI)
Carlos Tadeu dos Santos DiasI, ^{*} * Corresponding author < ctsdias@carpa.esalq.usp.br> ; Wojtek Janusz Krzanowski^{II}
^{I}USP/ESALQ  Depto. Ciências Exatas, Av. Pádua Dias 11, C.P. 09  13418900  Piracicaba, SP  Brasil
^{II}School of Mathematical Sciences, Laver Building, North Park Road, Exeter, EX4 4QE, UK
ABSTRACT
The additive main effect and multiplicative interaction (AMMI) models allows analysts to detect interactions between rows and columns in a twoway table. However, there are many methods proposed in the literature to determine the number of multiplicative components to include in the AMMI model. These methods typically give different results for any particular data set, so the user needs some guidance as to which methods to use. In this paper we compare four commonly used methods using simulated data based on real experiments, and provide some general recommendations.
Key words: twoway table, principal components, simulation
RESUMO
Os modelos de efeitos principais aditivos e interação significativa (AMMI) permitem ao analista detectar interações entre linhas e colunas em uma tabela de dupla entrada. Entretanto, existem muitos métodos propostos na literatura para determinar o número de componentes multiplicativos a incluir nos modelos AMMI. Esses métodos fornecem diferentes resultados para um particular conjunto de dados. Assim, o usuário necessita de orientação sobre quais métodos usar. Nesse trabalho investigamos quatro métodos comumente utilizados em um amplo espectro de condições usando dados simulados baseados em dados de ensaios reais. O método de EastmentKrzanowski apresentou melhor desempenho; cada um dos outros métodos exibiram baixo desempenho em alguma das situações estudadas.
Palavraschave: tabela de dupla entrada, componentes principais, simulação
INTRODUCTION
There are many tests in the literature for deciding how many components should be retained in additive and multiplicative interaction models. They include those by Gollob (1968), Mandel (1969; 1971), Corsten & van Eijnsbergen (1972), Johnson & Graybill (1972), Hegemann & Johnson (1976), Yochmowitz & Cornell (1978), and Cornelius (1993), among others. Data for analysis are assumed to come in the form of a twoway table. The classical model for such a table is the ANOVA model
i=1, ,g, j=1, ,e, where y_{ij} is the observation in the ith row and jth column of the table, µ is the overall mean effect, g_{i} represents the ith row effect, and e_{j} the jth column effect. The terms d_{ij} can either be seen as residuals or as interaction terms between rows and columns. The above model is called the additive model. However, it is quite possible that the interaction terms d_{ij} still contain some structure, and therefore one can model them by a set of multiplicative components plus residual error:
The number of components m < min {g1, e1} should be chosen in such a way that the residual e_{ij} represents white noise. Combining the last two expressions above yields an additive main effect and multiplicative interaction (AMMI) model for the twoway table:
To estimate the unknown parameters in the model, one usually first uses row/column means for the main effects and then performs a singular value decomposition of the residual matrix for the interaction parameters. This classical approach corresponds essentially to a least squares fit of the full model.
Nonadditive effects are frequently observed in twoway tables, and, as Daniel (1976) points out, the nonadditivity is often associated with just a few rows or columns of the table. As observed by Snee (1982), the interpretation of nonadditivity is less of a problem if replicate observations are present for each of the cells of the table. However when there is only one observation per cell, it is not possible to distinguish between row or columnrelated nonhomogeneous variance, and interaction from the observed data alone. In this situation, the use of a model that will detect a variety of different types of nonadditivity in twoway tables can be very helpful.
For this reason, a primary objective is to determine the optimal number of multiplicative terms in the AMMI model. Also, good prediction of the true trait response in each cell of the twoway table can be achieved by truncating the AMMI model, so criteria for determining the number of components needed to explain the patterns of interactions have been the objects of some research.
Dias & Krzanowski (2003) described two 'leaveoneout' methods to determine the number m of components to be retained in the AMMI model, and two other methods based on Ftests: the Gollob and Cornelius tests. However, these methods gave different results when applied to several data sets, so some guidance is needed regarding their properties if a choice between them is to be made in practice. In the present paper we conduct a comparison of the methods, to establish conditions under which each method might be expected to perform satisfactorily, by means of a simulation study built round real experimental data.
MATERIAL AND METHODS
Data and Model
Model and examples are described in terms of a set of g genotypes that have been tested experimentally in e environments. The mean yield for each combination of genotype and environment, obtained from n replications of an experiment (a balanced set of data) can be represented by the array Y_{(gxe)}=[y_{ij}].
The AMMI model postulates additive components for the main effects of genotypes (g_{i}) and environments (e_{j}), and multiplicative components for the effect of the interaction (ge)_{ij}, referred to as genotype x environment interaction (GEI). Thus, the mean response of genotype i in an environment j is modelled by:
in which (ge)_{ij} is represented by:
under the restrictions:
Estimates of the overall mean (µ) and the main effects (g_{i} and e_{j}) are obtained in the context of a simple twoway ANOVA of the array of means Y_{(gxe)}. The residuals from this array then constitute the array of interactions GE_{(gxe)}=[(ge)_{ij}], and the multiplicative interaction terms are estimated from the singular value decomposition (SVD) of this array. Thus, l_{k} is estimated by the kth singular value of GE, a_{ik} is estimated by ith element of the left singular vector a_{k(gx1)}, and g_{jk} is estimated by jth element of the right singular vector g^{T}_{k(1xe)} associated with l_{k} (Good, 1969; Mandel, 1971; Piepho, 1995). If the interaction contains one or more terms (m > 1), then it is assumed that the errors have a homogeneous variance after the nonadditivity effects have been taken into account.
Snee (1982) points out that nonadditivity in a twoway table, with one observation per cell, may be due either to row or columnrelated nonhomogeneous variance, or to interaction between the row and column factors. He discusses one example involving the effect of four varieties of wheat in 13 locations, and shows that the varietyrelated nonhomogeneous variance detected by several authors is more usefully interpreted as a leveldependent interaction due to the differences among the yields of different varieties of wheat being larger in those locations that have higher yields. Also, he comments that experience has indicated that frequently only one multiplicative term is needed to describe the nonadditivity.
Methods for selecting number of interaction components
Dias & Krzanowski (2003) present two methods based on a full 'leaveoneout' procedure that optimizes the crossvalidation process (i.e. maximises the number of data points left in the set at each iteration without incurring bias due to resubstitution), and two methods based on Ftests. The Eastment & Krzanowski (1982) and Gabriel (2002) crossvalidation methods assume that we wish to predict the elements x_{ij} of the matrix X by means of the model: . The methods are those outlined by Krzanowski (1987) and Gabriel (2002) respectively, in which we predict the value of x_{ij} (i=1, ,g; j=1, ,e) for each possible choice of m (the number of components), and measure the discrepancy between actual and predicted values as
However, to avoid bias in the prediction, the data point x_{ij} must not be used in the calculation of for each i and j. Hence, appeal to some form of crossvalidation is indicated, and the two approaches differ in the way that they handle this. However, both assume that the SVD of X can be written as X=UDV^{T}. The EastmentKrzanowski method utilizes the following statistic to determine the number of components in the model:
where D_{m} is the number of degrees of freedom required to fit the mth component and D_{r} is the number of degrees of freedom remaining after fitting the mth component. The standard crossvalidation procedure is to subdivide X into a number of groups, delete each group in turn from the data, evaluate the parameters of the predictor from the remaining data, and predict the deleted values (Wold, 1976; 1978). The most precise prediction results when each deleted group is as small as possible, and in the present instance that means a single element of X (Krzanowski, 1987). Denoting by X^{(i)} the result of deleting the ith row of X and meancentering the columns, and by X_{(j)} the result of deleting the jth column of X and mean centering the rows, following the scheme given by Eastment & Krzanowski (1982), one can write
and
Now consider the predictor
Each element on the righthand side of this equation is obtained from SVD of X, meancentered after omitting either the ith row or the jth column. Thus, the value x_{ij} has nowhere been used in calculating the prediction, and maximum use has been made of the other elements of X. Calculations here are exact, so there is no problem with convergence as opposed to expectation maximization approaches that have also been applied to AMMI, but are not guaranteed to converge.
On the other hand, the Gabriel method takes a mixture of regression and lowerrank approximation and utilizes PRESS(m) and to achieve the same goal. The proposed algorithm for crossvalidation of lower rank approximations is as follows.
For given (GEI) matrix X, use the partition
and approximate the submatrix X_{\11} by its rank m fit using the SVD
where U = [u_{1},..., u_{m}], V = [v_{1},..., v_{m}], and D = diag[d_{1},..., d_{m}].
Then predict x_{11} by
and obtain the crossvalidation residual .
Similarly obtain the crossvalidation fitted values and residuals for all other elements x_{ij}, i=1, ,g, j=1, ,m; (i,j) ¹ (1,1). Each will require a different partition of X. These residuals and fitted values can be summarized by PRESS(m) and PRECORR(m) respectively.
With either method, choice of m can be based on some suitable function of PRESS(m). However, the features of this statistic differ for the two methods. Gabriel's approach yields values that first decrease and then (usually) increase with m. He therefore suggests that the optimum value of m is the one that yields the minimum of the function. The EastmentKrzanowski approach produces (generally) a set of values that is monotonically decreasing with m. The optimum value for m is the highest value of m at which the statistic W is greater than unity.
The Gollob (1968) approximate Ftest assumes that is distributed as a chisquare variable, and he suggests using the statistic
against an Fdistribution with f_{1}=g+e1(2m) and ge(n1) degrees of freedom, or g(e1)(n1) degrees of freedom if blocks are present, to test the mth multiplicative term of the model for significance. By contrast, Cornelius et al. (1992) uses
against an Fdistribution with f_{2}=(g1m)(e1m) and ge(n1) degrees of freedom, or g(e1)(n1) degrees of freedom if blocks are present, to test significance of lack of fit of the mterm model.
Piepho (1995) investigated the robustness (to the assumptions of homogeneity and normality of the errors) of some alternative tests to select an AMMI model. He comments that F tests applied in accordance with Gollob's (1968) criterion are liberal in that they select too many multiplicative terms. Of the four methods he studied, including that of Gollob (1968), the test proposed by Cornelius et al. (1992) was the most robust. The author thus recommends that preliminary evaluations should be conducted to verify the validity of the assumptions if one of the other tests is to be used. In our study, the tests of Gollob and Cornelius were performed at a = 0.05. We note that other possible tests of both PRESS and F form exist (e.g. the PRESS test in Cornelius et al., 1996, and the F test based on the GoodmanHaberman theorem described by Cornelius, 1993), but we have not included them in the comparisons to keep the investigation to a manageable size.
Simulation Study
We studied by simulation the behaviour of the above methods of selecting the number of interaction components in AMMI models. However, to make the simulations realistic, we based them around real data sets for which the model is appropriate, according to the following scheme.
Suppose that Y = [y_{ij}] is the twoway table for a chosen data set, and that we fit the AMMI model
for a given value of m. Let an overbar denote average, with '.' indicating the subscripts over which summation takes place. Therefore,, and a_{ik}l_{k}g_{jk} is estimated from the kth component of the SVD of . The model assumes e_{ij} ~ N(0,s^{2}), and s^{2} is estimated from the error sum of squares in the ANOVA of the data set.
To generate a data set in the simulation study for a chosen value of m, we simply form the values
where the are independently drawn from a distribution. By these means it is possible to generate repeated data sets which mimic the patterns of the original data, but which have built in to them a specified number of interaction components (governed by the initial choice of m). Also, by varying some of the parameter values, it is possible to force particular structures on the data. For example, by changing the values of we can force some components to have a greater or smaller contribution to the interaction. We used both of these devices in the simulations reported below.
Wheat Yield Data
Our main study uses as a basis the wheat yield data, Trial 3, given by Cornelius & Crossa (1999) who describe five multienvironment (MET) international cultivar trials, all set in randomized block design. Trial 3 concerns maize (Zea mays L.), with 9 cultivars, 20 sites and n = 4.
The first series of simulation experiments was conducted to investigate the sensitivity of the methods to increasing numbers of interaction components in the data. We set the number of interaction components m successively equal to 0, 1, 2, ,9 and then used all estimated parameters from the original data, and the scheme outlined above, to produce the simulated data. We generated 100 repetitions for each of these 10 experiments, and used the four methods described above to determine in each case how many interaction components should be included in the model. Table 1 reports the value chosen most frequently for each method and each experiment.
We then conducted three further series of experiments, constraining the emphasis of components in each case as suggested in the previous section. In the first series of experiments, we arranged for each of the nine interaction components to have approximately equal weight. To do this we arranged that values had the following pattern: but with c chosen such that equalled the original GEI sum of squares and Once again 100 repetitions were provided for each experiment, and modal choices of number of interaction components are presented for each method and each experiment in Table 2.
In the second series of experiments, we inflated the first three to have high values and reduced the remaining , but in proportion to their actual values so that was again equal to the original sum of squares of GEI. In the third series of experiments, the first three values were again inflated while the remaining values were reduced, but the latter were now set to be approximately constant. Results from these two series are given in Tables 3 and 4 respectively.
As a final study we used the data of Trial 5 from Cornelius & Crossa (1999) as the basis for our simulations. This was a simple trial, and considered by the original author to be relatively free of interaction effects. Here we have 8 cultivars and 59 sites. We therefore only conducted the unconstrained form of experiments (i.e. as in Table 1), and show the results in Table 5.
DISCUSSION
We admit at the outset that our investigations are by no means exhaustive (a comparison of just four methods on a relatively small set of conditions), and are not designed in such a way that 'optimal' methods can be determined with respect to some welldefined mathematical objective function. Instead, we start from the view that there may often be conflict in practice among the methods, and that what is required by the practitioner is some reassurance that a chosen method is likely to perform 'reasonably'. So our conclusions in this discussion are essentially qualitative in nature, measured against what we believe to be reasonable interpretations of each studied situation.
The four methods of selecting the number of multiplicative interaction components have indeed yielded different results on the simulated data sets used here. Moreover, the results for these tests give conflicting recommendations on a particular m value for a data set. This feature can be seen clearly in Tables 15.
In general terms, the two methods based on Ftests are conservative and are prone to select a fairly high value of m. The EastmentKrzanowski method, on the other hand, is clearly the most parsimonious: it never selects a value of m greater than 3, and usually recommends either m = 1 or m = 2. We note that this accords closely with the view of Snee (1982) as regards an adequate description of nonadditivity.
Passing to particular results, consider the first experiment. The number of interaction components chosen by each method from a single analysis of the original data (as reported by Dias & Krzanowski, 2003) is given at the foot of Table 1: all methods suggest two components, except Gabriel's method which suggests seven components. Table 1 then shows the modal number of components chosen by each method as progressively more interaction terms are included in the simulated data. Since there appear to be only about two interaction components in the original data, and since the components are added according to their ranking in order of decreasing values (last column of Table 1), we would expect the addition of the third and succeeding components to have little effect. The four methods show different patterns of behaviour. The EastmentKrzanowski method is very stable, choosing an m value of 2 for the original data compared with 1 or 2 for all simulations. Gabriel's and Gollob's methods are very volatile, the former contrasting m = 7 for the original data with values 1, 2, 3 or 4 for the simulations, and the latter having all values between one and five in the simulations, but m = 2 for the original data. Cornelius' method is less volatile but still exhibits a range of values from 1 to 4.
Consider next the simulations with constrained interaction components. In Tables 3 and 4, where the interaction components have been constrained by heavily weighting the first three components, it is probably reasonable to say that no method should choose a value greater than 3. The methods of EastmentKrzanowski and Cornelius accord with this, while those of Gabriel and Gollob recommend a value of 4 for a substantial proportion of cases. However, in both cases we would expect the highest value reached to then stay constant as more components are added. Here the EastmentKrzanowski method fails in Table 4, where m drops back down from 3 to 1 from M = 5 onwards.
In Table 2 the components have a fairly even weighting, so perhaps we would expect there to be a monotonically increasing value of m as M rises. All methods apart from that of EastmentKrzanowski do indeed exhibit such an increasing set of values, whereas this method again drops back from 3 to 1 at M = 4 onwards.
Finally, Table 5 reports a case where terms are successively added in the simulations as m increases, but where the originator of the data believed there to be little interaction. The very even distribution and low values of the squared eigenvalues (last column of Table 5) support this view. The EastmentKrzanowski method suggests just one component (the minimum value possible) for all M values, while all other methods exhibit an upward trend in number of components chosen as m increases, and select very large values for at least some of these cases. Tables 15 do not present a modal number of AMMI multiplicative terms for the Cornelius method from m equal to eight because of negative degrees of freedom.
In summary, therefore, it is a rather mixed picture. The EastmentKrzanowski method is stable and behaves appropriately in the presence of a relatively small number of 'important' components, but generally underestimates when there is a larger number. The Gabriel method is very volatile and tends to choose slightly too many components in many situations. Of the Ftest methods, the one by Cornelius exhibits appropriate behaviour in the presence of 'important' components, but is less stable than EastmentKrzanowski. Finally, the Gollob version is broadly similar to the Cornelius method but with slightly worse stability, and is likely to choose more components in some situations. So, finally, the EastmentKrzanowski method appears to be the preferable crossvalidatory one while the Cornelius method would be the recommended Ftest method. The former is preferable to the latter if parsimony is a major consideration, but the latter is preferable if the presence of a large number of interaction components is suspected. Given the different performances of the two methods, however, both should perhaps be used before making a decision in a practical application. These conclusions support those arrived at by Dias & Krzanowski (2003) from the analysis of real data.
ACKNOWLEDGEMENTS
We are grateful to the anonymous reviewers, whose very thoughtful comments have considerably improved the manuscript. This research was financially supported by a CNPq PD grant to the first author.
Received February 18, 2005
Accepted February 28, 2006
 CORNELIUS, P.L. Statistical tests and retention of terms in the additive main effects and multiplicative interaction model for cultivar trials. Crop Science, v.33, p.11861193, 1993.
 CORNELIUS, P.L.; CROSSA, J Prediction assessment of shrinkage estimators of multiplicative model for multienvironment cultivar trials. Crop Science, v.39, p.9981009, 1999.
 CORNELIUS, P.L.; SEYEDSADR, M.; CROSSA, J. Using the shifted multiplicative model to search for "separability" in crop cultivar trials. Theoretical and Applied Genetics, v.84, p.161172, 1992.
 CORNELIUS, P.L.; CROSSA, J. SEYEDSADR, M. Statistical tests and estimators of multiplicative models for genotypebyenvironment interaction. In: KANG, M.S.; GAUCH, H.G. (Ed.) Genotypebyenvironment interaction. Boca Raton: CRC Press, 1996. p.199234.
 CORSTEN, L.C.A.; VAN EIJNSBERGEN, C.A. Multiplicative effects in twoway analysis of variance. Statistica Neerlandica, v.26, p.6168, 1972.
 DANIEL, C. Application of statistics to industrial experimentation New York: John Wiley, 1976.
 DIAS, C.T.S.; KRZANOWSKI, W.J. Model selection and crossvalidation in additive main effect and multiplicative interaction (AMMI) models. Crop Science, v.43, p.865873, 2003.
 EASTMENT, H.T.; KRZANOWSKI, W.J. Crossvalidatory choice of the number of components from a principal component analysis. Technometrics, v.24, p.7377, 1982.
 GABRIEL, K.R. Le biplot outil d'exploration de données multidimensionelles. Journal de la Societe Francaise de Statistique, v.143, p.555, 2002.
 GOLLOB, H.F. A statistical model which combines features of factor analytic and analysis of variance techniques. Psychometrika, v.33, p.73115, 1968.
 GOOD, I.J. Some applications of the singular decomposition of a matrix. Technometrics, v.11, p.823831, 1969.
 HEGEMANN, V.; JOHNSON, D.E. On analyzing twoway AOV with interaction. Technometrics, v.18, p.273281, 1976.
 JOHNSON, D.E.; GRAYBILL, F.A. An analysis of a twoway model with interaction and no replication. Journal of the American Statistical Association, v.67, p.862868, 1972.
 KRZANOWSKI, W.J. Crossvalidation in principal component analysis. Biometrics, v.43, p.575584, 1987.
 MANDEL, J. The partitioning of interaction in analysis of variance. Journal of Research of the National Bureau of Standards B. Mathematical Sciences, v.73, p.309328, 1969.
 MANDEL, J. A new analysis of variance model for nonadditive data. Technometrics, v.13, p.118, 1971.
 PIEPHO, H.P. Robustness of statistical test for multiplicative terms in additive main effects and multiplicative interaction model for cultivar trial. Theoretical and Applied Genetics, v.90, p.438443, 1995.
 SNEE, R.D. Nonadditivity in a twoway classification: is it interaction or nonhomogeneous variance? Journal of American Statistical Association, v.77, p.515519, 1982.
 WOLD, S. Pattern recognition by means of disjoint principal component models. Pattern Recognition, v.8, p.127139, 1976.
 WOLD, S. Crossvalidatory estimation of the number of components in factor and principal component models. Technometrics, v.20, p.397405, 1978.
 YOCHMOWITZ, M.G.; CORNELL, R.G. Stepwise tests for multiplicative components of interaction. Technometrics, v.20, p.7984, 1978.
Publication Dates

Publication in this collection
17 Apr 2006 
Date of issue
Apr 2006
History

Accepted
28 Feb 2006 
Received
18 Feb 2005