Skip to main content
  • Methodology article
  • Open access
  • Published:

Weighted analysis of general microarray experiments

Abstract

Background

In DNA microarray experiments, measurements from different biological samples are often assumed to be independent and to have identical variance. For many datasets these assumptions have been shown to be invalid and typically lead to too optimistic p-values. A method called WAME has been proposed where a variance is estimated for each sample and a covariance is estimated for each pair of samples. The current version of WAME is, however, limited to experiments with paired design, e.g. two-channel microarrays.

Results

The WAME procedure is extended to general microarray experiments, making it capable of handling both one- and two-channel datasets. Two public one-channel datasets are analysed and WAME detects both unequal variances and correlations. WAME is compared to other common methods: fold-change ranking, ordinary linear model with t-tests, LIMMA and weighted LIMMA. The p-value distributions are shown to differ greatly between the examined methods. In a resampling-based simulation study, the p-values generated by WAME are found to be substantially more correct than the alternatives when a relatively small proportion of the genes is regulated. WAME is also shown to have higher power than the other methods. WAME is available as an R-package.

Conclusion

The WAME procedure is generalized and the limitation to paired-design microarray datasets is removed. The examined other methods produce invalid p-values in many cases, while WAME is shown to produce essentially valid p-values when a relatively small proportion of genes is regulated. WAME is also shown to have higher power than the examined alternative methods.

Background

The DNA microarray technique involves a series of steps, from the harvesting of cells or biopsies to the preprocessing of the scanned arrays, before analysable data are obtained. During several of these steps the quality can be affected by random factors. For instance, depending on the handling of a biological sample the mRNA can be more or less degraded [1], and the cell-type composition of a biopsy can be more or less representative for the tissue in question. When arrays share sources of variation the deviations from the nominal value will be correlated. For example, two arrays from sources with degraded RNA will both tend to underestimate the expression of easily degradable genes, and two biopsies with a similar and non-representative cell-type composition will deviate in a similar fashion from the average expression for the ideal cell-type composition.

The procedure Weighted Analysis of Microarray Experiments (WAME) [2, 3] introduced a model where a covariance-structure matrix common for all genes aims at catching differences in quality by differences in variances and covarying deviations by correlations between arrays. For computations of test statistics and estimators this resulted in weighting of observations according to the estimated covariance-structure matrix, giving lower weight to imprecise or positively correlated arrays.

In order for the estimation of the covariance matrix to work in the current WAME method, the measurements of most genes must only measure noise, i.e. have an expected value of zero. This is the case in experiments where pair-wise log-ratios are observed and where few genes are differentially expressed between any of the pairwise measured conditions. In the present paper, this crucial constraint will be relaxed to only require that most genes are non-differentially expressed between the conditions actually being compared. Thus, non-paired experiments can be analysed, e.g. many additional ones based on one-channel microarray data. The relaxation is realised by transforming the data to remove irrelevant information in a manner yielding transformed data with expectation zero for non-differentially expressed genes, after which the current WAME method is applied. The transformed data are shown to give equivalent tests and estimates to those of the original data, given the corresponding covariance-structure matrices.

Problem formulation and current methods

Given a microarray experiment with n arrays and m genes, we observe for each gene g an n-dimensional vector X g of log2 transformed values measuring mRNA abundance. In WAME the vector X g is assumed to have expectation μ g described by a design matrix D and a gene-specific parameter vector γ g , typically having one dimension per studied condition. A covariance-structure matrix Σ, common for all genes, is used to model differences in quality between arrays as different variances and shared sources of variation between arrays as correlations. A gene-specific variance-scaling factor c g is assumed to have inverse gamma prior distribution with a global shape parameter α. Conditional on c g the vector X g is assumed to have a normal distribution with covariance matrix c g Σ. A matrix C specifies the differential expression vector δ g , describing the linear combinations of the parameters that are of main interest. Formally,

μ g = D γ g , X g | c g ~ N ( μ g , c g Σ ) , c g ~ Γ 1 ( α , 1 ) , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqadeWadaaabaaccmGae8hVd02aaSbaaSqaaiabdEgaNbqabaaakeaacqGH9aqpaeaacqWGebarcqWFZoWzdaWgaaWcbaGaem4zaCgabeaakiabcYcaSaqaaGqabiab+HfaynaaBaaaleaacqWGNbWzaeqaaOGaeiiFaWNaem4yam2aaSbaaSqaaiabdEgaNbqabaaakeaacqGG+bGFaeaacqqGobGtcqGGOaakcqWF8oqBdaWgaaWcbaGaem4zaCgabeaakiabcYcaSiabdogaJnaaBaaaleaacqWGNbWzaeqaaOGaeu4OdmLaeiykaKIaeiilaWcabaGaem4yam2aaSbaaSqaaiabdEgaNbqabaaakeaacqGG+bGFaeaacqqHtoWrdaahaaWcbeqaaiabgkHiTiabigdaXaaakiabcIcaOGGaciab9f7aHjabcYcaSiabigdaXiabcMcaPiabcYcaSaaaaaa@58F6@
(1)

and variables corresponding to different genes are assumed independent. We want to estimate the differential expressionδ g = C γ g

or we want to test for differential expression

H 0 : δ g = 0 H A : δ g 0 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqadeGadaaabaGaemisaG0aaSbaaSqaaiabicdaWaqabaGccqGG6aGoiiWacqWF0oazdaWgaaWcbaGaem4zaCgabeaaaOqaaiabg2da9aqaaGqabiab+bdaWaqaaiabdIeainaaBaaaleaacqWGbbqqaeqaaOGaeiOoaOJae8hTdq2aaSbaaSqaaiabdEgaNbqabaaakeaacqGHGjsUaeaacqGFWaamcqGGUaGlaaaaaa@3F3A@
(3)

In the current version of WAME [2, 3] the estimation of the covariance-structure matrix Σ is based on a temporary assumption of expectation zero, μ g = 0, for all genes, which is shown to give reasonable results if the expectation is close to zero for most genes. Thus, this is a suitable assumption for data with paired observations and few regulated genes between the pair-wise measured conditions.

The WAME model can be compared with the ordinary linear model (OLM) [4],X g ~ N(μ g , c g I)

which gives rise to the ordinary t- or F-tests, and with a widely used empirical Bayes model proposed in [5] and implemented in the LIMMA package [6],

X g | c g ~ N ( μ g , c g I ) , c g ~ Γ 1 ( α , β ) . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqadeGadaaabaacbeGae8hwaG1aaSbaaSqaaiabdEgaNbqabaGccqGG8baFcqWGJbWydaWgaaWcbaGaem4zaCgabeaaaOqaaiabc6ha+bqaaiabb6eaojabcIcaOGGadiab+X7aTnaaBaaaleaacqWGNbWzaeqaaOGaeiilaWIaem4yam2aaSbaaSqaaiabdEgaNbqabaGccqWGjbqscqGGPaqkcqGGSaalaeaacqWGJbWydaWgaaWcbaGaem4zaCgabeaaaOqaaiabc6ha+bqaaiabfo5ahnaaCaaaleqabaGaeyOeI0IaeGymaedaaOGaeiikaGccciGae0xSdeMaeiilaWIae0NSdiMaeiykaKIaeiOla4caaaaa@4FD3@
(5)

The novel feature of WAME was thus the introduction of the quality modelling covariance-structure matrix Σ.

After the introduction of WAME, a weighted version of LIMMA was proposed [7], which we will refer to as wLIMMA. There, a model with array-wise variance scales but no correlations is used,

X g | c g ~ N ( μ g , c g diag ( σ 1 2 , ... , σ n 2 ) ) , c g ~ Γ 1 ( α , β ) . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqadeGadaaabaacbeGae8hwaG1aaSbaaSqaaiabdEgaNbqabaGccqGG8baFcqWGJbWydaWgaaWcbaGaem4zaCgabeaaaOqaaiabc6ha+bqaaiabb6eaojabcIcaOGGadiab+X7aTnaaBaaaleaacqWGNbWzaeqaaOGaeiilaWIaem4yam2aaSbaaSqaaiabdEgaNbqabaGccqqGKbazcqqGPbqAcqqGHbqycqqGNbWzcqGGOaakiiGacqqFdpWCdaqhaaWcbaGaeGymaedabaGaeGOmaidaaOGaeiilaWIaeiOla4IaeiOla4IaeiOla4IaeiilaWIae03Wdm3aa0baaSqaaiabd6gaUbqaaiabikdaYaaakiabcMcaPiabcMcaPiabcYcaSaqaaiabdogaJnaaBaaaleaacqWGNbWzaeqaaaGcbaGaeiOFa4habaGaeu4KdC0aaWbaaSqabeaacqGHsislcqaIXaqmaaGccqGGOaakcqqFXoqycqGGSaalcqqFYoGycqGGPaqkcqGGUaGlaaaaaa@623B@
(6)

The parameters are estimated using a restricted maximum-likelihood (REML) approach.

A widely used approach is to only consider the ordinary least-squares estimated differential expression, often referred to as the log fold-change, here abbreviated as FC, or as the average M-value. In the present paper, the ranking of the genes imposed by this method will be included in comparisons, when applicable.

Results

The new version of WAME

In the current version of WAME [2, 3] the covariance-structure matrix Σ is estimated using a temporary assumption that μ g = 0 for most genes, i.e. that the measurements of most genes consist solely of biological and technical noise. In the new version of WAME we relax this to only assume that most genes are non-differentially expressed, i.e. δ g = 0. This allows a much larger class of experimental designs and design matrices D, most notably unpaired designs.

The trick used is to transform the data and consider

Y g = X g μ ˜ g 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaieqacqWFzbqwdaWgaaWcbaGaem4zaCgabeaakiabg2da9iab=HfaynaaBaaaleaacqWGNbWzaeqaaOGaeyOeI0cccmGaf4hVd0MbaGaadaqhaaWcbaGaem4zaCgabaGaeGimaadaaaaa@386D@
(7)

where μ ˜ g 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaDaaaleaacqWGNbWzaeaacqaIWaamaaaaaa@30EB@ is a suitable linear estimator of μ g which is unbiased under H0 and which preserves the estimability of the differential expression δ g , based on only the transformed data (see Methods for details). An example is (8) below where for each gene the mean value of all arrays is subtracted.

Since the transformed data contain only noise for non-differentially expressed genes by construction, the current version of WAME can essentially be applied to the transformed data Y g . As before, the covariance-structure matrix (now Σ Y ) and the hyperparameter α are first estimated under a provisional assumption (now δ g = 0). The maximum likelihood estimates of δ g and the likelihood ratio test statistics of (3) are then computed. The tests and estimators are in fact unchanged by the transformation (7), if the covariance-structure matrices for the transformed and untransformed data are known (details given in Methods). WAME is implemented as a package for the R language [8] and is available at [9].

Evaluation on real and resampled data

To investigate the properties of the new version of WAME, two real datasets are examined. Briefly, they are analysed both using WAME and the current methods described in Background. Array-specific weights, p-value distributions and rankings are produced showing clear differences between the procedures, most notably in the p-value distributions. To investigate the power of the different procedures and to look at p-value distributions in a controlled but realistic setting, we also analyse simulated data with real noise from the studied datasets and synthetic signal.

Description of the real datasets

Two public one-channel microarray datasets are analysed. The datasets are selected from the NCBI GEO database [10] with the criteria of having unpaired design and being sufficiently large to allow for the resample-based simulations in Resampled data below.

In the first dataset [11], biopsies were taken from the left atrium from 20 human hearts with normal sinus rythm and 10 hearts with permanent atrial fibrillation. It is here referred to as Atrium. In the second dataset [12], mechanisms in chronic obstructive pulmonary disease, COPD, were investigated by taking lung tissue biopsies from 12 smokers with mild or no emphysema and from 18 smokers with severe emphysema. In both datasets one Affymetrix HGU-133A array was used for each patient. In the present paper RMA [13] is used to obtain expression measures from the raw probe-wise intensities. The analyses are performed using the R language and the Bioconductor framework [14].

Analysis of the real datasets

A natural parameterisation of the included datasets is to have one parameter per condition, yielding design and hypothesis matrices

D = [ 1 0 1 0 0 1 0 1 ]  and  C = [ 1 1 ] . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGebarcqGH9aqpdaWadaqaauaabeqagiaaaaqaaiabigdaXaqaaiabicdaWaqaaiabl6Uinbqaaiabl6UinbqaaiabigdaXaqaaiabicdaWaqaaiabicdaWaqaaiabigdaXaqaaiabl6Uinbqaaiabl6UinbqaaiabicdaWaqaaiabigdaXaaaaiaawUfacaGLDbaacqqGGaaicqqGHbqycqqGUbGBcqqGKbazcqqGGaaicqWGdbWqcqGH9aqpdaWadaqaauaabeqabiaaaeaacqGHsislcqaIXaqmaeaacqaIXaqmaaaacaGLBbGaayzxaaGaeiOla4caaa@4D4F@

Under the null hypothesis, for each gene g and array i, an unbiased estimator of the expected value of the measurement X ig is obtained by the gene-wise mean value over all arrays from both groups. The transformation then becomes a subtraction of that mean value, cf. (7),

Y i g = X i g 1 n j = 1 n X j g . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGzbqwdaWgaaWcbaGaemyAaKMaem4zaCgabeaakiabg2da9iabdIfaynaaBaaaleaacqWGPbqAcqWGNbWzaeqaaOGaeyOeI0YaaSaaaeaacqaIXaqmaeaacqWGUbGBaaWaaabCaeaacqWGybawdaWgaaWcbaGaemOAaOMaem4zaCgabeaaaeaacqWGQbGAcqGH9aqpcqaIXaqmaeaacqWGUbGBa0GaeyyeIuoakiabc6caUaaa@453E@
(8)

Note how the transformation preserves the difference in mean value between the two groups of arrays.

If the elements in X g from the different arrays had in fact independent and identically distributed noise for each fixed gene g as assumed in OLM and unweighted LIMMA, the noise in Y g would have equal variances for all arrays. In Figure 1 array-wise density estimates for the transformed expression values are shown. For arrays from the same condition the distributions should be identical, reflecting the combined variability of signal and noise. For unregulated genes the expectation of Y g is zero, so if the assumption of few regulated genes holds the densities from all arrays should furthermore be essentially equal. Examination of Figure 1 reveals that neither of these statements are true, indicating that some variances are highly unequal.

Figure 1
figure 1

Density plots. Distribution of transformed expression values, Y, for the different arrays, in the two datasets. Colour-coding according to sample variance is used for increased clarity (blue for low variance, red for high variance). Differences in variability can be noted for both datasets.

Analogously, all pairs of arrays within each condition should have a common joint distribution and when few genes are regulated all pairs of arrays should essentially have a common joint distribution with a small negative correlation of -1/(n - 1). Examination of scatter plots for all pairs of arrays shows that this is clearly not the case (some obvious examples are shown in Figure 2, all pairs are included in Additional file 1 and Additional file 2).

Figure 2
figure 2

Pairwise plots. Transformed expression values, Y g , for selected pairs of arrays within the same group. Different pairs within the same group have distinctly different correlations. Upper triangle contains scatterplots. Lower triangle contains heatmaps of the corresponding two-dimensional kernel density estimates, where the majority of the genes are in the red portion of the plot, revealing important trends inside the black clouds. Diagonal red clouds in the heat maps reveal correlations between arrays. Off-diagonal numbers show estimated correlations from WAME. Diagonal boxes contain sample names and weights as well as estimated variances from WAME.

As expected from the observations above, unequal variances and non-zero correlations are estimated in the analyses with WAME, giving rise to highly unequal weights in the estimates of the differential expressions (shown in Table 1 and Table 2). In fact, the sign of the weight for some arrays even get switched compared to the sign of the weight of the other arrays from the same condition. This is an effect of strong correlations combined with unequal variances. It is an issue which is further addressed in Discussion.

Table 1 WAME weights for the Atrium dataset. Weights in percent from estimate of differential expression using WAME on the Atrium dataset.
Table 2 WAME weights for the COPD dataset. Weights in percent from estimate of differential expression using WAME on the COPD dataset.

The analysis methods described in Background are applied to the data and p-values and ranks computed. The respective probability plots are shown in Figure 3, demonstrating that there are substantial differences in the distribution of p-values between the different statistics. Since correlations and unequal variances are observed, the model assumptions of the alternative standard methods do not seem to hold. The p-values could thereby have become optimistic. On the other hand, it cannot be ruled out that the temporary assumption in WAME of no regulated genes makes its p-values conservative, which could also partly explain the differences. These problems are studied below by use of resampled data.

Figure 3
figure 3

Observed probability plots. Empirical distribution of p-values compared to the distribution expected for non-differentially expressed genes. The OLM and LIMMA curves largely coincide, as does the identity line and the WAME curve.

A common alternative to using the p-values as measures of significance is to consider the ranking of the genes, induced by the p-values or test statistics, and to select a fixed number of top ranked genes for further investigations. In Table 3 and Table 4 the concordance of the ranked lists are shown. The results from the included methods differ, for instance those from WAME compared to the other methods. This is not surprising since high correlations and highly unequal variances were identified by WAME, giving rise to highly unequal weights.

Table 3 Concordance of top lists in the Atrium dataset. Number of mutually included genes in the top-100 lists in the Atrium dataset as determined by the different methods.
Table 4 Concordance of top lists in the COPD dataset. Number of mutually included genes in the top-100 lists in the COPD dataset as determined by the different methods.

Resampled data

To examine closer the effect of violated assumptions of independence and identical distribution, we repeatedly selected two random subgroups of four arrays from within one group in the original data and performed tests between those groups. This was performed 100 times for the largest group in each of the two real datasets. Differentially expressed genes have unequal expected values in the two populations being sampled (cf. (2)). Since we now sample twice from the same condition, no differentially expressed genes exist.

Figure 4 shows the empirical p-value distributions for the resampled COPD data analysed with the four methods, together with the respective average empirical distribution,

Figure 4
figure 4

Probability plots. Empirical distributions of p-values for LIMMA, weighted LIMMA, OLM and WAME from tests on 100 resamples from the COPD dataset. Average empirical distribution indicated. Since no signal is added, the curves should ideally follow the diagonal.

F ¯ ( p ) = 1 100 i = 1 100 F i ( p ) , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacuWGgbGrgaqeaiabcIcaOiabdchaWjabcMcaPiabg2da9maalaaabaGaeGymaedabaGaeGymaeJaeGimaaJaeGimaadaamaaqahabaGaemOray0aaSbaaSqaaiabdMgaPbqabaGccqGGOaakcqWGWbaCcqGGPaqkaSqaaiabdMgaPjabg2da9iabigdaXaqaaiabigdaXiabicdaWiabicdaWaqdcqGHris5aOGaeiilaWcaaa@44D0@

where F i denotes the empirical CDF from the i th of the 100 resamples. For WAME, the p-value distributions are very close to the expected uniform. For OLM, LIMMA and weighted LIMMA there is a high variability between the p-value distributions and they are in many cases substantially different from the expected uniform. For WAME, OLM and LIMMA, the respective average empirical distribution is approximately correct, while for weighted LIMMA it is clearly optimistic. The results for the Atrium dataset (see Additional file 3) are very similar.

Evaluation of power

To evaluate the power of the tests in the studied datasets, a known regulation is added to randomly selected genes in one of the resampled groups, created according to the previous section. Thus, the noise is obtained from the real data and only the signal is synthetic. Ideally, the power can then be estimated by the proportion of differentially expressed genes that have a computed p-value less than a fixed level. However, valid p-values of the test statistics cannot be obtained from the respective models since, as demonstrated above, the corresponding assumptions are typically not valid. Ideally, the p-values would be determined by the true null distribution of the respective test statistics, given the array-wise quality deviations. In the simulation study, the critical value of the test statistics are therefore estimated from the empirical distribution of the test statistic for the unregulated genes. This is used to estimate the power of the different statistics (details are given in Methods).

The power estimates for the different methods are shown in Figure 5, for a level 0.1% test. The 0.1% level yields approximately 22 false positives if relatively few genes are in fact differentially expressed. For WAME, Σ is estimated both before and after adding a signal to 2228 genes (10%), thereby substantially affecting the estimate of Σ (cf. Figure 6). The powers of the two versions are nevertheless very similar (difference less than 0.003) and only the latter version is included in the plot.

Figure 5
figure 5

Estimated power. Estimated power in the simulated data for level 0.1% tests, based on resamples from the respective larger group in the Atrium and COPD datasets. Power is estimated at the marked points and spline interpolation is used in between.

Figure 6
figure 6

Average empirical p-value distribution for WAME under regulation. Average empirical p-value distribution of the unregulated genes, calculated using WAME, when 0%, 0.1%, 1%, 5% and 10% of the genes have a log2 differential expression of 1, i.e. a two-fold change. When genes are regulated the estimate of Σ is biased, leading to conservative, non-diagonal curves.

When the covariance-structure matrix Σ is estimated in WAME it is assumed that no genes are differentially expressed. Figure 6 includes the average empirical distribution for the p-values of the unregulated genes when different proportions of the genes have a log2 differential expression of 1. It is clear that the distributions are biased for high proportions, giving conservative p-values, which should be an effect of biased estimates of Σ.

The results from the studied datasets indicate (i) that WAME offers a relevant power increase compared to the included alternatives, (ii) that weighted LIMMA does not offer an advantage compared to LIMMA and (iii) that the moderated statistics (WAME, LIMMA and wLIMMA) are superior to the traditional methods of ranking by ordinary t-statistic (OLM) or estimated differential expression (FC).

Discussion

The WAME model and the simulations

The WAME model aims at catching quality deviations by one covariance-structure matrix common for all genes. This is certainly simplistic in some cases, e.g. when only certain physical parts of an array or certain types of mRNAs are of decreased quality. The estimated covariance structure can then only be expected to reflect a mixture of the qualities of the different genes. However, examining the simulations (Figure 5), we see a clear power gain in the WAME model compared to the other models. Also, WAME succeeds in catching enough of the quality deviations to make the p-value distributions more correct, thus providing increased usefulness of the p-values (Figure 3).

The models of LIMMA, weighted LIMMA and WAME are nested, where weighted LIMMA adds unequal variances and WAME adds unequal variances and correlations. Examination of Figure 1 shows that there are evident differences in variability between arrays. It is therefore interesting that we have not found a power increase of weighted LIMMA compared to LIMMA. Further, the p-values of weighted LIMMA turned out to be too optimistic (Figure 4). Comparison with the results of the WAME method, where the power increases and the p-value distributions get substantially more correct, suggests that the correlations are crucial in the model.

In the simulations, noise is taken from real data through resampling within a fixed group. This procedure provides data with fewer assumption on the noise structure compared to a fully parameterised simulation and should hopefully better reflect realistic situations. To evaluate the power of the different methods, a synthetic signal which is constant within each condition is added to the resample-based noise. This follows the assumption in the models of both WAME, OLM, LIMMA and weighted LIMMA, that the noise structure is equal for genes that are differentially expressed and non-differentially expressed. However, the biological variability of the expression of differentially expressed genes might be different under the different conditions due to the changed rôle of those genes. For complicated conditions such as complex diseases, the problem is more severe (cf. [1517]) since crucial genes might only be differentially expressed in a subset of the studied arrays. Further work is needed to evaluate the performance of WAME in such settings, as well as to possibly expand it to better fit these situations.

A relevant question regarding the modelling of quality deviations by the covariance-structure matrix Σ is whether biologically interesting features may be hidden by this model. In the present datasets, the question can partly be answered by examining the pairwise plots (cf. Figure 2) and noticing that a large proportion of the genes show similar deviations, which should speak against a specific interesting biological explanation. Also, the estimated covariance structure matrix Σ can be inspected with the goal of finding relevant correlations between arrays and thus highlighting interesting features in the data. Possible future work is to use such inspections to reveal unwanted features in normalisation or in preprocessing wet-lab steps that give rise to correlated errors for a large proportion of the genes.

Weights with switched signs

In the studied datasets, strong correlations combined with unequal variances make some weights within a group switch sign, in essence meaning that it is beneficial to partly subtract some arrays within a group in the estimate to be able to add more of the others in the same group (cf. Table 1 and Table 2). Since this might seem counter-intuitive, an elucidating example of possible mechanisms behind such weights follows.

Consider an example where two two-colour arrays are observed, X1 and X2. Let the two arrays have two sources of variation, one that is mutually independent (ε1, ε2) and one consisting of different proportions, a1 and a2, of one common source of variation η. Let ε1, ε2 and η be independent and normally distributed with expectation 0 and variances σ ε 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaaa@3137@ and σ η 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiGacqWFdpWCdaqhaaWcbaGae83TdGgabaGaeGOmaidaaaaa@313C@ , respectively. Furthermore, let μ be the parameter to be estimated. The model becomesX i = μ + a i η + ε i , i {1,2}.

Then, X1 gets a negative weight if and only if

a 1 > a 2 + σ ε 2 a 2 σ η 2 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGHbqydaWgaaWcbaGaeGymaedabeaakiabg6da+iabdggaHnaaBaaaleaacqaIYaGmaeqaaOGaey4kaSYaaSaaaeaaiiGacqWFdpWCdaqhaaWcbaGae8xTdugabaGaeGOmaidaaaGcbaGaemyyae2aaSbaaSqaaiabikdaYaqabaGccqWFdpWCdaqhaaWcbaGae83TdGgabaGaeGOmaidaaaaakiabcYcaSaaa@4000@

i.e. if array 1 includes a large enough contribution from the common source of variation. When a negative weight is allowed instead of removing the array, a smaller proportion of the common source of variation is included in the final estimate. Its precision is thus increased.

Validity of the p-values and derived entities

Varying quality of arrays and correlated errors were demonstrated in [2, 3] and in the present paper through examination of the data. These questions are typically neglected in microarray analyses, both when using parametric and when using non-parametric analysis methods, since independence and identical distribution or exchangeability are generally assumed under the null hypothesis. Thus, the validity is questionable of the corresponding p-values and their derived entities, e.g. false discovery rates and estimates of proportions of differentially expressed genes. This problem is obvious in the resample based simulations.

A number of experiments have been analysed (data not shown) in addition to those published in the present paper and in [2, 3]. In almost all cases relevant unequal variances and correlations have been identified, indicating that the problem is common.

In the resample based simulations with added signal, WAME is shown to be conservative, which is an effect of the biased estimate of Σ. Further work on an estimator of Σ with better characteristics under regulation is therefore needed. However, the simulations indicate (i) that the power of the test is basically unaffected by the bias and (ii) that hundreds of genes may be differentially expressed (two-fold) with only mildly conservative p-values as result.

Correlations between genes or between arrays?

It has recently been argued that the expression of different genes are highly dependent, making the law of large number normally inapplicable [18] and standard estimators of e.g. the false discovery rate (FDR) imprecise [19]. In [19], a latent FDR is introduced, being the conditional FDR given a random effect b that captures the correlation effects between genes. The FDR is then the marginal latent FDR, that is the average over the random effect b.

For the datasets examined in the present paper, the model assumptions of e.g. the ordinary linear model are shown not to hold (cf. Figure 1 and Figure 2). This can be expected to result in invalid p-values, which is indeed observed in Figure 4. Interestingly, the p-value distribution seem to be valid marginally, i.e. on average over the different resamples, which would yield valid but imprecise estimates of the FDR. This type of failed model assumptions is not taken into account in e.g. [18, 19]. Since for a performed experiment, the p-values from the ordinary t-statistic (OLM) share a common bias conditional on the experiment (see Figure 4), the different p-values may be highly dependent. However, this dependency is due to failure of taking array-wide quality deviations into account in the model and not due to the nature of microarray data per se, e.g. through substantial long-range gene-gene interactions.

Consequently, the strong observed dependencies between statistics from different genes might largely be explainable by quality deviations between the arrays in the experiment, e.g. correlations between arrays. Since WAME models these deviations such that the p-values are essentially correctly distributed when few genes are differentially expressed in the studied datasets, the dependency between genes should be greatly decreased. The covariance structure matrix Σ is therefore in a sense a parallel to the random factor b in [19]. It remains as future work to evaluate the gene-gene dependencies and estimates of e.g. the FDR in the context of the WAME model.

In the WAME model, the data from different genes are assumed independent, which is unrealistic, e.g. since genes act together in pathways. However, this is only used in the derivation of the maximum likelihood estimaties of the covariance structure matrix Σ and the shape parameter α. The assumption could thus be relaxed to a dependence between the different genes that is weak enough that the estimates of Σ and α become precise, and accurate under H0. This holds if the law of large numbers is applicable for averages of certain functions of the gene-wise observed data (cf. the likelihood functions in [2, 3]). Given the large number of genes and the observed p-value distributions in Figure 4, this relaxed assumption seems plausible.

It can be noted that for the studied data, WAME has higher power and considerably more valid p-values than weighted LIMMA. Since the difference between the weighted LIMMA and WAME models is the inclusion of correlations between arrays, this emphasises the importance of the correlations in the model.

Conclusion

Statistical methods in microarray analysis are typically based on the often erroneous assumption that the data from different arrays are independent and identically distributed. An exception is Weighted Analysis of Microarray Experiment (WAME) where heteroscedasticity and correlations between arrays are modelled by a covariance-structure common for all genes. In the present paper, WAME has been extended to handle datasets without a natural pairing, e.g. data from one-channel microarrays, and corresponding estimates and test statistics have been derived. In the examined one-channel microarray datasets WAME detected unequal variances and nonzero correlations.

WAME was compared with four other common methods: an ordinary linear model with t-tests, LIMMA, weighted LIMMA, and fold-change ranking. The comparison was performed using resampling of the different arrays within the datasets. Here, WAME had the highest power. When a relatively small proportion of the genes are regulated, WAME also generates close to correct p-value distributions while the p-value distributions from the other methods are highly variable. However, when many genes are differentially expressed, the p-values from WAME tend to be conservative.

In conclusion, p-values from the standard methods for microarray analysis should in general not be trusted and any result based on p-values, e.g. estimates of the number of regulated genes and false discovery rates, should be interpreted with care. The analyses of the examined datasets showed that WAME gives a powerful approach for finding differentially expressed genes and that it produces more trustworthy p-values when a relatively small proportion of genes are differentially expressed.

Methods

Details on the new version of WAME

Model Framework

For g = 1,..., m, let X g be an n-dimensional vector with expectation μ g = D γ g , where D is the design matrix, having rank k, and γ g qis the parameter vector. Furthermore, letX g | c g ~ N(μ g , c g Σ),c g ~ Γ-1(α, 1),

where Σ is the non-singular covariance-structure matrix, c g is the variance-scaling factor, α is the shape parameter for c g and (c1, X1),...,(c m , X m ) are assumed independent. The differential expression vector is defined asδ g = C γ g ,

where C is a matrix of rank p such that δ g is estimable. Here, an estimator of δ g and a test for

H 0 : δ g = 0 H A : δ g 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqadeGadaaabaGaemisaG0aaSbaaSqaaiabicdaWaqabaGccqGG6aGoiiWacqWF0oazdaWgaaWcbaGaem4zaCgabeaaaOqaaiabg2da9aqaaGqabiab+bdaWaqaaiabdIeainaaBaaaleaacqWGbbqqaeqaaOGaeiOoaOJae8hTdq2aaSbaaSqaaiabdEgaNbqabaaakeaacqGHGjsUaeaacqGFWaamaaaaaa@3E56@
(9)

are in focus.

As mentioned in Background, one main obstacle is that Σ is hard to estimate. In fact, Σ and δ g cannot be maximum likelihood estimated simultaneously, since there are trivial infinite suprema of the likelihood, e.g. when the variance of one observation is set to zero and the corresponding mean is selected so that it equals that observation.

The current WAME method

In the current version of WAME [3], Σ is estimated as follows. First, temporarily assume that μ g = 0 for all genes, which is reasonable for paired experimental designs with few differentially expressed genes between any pairwise measured conditions. For each gene, the variance scaling factor c g is removed by dividing the n measurements with the first measurement, yielding a random vector distributed according to a multivariate generalisation of the Cauchy distribution. A scaled version of Σ is then maximum likelihood estimated numerically. Second, the unknown scale and the hyperparameter α of the prior distribution of c g are maximum likelihood estimated numerically without the assumption of μ g = 0. The parameters Σ and α are subsequently treated as known in the maximum-likelihood estimates and likelihood-ratio tests for the different genes.

The new WAME method

The new version of WAME relaxes the assumption from μ g = 0 to δ g = 0, which incorporates many designs without a natural pairing. This is performed by subtracting an arbitrary estimator μ ˜ g 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaDaaaleaacqWGNbWzaeaacqaIWaamaaaaaa@30EB@ of μ g , which is unbiased under H0 and has as image the space V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 of possible values for μ g under H0,

Y g = X g μ ˜ g 0 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaieqacqWFzbqwdaWgaaWcbaGaem4zaCgabeaakiabg2da9iab=HfaynaaBaaaleaacqWGNbWzaeqaaOGaeyOeI0cccmGaf4hVd0MbaGaadaqhaaWcbaGaem4zaCgabaGaeGimaadaaOGaeiOla4caaa@395B@
(10)

It can be shown that this transformation preserves the estimability of δ g .

By construction, the transformed data Y g will have expectation zero for non-differentially expressed genes and the current WAME method can be applied on Y g , including the estimation of the covariance-structure matrix Σ Y for Y g . It will now be proved that the likelihood ratio tests of (9) and the maximum likelihood estimates of δ g based on X g or Y g are identical, if α and Σ or Σ Y respectively are considered known.

We shall henceforth consider a fixed gene g and drop the g index.

Equality of tests and estimators

Before beginning, some further definitions are needed. Define the Mahalanobis inner product corresponding to a symmetric n by n matrix A as

x 1 , x 2 A = x 1 T A x 2 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqGHPms4ieqacqWF4baEdaWgaaWcbaGaeGymaedabeaakiabcYcaSiab=Hha4naaBaaaleaacqaIYaGmaeqaaOGaeyOkJe=aaSbaaSqaaiabdgeabbqabaGccqGH9aqpcqWF4baEdaqhaaWcbaGaeGymaedabaGaeeivaqfaaOGaemyqae0aaWbaaSqabeaacqGHsislaaGccqWF4baEdaWgaaWcbaGaeGOmaidabeaakiabcYcaSaaa@420F@
(11)

and the norm · A as

x A 2 = x , x = x T A x , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaadaqbdaqaaGqabiab=Hha4bGaayzcSlaawQa7amaaDaaaleaacqWGbbqqaeaacqaIYaGmaaGccqGH9aqpcqGHPms4cqWF4baEcqGGSaalcqWF4baEcqGHQms8cqGH9aqpcqWF4baEdaahaaWcbeqaaiabbsfaubaakiabdgeabnaaCaaaleqabaGaeyOeI0caaOGae8hEaGNaeiilaWcaaa@443E@

where x, x1, x2 lies in the rowspace of A and the generalised inverse A-is any matrix satisfying AA-A = A. Let χ denote the n-dimensional inner product space with ·,·Σ as inner product. Define V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ χ as the space of possible values for μ g ,

V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@

= {μ: μ= D γ, γ q}

and let V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 χ denote the corresponding space restricted by the null hypothesis,

V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@

0 = {μ: μ= D γ, C γ = 0, γ q}.

Proposition Let μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@ be an arbitrary linear estimator of μ, which is unbiased under H0 and which has image V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0. Let

Y = X μ ˜ 0 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaieqacqWFzbqwcqGH9aqpcqWFybawcqGHsisliiWacuGF8oqBgaacamaaCaaaleqabaGaeGimaadaaOGaeiilaWcaaa@34E6@

and let Σ Y be the covariance-structure matrix of Y. Then the likelihood ratio test of (9) and the maximum likelihood estimate of δ based on X with Σ and α known are identical to the ones based on Y with Σ Y and α known.

Proof of the Proposition

The proof is divided into two steps which combined conclude the proof.

  1. 1.

    The likelihood ratio test (LRT) of (9) and the maximum likelihood estimator (MLE) of δ does not depend on the choice of μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@ .

  2. 2.

    The proposition holds when μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@ is the MLE of μ under H0.

Proof of step 1

Let μ' and μ'' be two valid choices of μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@ , i.e. they are both unbiased estimators of μ under H0 and have V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 as image. Let Y' = X - μ' and Y'' = X - μ''. Recall that a matrix P is a projection matrix projecting on V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 if and only if for all x n, Px V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 and for all x0 V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0, P x0 = x0. It can be shown that μ' and μ'' can be written as μ'= P' X and μ''= P'' X for some projection matrices P' and P'' projecting on V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0. Since P' and P'' project on the same space it follows that P' P'' = P'' and P'' P' = P', and thus (I - P') Y'' = Y' and (I - P'') Y' = Y''. Hence there is an invertible map between Y' and Y'' and likelihood methods based on Y' and Y'' respectively will give equal results. Consequently, the MLE of (9) and the LRT of δ will not depend on the choice of μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@

Proof of step 2

Since δ is estimable based on X, there exist a matrix A such that C = AD and thus δ= A μ. The likelihood of μ can therefore be examined instead of the likelihood of δ.

The likelihood of μ based on X can be shown to be

L ( μ | X ) = 0 f ( X | μ , c ) f ( c ) d c [ X μ Σ 2 / 2 + 1 ] n / 2 α , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGadaaabaGaemitaWKaeiikaGcccmGae8hVd0MaeiiFaWhcbeGae4hwaGLaeiykaKcabaGaeyypa0dabaWaa8qmaeaacqWGMbGzcqGGOaakcqGFybawcqGG8baFcqWF8oqBcqGGSaalcqWGJbWycqGGPaqkcqGHflY1cqWGMbGzcqGGOaakcqWGJbWycqGGPaqkcqWGKbazcqWGJbWyaSqaaiabicdaWaqaaiabg6HiLcqdcqGHRiI8aaGcbaaabaGaeyyhIulabaWaamWaaeaadaqbdaqaaiab+HfayjabgkHiTiab=X7aTbGaayzcSlaawQa7amaaDaaaleaacqqHJoWuaeaacqaIYaGmaaGccqGGVaWlcqaIYaGmcqGHRaWkcqaIXaqmaiaawUfacaGLDbaadaahaaWcbeqaaiabgkHiTiabd6gaUjabc+caViabikdaYiabgkHiTGGaciab9f7aHbaakiabcYcaSaaaaaa@6514@
(12)

where denotes proportionality. Using the Projection Theorem [20], the MLE of μ is the orthogonal projection of X on V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ ,

μ ^ = P V X , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaqcaiabg2da9iabbcfaqnaaBaaaleaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqGFveVvaeqaaGqabOGae0hwaGLaeiilaWcaaa@3E97@

where the orthogonality is according to the inner product of χ. When H0 is true, μ is restricted to V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 and thus the MLE of μ becomes

μ ^ 0 = P V 0 X . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaqcamaaCaaaleqabaGaeGimaadaaOGaeyypa0Jaeeiuaa1aaSbaaSqaamrtHrhAL1wy0L2yHvtyaeHbnfgDOvwBHrxAJfwnaGabaiab+vr8wnaaBaaameaacqaIWaamaeqaaaWcbeaaieqakiab9HfayHqaaiab85caUaaa@40E9@

Note that μ ^ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaqcamaaCaaaleqabaGaeGimaadaaaaa@2F95@ is a valid choice for μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@ , i.e. μ ^ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaqcamaaCaaaleqabaGaeGimaadaaaaa@2F95@ is unbiased under H0 and has V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 as image. Let

Z = X μ ^ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaieqacqWFAbGwcqGH9aqpcqWFybawcqGHsisliiWacuGF8oqBgaqcamaaCaaaleqabaGaeGimaadaaaaa@33FF@

which gives Z = P V 0 X MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaieqacqWFAbGwcqGH9aqpcqqGqbaudaWgaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae4xfXB1aa0baaWqaaiabicdaWaqaambvLH1qnnvETj2BSbqehmwBZLxmWaacfaGae0xPI8daaaWcbeaakiab=Hfaybaa@4640@ , where V 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvdaqhaaWcbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaaa@415B@ denotes the orthogonal complement of V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 in χ. Standard properties of the normal distribution givesM

Z|c ~ N(μ z , c Σ z ),

where μ z = D z γ with D z = P V 0 D MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqGqbaudaWgaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae8xfXB1aa0baaWqaaiabicdaWaqaambvLH1qnnvETj2BSbqehmwBZLxmWaacfaGae4xPI8daaaWcbeaakiabdseaebaa@43D5@ , and where Σ z = P V 0 Σ P V 0 T MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqHJoWudaWgaaWcbaGaemOEaOhabeaakiabg2da9iabbcfaqnaaBaaaleaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvdaqhaaadbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaleqaaOGaeu4OdmLaeeiuaa1aa0baaSqaaiab=vr8wnaaDaaameaacqaIWaamaeaacqGFLkYpaaaaleaacqqGubavaaaaaa@4F9A@ .

The likelihood function of μ z (with respect to the Lebesgue measure on the space of possible values of Z spanned by the column vectors of Σ z ) can be written as

L ( μ | Z ) [ Z μ z Σ z 2 / 2 + 1 ] n / 2 α . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGmbatcqGGOaakiiWacqWF8oqBcqGG8baFieqacqGFAbGwcqGGPaqkcqGHDisTdaWadaqaamaafmaabaGae4NwaOLaeyOeI0Iae8hVd02aaSbaaSqaaiabdQha6bqabaaakiaawMa7caGLkWoadaqhaaWcbaGaeu4Odm1aaSbaaWqaaiabdQha6bqabaaaleaacqaIYaGmaaGccqGGVaWlcqaIYaGmcqGHRaWkcqaIXaqmaiaawUfacaGLDbaadaahaaWcbeqaaiabgkHiTiabd6gaUjabc+caViabikdaYiabgkHiTGGaciab9f7aHbaakiabc6caUaaa@5003@

Since, δ is estimable based on Z, the likelihood of μ z can be examined instead of the likelihood of δ.

The likelihood ratio statistic of (9) based on X is defined by

T = sup μ V L ( μ | X ) sup μ V 0 L ( μ | X ) , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGubavcqGH9aqpdaWcaaqaamaaxababaGagi4CamNaeiyDauNaeiiCaahaleaaiiWacqWF8oqBcqGHiiIZt0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqGFveVvaeqaaOGaemitaWKaeiikaGIae8hVd0MaeiiFaWhcbeGae0hwaGLaeiykaKcabaWaaCbeaeaacyGGZbWCcqGG1bqDcqGGWbaCaSqaaiab=X7aTjabgIGiolab+vr8wnaaBaaameaacqaIWaamaeqaaaWcbeaakiabdYeamjabcIcaOiab=X7aTjabcYha8jab9HfayjabcMcaPaaacqGGSaalaaa@5C7F@

which can be rewritten (cf. [3]) as a strictly increasing function of

T = n p + 2 α k P V X P V 0 X Σ 2 X P V X Σ 2 + 2 = n p + 2 α k P V V 0 X Σ 2 P V X Σ 2 + 2 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqadeGabaaabaGafmivaqLbauaacqGH9aqpdaWcaaqaaiabd6gaUjabgkHiTiabdchaWjabgUcaRiabikdaYGGaciab=f7aHbqaaiabdUgaRbaadaWcaaqaamaafmaabaGaeeiuaa1aaSbaaSqaamrtHrhAL1wy0L2yHvtyaeHbnfgDOvwBHrxAJfwnaGabaiab+vr8wbqabaacbeGccqqFybawcqGHsislcqqGqbaudaWgaaWcbaGae4xfXB1aaSbaaWqaaiabicdaWaqabaaaleqaaOGae0hwaGfacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaaaOqaamaafmaabaGae0hwaGLaeyOeI0Iaeeiuaa1aaSbaaSqaaiab+vr8wbqabaGccqqFybawaiaawMa7caGLkWoadaqhaaWcbaGaeu4OdmfabaGaeGOmaidaaOGaey4kaSIaeGOmaidaaaqaaiabg2da9maalaaabaGaemOBa4MaeyOeI0IaemiCaaNaey4kaSIaeGOmaiJae8xSdegabaGaem4AaSgaamaalaaabaWaauWaaeaacqqGqbaudaWgaaWcbaGae4xfXBLaeyykICSae4xfXB1aa0baaWqaaiabbcdaWaqaambvLH1qnnvETj2BSbqehmwBZLxmWaacfaGaeWxPI8daaaWcbeaakiab9HfaybGaayzcSlaawQa7amaaDaaaleaacqqHJoWuaeaacqaIYaGmaaaakeaadaqbdaqaaiabbcfaqnaaBaaaleaacqGFveVvdaahaaadbeqaaiab8vQi=aaaaSqabaGccqqFybawaiaawMa7caGLkWoadaqhaaWcbaGaeu4OdmfabaGaeGOmaidaaOGaey4kaSIaeGOmaidaaiabcYcaSaaaaaa@8FB5@
(13)

where V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ and V 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvdaqhaaWcbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaaa@415B@ are the orthogonal complements of V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ and V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0 respectively.

Note that the space of possible values for μ z is V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ V 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvdaqhaaWcbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaaa@415B@ and that μ z = 0 under H0. Let P z MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqGqbaudaahaaWcbeqaaiabdQha6baaaaa@2F7D@ denote the orthogonal projection according to , Σ z MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqGHPms4cqGHflY1cqGGSaalcqGHflY1cqGHQms8daWgaaWcbaGaeu4Odm1aaSbaaWqaaiabdQha6bqabaaaleqaaaaa@3908@ . Then, the likelihood ratio statistic of (9) based on Z can in analogy with (13) be shown to be a strictly increasing function of

T z = n p + 2 α k P V V 0 z Z Σ z 2 Z P V V 0 z Z Σ z 2 + 2 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacuWGubavgaqbamaaBaaaleaacqWG6bGEaeqaaOGaeyypa0ZaaSaaaeaacqWGUbGBcqGHsislcqWGWbaCcqGHRaWkcqaIYaGmiiGacqWFXoqyaeaacqWGRbWAaaWaaSaaaeaadaqbdaqaaiabbcfaqnaaDaaaleaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqGFveVvcqGHPiYXcqGFveVvdaqhaaadbaGaeeimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqqFLkYpaaaaleaaieGacqaF6bGEaaacbeGccqWEAbGwaiaawMa7caGLkWoadaqhaaWcbaGaeu4Odm1aaSbaaWqaaiabdQha6bqabaaaleaacqaIYaGmaaaakeaadaqbdaqaaiab7PfaAjab=jHiTiabbcfaqnaaDaaaleaacqGFveVvcqGHPiYXcqGFveVvdaqhaaadbaGaeeimaadabaGae0xPI8daaaWcbaGaeWNEaOhaaOGaeSNwaOfacaGLjWUaayPcSdWaa0baaSqaaiabfo6atnaaBaaameaacqWG6bGEaeqaaaWcbaGaeGOmaidaaOGaey4kaSIaeGOmaidaaiabc6caUaaa@7599@
(14)

The Lemma below yields that for all W V 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFwe=vcqGHgksZcqWFveVvdaqhaaWcbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaaa@453F@ and all z V 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvdaqhaaWcbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaaa@415B@ , z Σ z 2 = z Σ 2 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaadaqbdaqaaGqabiab=Pha6bGaayzcSlaawQa7amaaDaaaleaacqqHJoWudaWgaaadbaGaemOEaOhabeaaaSqaaiabikdaYaaakiabg2da9maafmaabaGae8NEaOhacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaaaaa@3E01@ and P W z z = P W z MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqGqbaudaqhaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae8NfXFfabaGaemOEaOhaaGqabOGae4NEaONaeyypa0Jaeeiuaa1aaSbaaSqaaiab=zr8xbqabaGccqGF6bGEaaa@4267@ . The equivalence of the test statistics (13) and (14) is now straight-forward,

T z = n p + 2 α k P V V 0 z Z Σ z 2 Z v V V 0 z Z Σ z 2 + 2 = n p + 2 α k P V V 0 P V 0 X Σ 2 ( P V + P V ) ( P V 0 X P V V 0 P V 0 X ) Σ 2 + 2 = n p + 2 α k P V V 0 X Σ 2 P V X Σ 2 + 2 = T . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeWadaaabaGafmivaqLbauaadaWgaaWcbaGaemOEaOhabeaaaOqaaiabg2da9aqaamaalaaabaGaemOBa4MaeyOeI0IaemiCaaNaey4kaSIaeGOmaidcciGae8xSdegabaGaem4AaSgaamaalaaabaWaauWaaeaacqqGqbaudaqhaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae4xfXBLaeyykICSae4xfXB1aa0baaWqaaiabbcdaWaqaambvLH1qnnvETj2BSbqehmwBZLxmWaacfaGae0xPI8daaaWcbaGaemOEaOhaaGqabOGaeWNwaOfacaGLjWUaayPcSdWaa0baaSqaaiabfo6atnaaBaaameaacqWG6bGEaeqaaaWcbaGaeGOmaidaaaGcbaWaauWaaeaacqaFAbGwcqGHsislcqGFVeVDdaqhaaWcbaGae4xfXBLaeyykICSae4xfXB1aa0baaWqaaiabbcdaWaqaaiab9vQi=aaaaSqaaiabdQha6baakiab8PfaAbGaayzcSlaawQa7amaaDaaaleaacqqHJoWudaWgaaadbaGaemOEaOhabeaaaSqaaiabikdaYaaakiabgUcaRiabikdaYaaaaeaaaeaacqGH9aqpaeaadaWcaaqaaiabd6gaUjabgkHiTiabdchaWjabgUcaRiabikdaYiab=f7aHbqaaiabdUgaRbaadaWcaaqaamaafmaabaGaeeiuaa1aaSbaaSqaaiab+vr8wjabgMIihlab+vr8wnaaDaaameaacqaIWaamaeaacqqFLkYpaaaaleqaaOGaeeiuaa1aaSbaaSqaaiab+vr8wnaaDaaameaacqaIWaamaeaacqqFLkYpaaaaleqaaOGaeWhwaGfacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaaaOqaamaafmaabaGaeiikaGIaeeiuaa1aaSbaaSqaaiab+vr8wbqabaGccqGHRaWkcqqGqbaudaWgaaWcbaGae4xfXB1aaWbaaWqabeaacqqFLkYpaaaaleqaaOGaeiykaKIaeiikaGIaeeiuaa1aaSbaaSqaaiab+vr8wnaaDaaameaacqaIWaamaeaacqqFLkYpaaaaleqaaOGaeWhwaGLaeyOeI0Iaeeiuaa1aaSbaaSqaaiab+vr8wjabgMIihlab+vr8wnaaDaaameaacqaIWaamaeaacqqFLkYpaaaaleqaaOGaeeiuaa1aaSbaaSqaaiab+vr8wnaaDaaameaacqaIWaamaeaacqqFLkYpaaaaleqaaOGaeWhwaGLaeiykaKcacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaakiabgUcaRiabikdaYaaaaeaaaeaacqGH9aqpaeaadaWcaaqaaiabd6gaUjabgkHiTiabdchaWjabgUcaRiabikdaYiab=f7aHbqaaiabdUgaRbaadaWcaaqaamaafmaabaGaeeiuaa1aaSbaaSqaaiab+vr8wjabgMIihlab+vr8wnaaDaaameaacqaIWaamaeaacqqFLkYpaaaaleqaaOGaeWhwaGfacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaaaOqaamaafmaabaGaeeiuaa1aaSbaaSqaaiab+vr8wnaaCaaameqabaGae0xPI8daaaWcbeaakiab8HfaybGaayzcSlaawQa7amaaDaaaleaacqqHJoWuaeaacqaIYaGmaaGccqGHRaWkcqaIYaGmaaGaeyypa0JafmivaqLbauaacqGGUaGlaaaaaa@EDB6@
(15)

Lemma Let W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFwe=vaaa@384C@ be a subspace of χ and let P W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqGqbaudaWgaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae8NfXFfabeaaaaa@399F@ be the orthogonal projection from χ onto W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFwe=vaaa@384C@ . Then for any x1, x2 W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFwe=vaaa@384C@ ,

x 1 , x 2 Σ = x 1 , x 2 Σ W , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqGHPms4ieqacqWF4baEdaWgaaWcbaGaeGymaedabeaakiabcYcaSiab=Hha4naaBaaaleaacqaIYaGmaeqaaOGaeyOkJe=aaSbaaSqaaiabfo6atbqabaGccqGH9aqpcqGHPms4cqWF4baEdaWgaaWcbaGaeGymaedabeaakiabcYcaSiab=Hha4naaBaaaleaacqaIYaGmaeqaaOGaeyOkJe=aaSbaaSqaaiabfo6atnaaBaaameaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqGFwe=vaeqaaaWcbeaakiabcYcaSaaa@511D@

where Σ W = P W Σ P W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqHJoWudaWgaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae8NfXFfabeaakiabg2da9iabbcfaqnaaBaaaleaacqWFwe=vaeqaaOGaeu4OdmLaeeiuaa1aaSbaaSqaaiab=zr8xbqabaaaaa@4306@ .

Proof Let A be a matrix of a change of basis [20] from the standard basis to an orthonormal basis of χ such that the first l basis vectors span W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFwe=vaaa@384C@ . Let I(l)denote the identity matrix with all but the l top left diagonal elements set to zero. It follows that AT A = Σ-1 and A P W MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqqGqbaudaWgaaWcbaWenfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae8NfXFfabeaaaaa@399F@ = I(l)A and therefore,

x 1 , x 2 Σ = x 1 T Σ 1 x 2 = x 1 T A T A P W x 2 = x 1 T A T ( I ( l ) ) A x 2 = x 1 T A T ( I ( l ) A Σ A T I ( l ) ) A x 2 = x 1 T A T ( A P W Σ P W T A T ) A x 2 = x 1 T ( P W Σ P W T ) x 2 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGbdaaaaeaacqGHPms4ieqacqWF4baEdaWgaaWcbaGaeGymaedabeaakiabcYcaSiab=Hha4naaBaaaleaacqaIYaGmaeqaaOGaeyOkJe=aaSbaaSqaaiabfo6atbqabaaakeaacqGH9aqpaeaacqWF4baEdaqhaaWcbaGaeGymaedabaGaeeivaqfaaOGaeu4Odm1aaWbaaSqabeaacqGHsislcqaIXaqmaaGccqWF4baEdaWgaaWcbaGaeGOmaidabeaaaOqaaaqaaiabg2da9aqaaiab=Hha4naaDaaaleaacqaIXaqmaeaacqqGubavaaGccqWGbbqqdaahaaWcbeqaaiabbsfaubaakiabdgeabjabbcfaqnaaBaaaleaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqGFwe=vaeqaaOGae8hEaG3aaSbaaSqaaiabikdaYaqabaaakeaaaeaacqGH9aqpaeaacqWF4baEdaqhaaWcbaGaeGymaedabaGaeeivaqfaaOGaemyqae0aaWbaaSqabeaacqqGubavaaGccqGGOaakcqWGjbqsdaWgaaWcbaGaeiikaGIaemiBaWMaeiykaKcabeaakiabcMcaPmaaCaaaleqabaGaeyOeI0caaOGaemyqaeKae8hEaG3aaSbaaSqaaiabikdaYaqabaaakeaaaeaacqGH9aqpaeaacqWF4baEdaqhaaWcbaGaeGymaedabaGaeeivaqfaaOGaemyqae0aaWbaaSqabeaacqqGubavaaGccqGGOaakcqWGjbqsdaWgaaWcbaGaeiikaGIaemiBaWMaeiykaKcabeaakiabdgeabjabfo6atjabdgeabnaaCaaaleqabaGaeeivaqfaaOGaemysaK0aaSbaaSqaaiabcIcaOiabdYgaSjabcMcaPaqabaGccqGGPaqkdaahaaWcbeqaaiabgkHiTaaakiabdgeabjab=Hha4naaBaaaleaacqaIYaGmaeqaaaGcbaaabaGaeyypa0dabaGae8hEaG3aa0baaSqaaiabigdaXaqaaiabbsfaubaakiabdgeabnaaCaaaleqabaGaeeivaqfaaOGaeiikaGIaemyqaeKaeeiuaa1aaSbaaSqaaiab+zr8xbqabaGccqqHJoWucqqGqbaudaqhaaWcbaGae4NfXFfabaGaeeivaqfaaOGaemyqae0aaWbaaSqabeaacqqGubavaaGccqGGPaqkdaahaaWcbeqaaiabgkHiTaaakiabdgeabjab=Hha4naaBaaaleaacqaIYaGmaeqaaaGcbaaabaGaeyypa0dabaGae8hEaG3aa0baaSqaaiabigdaXaqaaiabbsfaubaakiabcIcaOiabbcfaqnaaBaaaleaacqGFwe=vaeqaaOGaeu4OdmLaeeiuaa1aa0baaSqaaiab+zr8xbqaaiabbsfaubaakiabcMcaPmaaCaaaleqabaGaeyOeI0caaOGae8hEaG3aaSbaaSqaaiabikdaYaqabaGccqGGSaalaaaaaa@B789@

where the last equality uses the fact that (AB)- = B- A-1 when A is invertible □.

The next step is to show that the MLE of δ when X is observed is identical to the MLE of δ when Z is observed. The former is defined by

δ ^ = C γ ^ = C arg min γ X D γ Σ 2 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF0oazgaqcaiabg2da9iabdoeadjqb=n7aNzaajaGaeyypa0Jaem4qam0aaCbeaeaacyGGHbqycqGGYbGCcqGGNbWzcyGGTbqBcqGGPbqAcqGGUbGBaSqaaiab=n7aNbqabaGcdaqbdaqaaGqabiab+HfayjabgkHiTiabdseaejab=n7aNbGaayzcSlaawQa7amaaDaaaleaacqqHJoWuaeaacqaIYaGmaaGccqGGUaGlaaa@49F0@

Define G MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaaliab=zq8hbaa@3837@ 0 = {γ: D γ V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ 0} and G MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaaliab=zq8hbaa@3837@ 1 = {γ: D γ V 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvdaqhaaWcbaGaeGimaadabaWeuvgwd10u51MyVXgarCWyTnxEXadaiuaacqGFLkYpaaaaaa@415B@ } and note that for any γ there exist γ0 G MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaaliab=zq8hbaa@3837@ 0 and γ1 G MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaaliab=zq8hbaa@3837@ 1 such that γ= γ0 + γ1. Thus,

δ ^ = C arg min γ 0 + γ 1 : γ 0 G 0 , γ 1 G 1 X D ( γ 0 + γ 1 ) Σ 2 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF0oazgaqcaiabg2da9iabdoeadnaaxababaGagiyyaeMaeiOCaiNaei4zaCMagiyBa0MaeiyAaKMaeiOBa4galeaacqWFZoWzdaWgaaadbaGaeGimaadabeaaliabgUcaRiab=n7aNnaaBaaameaacqaIXaqmaeqaaSGaeiOoaOJae83SdC2aaSbaaWqaaiabicdaWaqabaWccqGHiiIZt0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqGFge=rdaWgaaadbaGaeGimaadabeaaliabcYcaSiab=n7aNnaaBaaameaacqaIXaqmaeqaaSGaeyicI4Sae4NbXF0aaSbaaWqaaiabigdaXaqabaaaleqaaOWaauWaaeaaieqacqqFybawcqGHsislcqWGebarcqGGOaakcqWFZoWzdaWgaaWcbaGaeGimaadabeaakiabgUcaRiab=n7aNnaaBaaaleaacqaIXaqmaeqaaOGaeiykaKcacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaakiabc6caUaaa@6B7F@

Now, since P V 0 + P V 0 = I MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaliabbcfaqPWaaSbaaSqaamrtHrhAL1wy0L2yHvtyaeHbnfgDOvwBHrxAJfwnaGabaiab=vr8wnaaDaaameaacqaIWaamaeaatqvzynuttLxBI9gBaeXbJ12C5fdmaGqbaiab+vQi=aaaaSqabaGccqGHRaWkliabbcfaqPWaaSbaaSqaaiab=vr8wnaaBaaameaacqaIWaamaeqaaaWcbeaakiabg2da9iabdMeajbaa@4A55@ ,

δ ^ = C arg min γ 0 + γ 1 : γ 0 G 0 , γ 1 G 1 ( P V 0 ( X D ( γ 0 + γ 1 ) ) + P V 0 ( X D ( γ 0 + γ 1 ) ) Σ 2 ) = C arg min γ 0 + γ 1 : γ 0 G 0 , γ 1 G 1 ( P V 0 ( X D γ 0 ) Σ 2 + P V 0 ( X D γ 1 ) Σ z 2 ) , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGadaaabaaccmGaf8hTdqMbaKaaaeaacqGH9aqpaeaacqWGdbWqdaWfqaqaaiGbcggaHjabckhaYjabcEgaNjGbc2gaTjabcMgaPjabc6gaUbWcbaGae83SdC2aaSbaaWqaaiabicdaWaqabaWccqGHRaWkcqWFZoWzdaWgaaadbaGaeGymaedabeaaliabcQda6iab=n7aNnaaBaaameaacqaIWaamaeqaaSGaeyicI48enfgDOvwBHrxAJfwnHbqeg0uy0HwzTfgDPnwy1aaceaGae4NbXF0aaSbaaWqaaiabicdaWaqabaWccqGGSaalcqWFZoWzdaWgaaadbaGaeGymaedabeaaliabgIGiolab+zq8hnaaBaaameaacqaIXaqmaeqaaaWcbeaakmaabmaabaWaauWaaeaaliabbcfaqPWaaSbaaSqaaiab+vr8wnaaBaaameaacqaIWaamaeqaaaWcbeaakiabcIcaOGqabiab9HfayjabgkHiTiabdseaejabcIcaOiab=n7aNnaaBaaaleaacqaIWaamaeqaaOGaey4kaSIae83SdC2aaSbaaSqaaiabigdaXaqabaGccqGGPaqkcqGGPaqkcqGHRaWkliabbcfaqPWaaSbaaSqaaiabbAfawnaaDaaameaacqaIWaamaeaatqvzynuttLxBI9gBaeXbJ12C5fdmaGqbaiab8vQi=aaaaSqabaGccqGGOaakcqqFybawcqGHsislcqWGebarcqGGOaakcqWFZoWzdaWgaaWcbaGaeGimaadabeaakiabgUcaRiab=n7aNnaaBaaaleaacqaIXaqmaeqaaOGaeiykaKIaeiykaKcacaGLjWUaayPcSdWaa0baaSqaaiabfo6atbqaaiabikdaYaaaaOGaayjkaiaawMcaaaqaaaqaaiabg2da9aqaaiabdoeadnaaxababaGagiyyaeMaeiOCaiNaei4zaCMagiyBa0MaeiyAaKMaeiOBa4galeaacqWFZoWzdaWgaaadbaGaeGimaadabeaaliabgUcaRiab=n7aNnaaBaaameaacqaIXaqmaeqaaSGaeiOoaOJae83SdC2aaSbaaWqaaiabicdaWaqabaWccqGHiiIZcqGFge=rdaWgaaadbaGaeGimaadabeaaliabcYcaSiab=n7aNnaaBaaameaacqaIXaqmaeqaaSGaeyicI4Sae4NbXF0aaSbaaWqaaiabigdaXaqabaaaleqaaOWaaeWaaeaadaqbdaqaaSGaeeiuaaLcdaWgaaWcbaGae4xfXB1aaSbaaWqaaiabicdaWaqabaaaleqaaOGaeiikaGIae0hwaGLaeyOeI0IaemiraqKae83SdC2aaSbaaSqaaiabicdaWaqabaGccqGGPaqkaiaawMa7caGLkWoadaqhaaWcbaGaeu4OdmfabaGaeGOmaidaaOGaey4kaSYaauWaaeaaliabbcfaqPWaaSbaaSqaaiab+vr8wnaaDaaameaacqaIWaamaeaacqaFLkYpaaaaleqaaOGaeiikaGIae0hwaGLaeyOeI0IaemiraqKae83SdC2aaSbaaSqaaiabigdaXaqabaGccqGGPaqkaiaawMa7caGLkWoadaqhaaWcbaGaeu4Odm1aaSbaaWqaaiabdQha6bqabaaaleaacqaIYaGmaaaakiaawIcacaGLPaaacqGGSaalaaaaaa@D7C0@

where the second equality follows from the generalised Theorem of Pythagoras [20], the Lemma, and the fact that P V 0 D γ 0 = 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaliabbcfaqPWaaSbaaSqaamrtHrhAL1wy0L2yHvtyaeHbnfgDOvwBHrxAJfwnaGabaiab=vr8wnaaDaaameaacqaIWaamaeaatqvzynuttLxBI9gBaeXbJ12C5fdmaGqbaiab+vQi=aaaaSqabaGccqWGebariiWacqqFZoWzdaWgaaWcbaGaeGimaadabeaakiabg2da9iabicdaWaaa@48AF@ and P V 0 D γ 1 = 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaliabbcfaqPWaaSbaaSqaamrtHrhAL1wy0L2yHvtyaeHbnfgDOvwBHrxAJfwnaGabaiab=vr8wnaaBaaameaacqaIWaamaeqaaaWcbeaakiabdseaeHGadiab+n7aNnaaBaaaleaacqaIXaqmaeqaaOGaeyypa0JaeGimaadaaa@40BB@ . Now since γ0 and γ1 minimise the expression independently of each other and since C γ0 = 0 by construction,

δ ^ = C ( arg min γ 0 G 0 P V 0 ( X D γ 0 ) Σ 2 + arg min γ 1 G 1 Z D z γ 1 Σ z 2 ) = C arg min γ 1 G 1 Z D z γ 1 Σ z 2 . MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaafaqaaeGadaaabaaccmGaf8hTdqMbaKaaaeaacqGH9aqpaeaacqWGdbWqdaqadaqaamaaxababaGagiyyaeMaeiOCaiNaei4zaCMagiyBa0MaeiyAaKMaeiOBa4galeaacqWFZoWzdaWgaaadbaGaeGimaadabeaaliabgIGioprtHrhAL1wy0L2yHvtyaeHbnfgDOvwBHrxAJfwnaGabaiab+zq8hnaaBaaameaacqaIWaamaeqaaaWcbeaakmaafmaabaWccqqGqbaukmaaBaaaleaacqGFveVvdaWgaaadbaGaeGimaadabeaaaSqabaGccqGGOaakieqacqqFybawcqGHsislcqWGebarcqWFZoWzdaWgaaWcbaGaeGimaadabeaakiabcMcaPaGaayzcSlaawQa7amaaDaaaleaacqqHJoWuaeaacqaIYaGmaaGccqGHRaWkdaWfqaqaaiGbcggaHjabckhaYjabcEgaNjGbc2gaTjabcMgaPjabc6gaUbWcbaGae83SdC2aaSbaaWqaaiabigdaXaqabaWccqGHiiIZcqGFge=rdaWgaaadbaGaeGymaedabeaaaSqabaGcdaqbdaqaaiab9PfaAjabgkHiTiabdseaenaaBaaaleaacqWG6bGEaeqaaOGae83SdC2aaSbaaSqaaiabigdaXaqabaaakiaawMa7caGLkWoadaqhaaWcbaGaeu4Odm1aaSbaaWqaaiabdQha6bqabaaaleaacqaIYaGmaaaakiaawIcacaGLPaaaaeaaaeaacqGH9aqpaeaacqWGdbWqdaWfqaqaaiGbcggaHjabckhaYjabcEgaNjGbc2gaTjabcMgaPjabc6gaUbWcbaGae83SdC2aaSbaaWqaaiabigdaXaqabaWccqGHiiIZcqGFge=rdaWgaaadbaGaeGymaedabeaaaSqabaGcdaqbdaqaaiab9PfaAjabgkHiTiabdseaenaaBaaaleaacqWG6bGEaeqaaOGae83SdC2aaSbaaSqaaiabigdaXaqabaaakiaawMa7caGLkWoadaqhaaWcbaGaeu4Odm1aaSbaaWqaaiabdQha6bqabaaaleaacqaIYaGmaaGccqGGUaGlaaaaaa@9F1F@

For all γ0 G MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaaliab=zq8hbaa@3837@ 0, C γ0 = 0 and D z γ0 = 0, so the area of minimisation can be extended,

δ ^ = C arg min γ Z D z γ Σ z 2 , MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF0oazgaqcaiabg2da9iabdoeadnaaxababaGagiyyaeMaeiOCaiNaei4zaCMagiyBa0MaeiyAaKMaeiOBa4galeaacqWFZoWzaeqaaOWaauWaaeaaieqacqGFAbGwcqGHsislcqWGebardaWgaaWcbaGaemOEaOhabeaakiab=n7aNbGaayzcSlaawQa7amaaDaaaleaacqqHJoWudaWgaaadbaGaemOEaOhabeaaaSqaaiabikdaYaaakiabcYcaSaaa@4991@

which is the MLE of δ based on Z by definition.    □

Remark 1 Using the invertible map between any two choices of Y, Y and Y', as defined in Step 1 above, the respective maximum likelihood estimates of α, Σ y and Σ y' can be shown to be identical based on Y or Y'. In this sense, the choice of μ ˜ 0 MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaaiiWacuWF8oqBgaacamaaCaaaleqabaGaeGimaadaaaaa@2F94@ is thus truly irrelevant.

Remark 2 Sometimes, additional linear combinations of γ can be assumed to be zero for most genes, C* γ= 0 for some matrix C* with rowspace being a superspace of the rowspace of C. Let P* be any projection matrix on the corresponding space V MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaat0uy0HwzTfgDPnwy1egaryqtHrhAL1wy0L2yHvdaiqaacqWFveVvaaa@384A@ * = {μ: μ= D γ, C* γ= 0, γ q} and let Y* = X - P* X. It is straight forward to show that a variant of the Proposition still holds, so given the covariance structure matrices the inference results concerning C γ will be identical, based on Y or Y* respectively. However, the estimates of the covariance structure matrices for Y and Y* might not be coherent and the results are expected to differ slightly.

The estimator of power

Consider a certain experimental design, a level 1-λ test and a differential expression δ. Let a realisation of the experiment be given, which e.g. results in certain quality deviations between arrays. The conditional power is defined as the probability of identifying a random gene in the current experiment, i.e. conditional on e.g. the quality deviations, when the gene has differential expression δ. The power is then defined as the average conditional power over all possible realisations of the experimental design. The power is thus related to an unperformed experiment while the conditional power is related to a specific performed experiment. Here, the test is assumed to be valid conditional on the experiment, i.e. to have conditional power λ when δ= 0.

In Evaluation of power, the aim is to estimate the power for a hypothetical experiment where the distribution of the data under the null hypothesis is obtained by resampling of real data. For a given resample, a constant differential expression is added to randomly selected genes and the statistics t g are computed. The estimate of the conditional critical value is computed so that a proportion λ of the unregulated genes satisfy |t g | ≥ . The conditional power is then estimated by the proportion of regulated genes satisfying |t g | ≥ . The power is finally estimated by averaging the estimated conditional power over the resamples.

References

  1. Auer H, Lyianarachchi S, Newsom D, Klisovic M, Marcucci G, Kornacker K: Chipping away at the chip bias: RNA degradation in microarray analysis. Nature Genetics. 2003, 35: 292-293. 10.1038/ng1203-292.

    Article  CAS  PubMed  Google Scholar 

  2. Kristiansson E, Sjögren A, Rudemo M, Nerman O: Weighted Analysis of Paired Microarray Experiments. Stat Appl Genet Mol Biol. 2005, 4: Article30-10.2202/1544-6115.1160.

    PubMed  Google Scholar 

  3. Kristiansson E, Sjögren A, Rudemo M, Nerman O: Quality Optimised Analysis of General Paired Microarray Experiments. Stat Appl Genet Mol Biol. 2006, 5: Article10-10.2202/1544-6115.1209.

    PubMed  Google Scholar 

  4. Arnold S: The Theory of Linear Models and Multivariate Analysis. 1980, John Wiley & Sons

    Google Scholar 

  5. Lönnstedt I, Speed T: Replicated microarray data. Statistica Sinica. 2002, 12: 31-46.

    Google Scholar 

  6. Smyth G: Linear models and empirical Bayes methods for assessing differential expression in microarray experiments. Stat Appl Genet Mol Biol. 2004, 3: Article3-

    PubMed  Google Scholar 

  7. Ritchie M, Diyagama D, Neilson J, van Laar R, Dobrovic A, Holloway A, Smyth G: Empirical array quality weights in the analysis of microarray data. BMC Bioinformatics. 2006, 7: 261-10.1186/1471-2105-7-261.

    Article  PubMed Central  PubMed  Google Scholar 

  8. R Development Core Team: R: A Language and Environment for Statistical Computing. 2006, R Foundation for Statistical Computing, Vienna, Austria, [ISBN 3-900051-07-0], [http://www.R-project.org]

    Google Scholar 

  9. The WAME homepage. [http://wame.math.chalmers.se]

  10. Edgar R, Domrachev M, Lash A: Gene Expression Omnibus: NCBI gene expression and hybridization array data repository. Nucleic Acids Research. 2002, 30: 207-210. 10.1093/nar/30.1.207.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  11. Barth A, Merk S, Arnoldi E, Zwermann L, Kloos P, Gebauer M, Steinmeyer K, Bleich M, Kaab S, Hinterseer M, Kartmann H, Kreuzer E, Dugas M, Steinbeck G, Nabauer M: Reprogramming of the Human Atrial Transcriptome in Permanent Atrial Fibrillation. Circulation Research. 2005, 96 (9): 1022-1029. 10.1161/01.RES.0000165480.82737.33.

    Article  CAS  PubMed  Google Scholar 

  12. Spira A, Beane J, Pinto-Plata V, Kadar A, Liu G, Shah V, Celli B, Brody J: Gene expression profiling of human lung tissue from smokers with severe emphysema. Am J Respir Cell Mol Biol. 2004, 31 (6): 601-10. 10.1165/rcmb.2004-0273OC.

    Article  CAS  PubMed  Google Scholar 

  13. Irizarry R, Hobbs B, Collin F, Beazer-Barclay Y, Antonellis K, Scherf U, Speed T: Exploration, Normalization, and Summaries of High Density Oligonucleotide Array Probe Level Data. Biostatistics. 2003, 4 (2): 249-64. 10.1093/biostatistics/4.2.249.

    Article  PubMed  Google Scholar 

  14. Gentleman R, Carey V, Bates D, Bolstad B, Dettling M, Dudoit S, Ellis B, Gautier L, Ge Y, Gentry J, Hornik K, Hothorn T, Huber W, Iacus S, Irizarry R, Friedrich L, Li C, Maechler M, Rossini A, Sawitzki G, Smith C, Smyth G, Tierney L, Yang J, Zhang J: Bioconductor: Open software development for computational biology and bioinformatics. Genome Biology. 2004, 5: R80-10.1186/gb-2004-5-10-r80.

    Article  PubMed Central  PubMed  Google Scholar 

  15. Tomlins S, Rhodes D, Perner S, Dhanasekaran S, Mehra R, Sun X, Varambally S, Cao X, Tchinda J, Kuefer R, Lee C, Montie J, Shah R, Pienta K, Rubin M, Chinnaiyan A: Recurrent Fusion of TMPRSS2 and ETS Transcription Factor Genes in Prostate Cancer. Science. 2005, 310: 644-648. 10.1126/science.1117679.

    Article  CAS  PubMed  Google Scholar 

  16. van Wieringen W, van de Wiel M, van der Vaart A: A Test for Partial Differential Expression. Tech Rep WS2006-4. 2006, Department of Mathematics, Vrije Universiteit

    Google Scholar 

  17. Tibshirani R, Hastie T: Outlier sums for differential gene expression analysis. Biostatistics. 2007, 8 (1): 2-8. 10.1093/biostatistics/kxl005.

    Article  PubMed  Google Scholar 

  18. Klebanov L, Yakovlev A: Treating Expression Levels of Different Genes as a Sample in Microarray Data Analysis: Is it Worth a Risk?. Stat Appl Genet Mol Biol. 2006, 5: Article9-

    PubMed  Google Scholar 

  19. Pawitan Y, Calza S, Ploner A: Estimation of the false discovery proportion under general dependence. Bioinformatics. 2006, 22 (24): 3025-3031. 10.1093/bioinformatics/btl527.

    Article  CAS  PubMed  Google Scholar 

  20. Anton H: Elementary Linear Algebra. 1991, Wiley, 6

    Google Scholar 

Download references

Acknowledgements

Lina Gunnarsson is acknowledged for providing feedback on the manuscript. AS and EK acknowledge the National Research School in Genomics and Bioinformatics for funding.

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Anders Sjögren.

Additional information

Competing interests

The author(s) declares that there are no competing interests.

Authors' contributions

ON provided initial ideas related to the generalisation. AS formulated the generalisation, performed the proofs in Methods, designed and programmed the analyses, simulations and plots. EK and ON helped refining the generalisation and Methods. AS, EK and MR drafted the manuscript. All authors continuously provided feedback on various parts of the work leading to the manuscript and approved the final version of the manuscript.

Electronic supplementary material

12859_2007_1759_MOESM1_ESM.png

Additional file 1: Pairwise plots of all arrays in the Atrium dataset. Transformed expression values for all arrays in the Atrium dataset. See legend of Figure 2 for details. (PNG 3 MB)

12859_2007_1759_MOESM2_ESM.png

Additional file 2: Pairwise plots of all arrays in the COPD dataset. Transformed expression values for all arrays in the COPD dataset. See legend of Figure 2 for details. (PNG 3 MB)

12859_2007_1759_MOESM3_ESM.pdf

Additional file 3: Probability plots for the Atrium dataset. Empirical distributions of p-values for LIMMA, weighted LIMMA, OLM and WAME from tests on 100 resamples from the Atrium dataset. Average empirical distribution indicated. Since no signal is added, the curves should ideally follow the diagonal. (PDF 487 KB)

Authors’ original submitted files for images

Rights and permissions

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and permissions

About this article

Cite this article

Sjögren, A., Kristiansson, E., Rudemo, M. et al. Weighted analysis of general microarray experiments. BMC Bioinformatics 8, 387 (2007). https://doi.org/10.1186/1471-2105-8-387

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/1471-2105-8-387

Keywords