Robustness and Comparison of Wilks’ Test Statistic for Two-Way MANOVA

The classical Wilks’ statistic is mostly used to test hypotheses in the two way multivariate analysis of variance (MANOVA), which is highly sensitive to the effects of outliers. The non-robustness of test statistics based on normal theory has lead many researchers to study different options. In this paper, we presented a robust version of the Wilks’ test statistic based on highly robust and efficient is reweighted minimum covariance determinant estimates (RMCD). Monte Carlo simulations are used to evaluate the performance of the test statistics under various distributions. In addition, type I error rate results and test power are considered as statistical tools for comparing test statistics .


. Introduction""
Two-way multivariate analysis of variance (MANOVA) deals with testing the effects of the two factors on the measured observations with or without the effects of interaction between the factors.To formalize the two-way multivariate model, it is assumed that two factors (row and column) R with r levels and C with c levels and i j1 i j 2 i jn , ,… , y y y are independent observations of size n have p -variate normal distribution with mean vector ij μ and equal covariance matrix Σ for all i = 1, 2 ,… ,r , j = 1, 2 ,… ,c .The fixed model for two-way MANOVA with interactions is  where μ is an overall effect, i α is the th i level effect of the row factor R , j β is the th j level effect of the factor C , ij γ is the interaction effect between the th i level effect and the th j level effect, and ε is an error random vector that is assumed under classical assumptions to be independent and identically distributed as Σ) , ( N p 0 for all ,, i j and k .Additionally, we have Then, the null and alternative hypotheses for testing the row effects i α , column effects j β , and interaction effects ij γ , can be written respectively as : Four test statistics are used for testing the hypothesis in multivariate case [12] : Wilks' statistic (1932), Lawley -Hotelling statistic (1938)(1939)(1940)(1941)(1942)(1943)(1944)(1945)(1946)(1947)(1948)(1949)(1950)(1951)), Pillai's statistic (1955), and Roy's statistic (1964).
The Wilks' test statistic Λ (1932) is mostly used, which is defined as: E Λ = .
E + H (1.5) where H and E "between" and " within" matrices.
So, the Wilks

ΛΛ
respectively where , and the tabulated values with level of significance α and the degrees of freedom . Alternatively, the two-way MANOVA model without interactions is Using the RC H instead of E in (1.6) and (1.7).The Wilks' statistics for testing the row effects i α and the column effects The null hypotheses of 0R H , and


the tabulated values with level of significance  and the degrees of freedom p , , and Most classical statistics are extremely sensitive to the effect of outliers (Beakman and Cook [3]).Several statistics that are robust against possible outliers in the data have been presented.In 1985, Nath and Pavur [9] presented an alternative statistic for the one-way MANOVA depend on the rank order of the data.In 2010, Todorov and Filzmoser [14] introduced a robust Wilks' statistic for the one-way MANOVA depend on RMCD estimator.Van Aelst and Willems (2011) [15] used S and MM-estimators to construct a robust Wilks' statistic for testing the hypotheses in the one-way MANOVA.Wilcox (2012) [16] suggested a one way MANOVA testing procedure based on trimmed means and Winzorized covariance matrices.Spangl (2018) [13] extended the approach of Todorov and Filzmoser (2010) to two-way MANOVA designs and used robust Wilks' statistics based on the minimum covariance determinant estimator.Ameen and Hadi (2019) [2] presented a robust version of the Wilks' statistic for the one-way MANOVA based on RMCD and constructed its approximate distribution.To perform the robust Wilks' test statistic that proposed by Spangl, it will take a lot of time .Therefore, in this the present study, the same robust Wilks' statistic, which proposed by Spangl is used but with different weight function.The novel of the study is that, construct another an approximate distribution for this statistic.The MCD estimator that proposed by Rousseeuw in (1985) [11] is a highly robust location and scatter estimator that is used for this reason.Reweighted MCD estimator (RMCD) can be used to increase efficiency while maintaining high robustness in section 2. The robust Wilks' statistic is reviewed in section 3.In section 4, we are constructing the approximation proposed and testing its accuracy.A simulation study is used to evaluate the statistical performance proposed and to compare the different test statistics in different distribution cases in terms of significance level, test power and robustness.Section 5 describes the simulation study and its results.

Minimum Covariance Determent (MCD) Estimator
The MCD estimator introduced by Rousseeuw (1985) looks for a subset of h observations with the lowest determinant of the sample covariance matrix, where the subset size h is selected between half and the full size of sample.The mean observations of the subset h represent the MCD location estimate T and a multiple of its covariance matrix is the MCD scatter estimate C .The effective algorithm for calculating the MCD estimates is found in most known statistical software packages such as SAS Plus S R , , 

, and
Matlab .For computing the MCD and the related estimators, the FAST-MCD algorithm of Rousseeuw and Van Driessen (1999) [10] will be used as implemented in the R package.To increase the efficiency of the MCD estimator, a reweighted version is used.
Several methods have been proposed to estimate the common covariance matrix.The method which was introduced by He and Fung (2000) [7] for S estimates and by Hubert and Van Driessen (2004) [8] for MCD estimates is used.In this method, the observations (2.2) In order to increase efficiency while retaining high robustness, one can apply the RMCD estimators.By using the final obtained estimates .ij μ ˆ and Z  we can calculate the Mahalanobis distances as: With these initial robust distances, we can define the weight ijk w for each observation ijk y computed by the Huber weight function defined as:

The Robust Wilks' Statistic
The effect of outliers on the quality of the hypothesis test, which is based on the classical statistics, makes most of the researchers using robust statistics instead of the classical ones.Spangl [13] proposed a robust version of the Wilks' test statistics for two-way MANOVA model with interactions in (1.1) based on the reweighted MCD estimator as: where , ) .

The proposed approximation distribution of Wilks' statistic
The distribution of classical Wilks' statistic  , which was considered by Anderson (1958)  [1] as a ratio of two Wishart distributions, is very complicated.Therefore, the asymptotic approximation based on F or 2  distributions are used .
Bartlett introduced a good approximation based on 2  distribution of the Wilks' statistic given by (see [4]): Todorov and Filzmoser are assumed for the robust Wilks' test statistic Λ R in one way MANOVA the following approximation: , where the multiplication factor d and the degrees of freedom q of the 2  distribution defined as: Spangl extends the approximation of Todorov and Filzmoser (2010) [14] for the robust Wilks' test statistic Λ R in two way MANOVA : The mean

Var
of the approximation L are not possible to obtain analytically.So, they are determined by simulation after repeated  times as follows: For a given dimension p , number of levels c r , , and sample size n in each factor combination group, samples For each sample the robust Wilks' test statistic Λ R based on the weighted MCD will be calculated.After performing 3000  m trials, the sample mean and variance of L will be obtained as: To perform the robust Wilks' statistic that proposed by Spangl , it will take a lot of time during simulations to find d and q for approximate distribution.Therefore, in this the present study, the same robust Wilks' test statistics, which proposed by Spangl is used with different weight function, denoted as the modified robust Wilks' statistic Λ M (i.e.

  R Λ M
).The novel of the study is that, is that, construct another an approximate distribution for this statistic.The matrices R H , C H , RC H , and E in (1.9)-(1.12)can be written as: Where Therefore, the degrees of freedom Analogously we can be written the matrices . where , and ,and , and   Campbell, (1980)   [5]) as:

Monte Carlo Simulation
Monte Carlo study is conducted to evaluate the performance of the proposed test statistics.The evaluation of the performance of any test statistics involves two measures: significance level (the type I error rate) and the power of the test.To compute the simulated : significance level and the power of the test for classical statistics and proposed robust statistics under various distributions( p -variate normal and contaminated distributions), let us consider several dimensions   the significance level α is set to 0.05.

Significance level
To study the attained significance level of the test statistics in two way MANOVA , it is assumed that the observations come from the same multivariate distributions under the null hypotheses of i α , j β and  () LT is the number of times of the test statistic rejects the hypothesis 0 H when 0 H is true, of the test statistics which are above the appropriate critical value of the corresponding approximate distribution are taken as an estimate of the true significance level.The p -value plots proposed by Davidson and MacKinnon (1998) [6] are used ,which give more complete picture of how test statistics follow the approximate distribution under the null hypothesis in the simulated samples.Figure 1 show plots of the empirical distribution functions of the p -values of the two-way MANOVA with and without interactions.we considered several dimensions p ∈ {2,6}, number of levels , rc∈ {2,3,5}, and number of samples n ∈{20,30,50} assuming an equal sample size in each factor combination group.It is seen that the test statistics  , and M  are close to the 45 ° line, and the robust Wilks' statistic R  is considerably below the 45 line for small sample sizes.

Power of test
In order to evaluate the power of the robust Wilks' test statistic, we will generate data samples under an alternative hypothesis.It will incorrectly examine the frequency of failures in rejecting the null hypothesis, i.e. the frequency of type II errors.The same combinations of dimensions p , number of levels ,, rc and sample sizes n in each factor combination group.We're going to differentiate between the two called MANOVA models.Data samples are created from a normal p -dimensional distribution for the two-way MANOVA with interactions, where each () KT the number of times of rejected the test statistic when the hypothesis 0 H is false) of the test statistics when the statistic exceeds its appropriate critical value it will be the estimate of the power for the specific configuration..The results for the sample size n in each factor combination groups are shown in Figure 2. It is clearly seen that the Size-power of two-way MANOVA with and without interaction the classical statistic  and the proposed statistic M  are close while the robust Wilks' statistic R  by Spangl is less.In Figure 5 show power plot of two-way MANOVA with and without interactions the proposed statistic M  and the classical Wilks' Lambda test  are close and while the test R  by Spangl is less.

Robustness comparisons
Now we're going to look at the robustness of the proposed MANOVA two-way test statistic.We will therefore construct data samples under the null and alternative, And by adding outliers, we can contaminate them.The same cases of dimensions , number of levels , rcand sample sizes n in each factor combination group will be used.

Significance level
Here, for each design in Table 1, data samples are generated under empty hypotheses for all combination groups of factors  In these Figure 5 show p -value plot of two-way MANOVA with interactions.Further, Figure 3 show p -value plots (actual size) of two-way MANOVA with and without interactions the p -value plots (actual size) based on the test statistics M  , is so close to the 45 ° line compared to the same of the test statistic R  , while the classical statistic  is very bad for all the different cases of dimension  and sample sizes.Figure 7 show observed type I error rates of the two-way MANOVA with and without interactions in the presence of outliers.Further.Thus, for the M  test, the type I error rates turn out to be quite robust and close to the nominal value.The test R  and the classical Wilks' Lambda test  are seen to be prone to outliers.Both yield very erroneous type I error rates.

Power of test
Under the alternative hypothesis the data samples will be generated from the following contamination model:

Conclusions
In this paper, we presented a robust version of the Wilks' statistics and constructed its approximate distribution.The pvalue and size -power of the new proposed statistics were compared with the classical and robust Wilks' of Spangl statistics in Monte Carlo studies.Various simulations were performed considering different dimensions, number of levels, and sample sizes of factor combination groups.Therefore it can be concluded that p-value plots and size-power curves for the proposed robust statistic are close to the classical and the robust Wilks' statistic of Spangl in case of normal distribution for the data set, while in case of contaminated distribution the proposed robust statistic is the best especially with small sample sizes.Although only a selection of the results is presented in the paper these are typical representative outcomes.


is estimated as the RMCD covariance matrix of the centered observations Z .Finally, the location estimate Z μ of Z is used to adjust the group meansZ Alternatively, the robust Wilks' test statistics for two-way MANOVA model without interactions in (1.13) are:

1 is
observation ijk y computed by the Hample weight function defined (see

2 
approximation of the classical Wilks' statistic Λ in (4.1), we can define for Λ M the following approximations: factor combination group, ij y has a different mean ij  and all of them have the same covariance matrix P I , the data samples are created from a normal p -dimensional distribution for the two-way MANOVA without interactions, where each factor combination group, ij y has a different mean ij  and all of them have the same covariance matrix d again takes the following values: 0.0,0.2,0.5,0.7,1.0,1.5, 2.0 d  .Thus, the classical Wilks' test statistics  , the robust Wilks' test statistics  R ,and the proposed robust Wilks' test statistics M are calculated .This is repeated m = 3000 and the power of test ˆ( Figures 6 show the power plot of two-way MANOVA with and without interaction .It is clearly seen that the proposed robust Wilks' R  statistic is the best compared to the other .r=2 , c=2 , p=2 , n=10

Figure 1 :Figure 2 :Figure 3 :Figure 4 :
Figure 1: p -values plot for test statistics ( ), ( ), ( ) red line R black line and M blue line    of two-way MANOVA with and Figure5: power curves for test statistics ( ), ( ), ( ) red line R black line and M blue line    of two-way MANOVA with and

Figure 6 :Figure 7 :
Figure 7: Type I error for test statistics ( ), ( ), ( ) red line R black line and M blue line    of two-way MANOVA with and without interaction with the same combinations of dimensions p , number of levels , rc , and sample sizes n in each factor combination group of multivariate contaminated distribution .The 45 line is given too.

Table 1 : Selected designs for the simulation study
sample size in each factor combination group.The different designs for two-way MANOVA with and without interactions are given in Table1.