,<':-
4"
).7',
STATISTICS& PROBABIUTY LEI"~RS ELSEVIER
Statistics & ProbabilityLetters 33 (1997) 333-339
A multivariate approach for estimating the random effects variance component in one-way random effects model Brajendra C. Sutradhar Department of Mathematics and Statistics, Memorial University of Newfoundland, St. John's, Nfld, Canada AIC 5S7 Received February 1996; revised April 1996
Abstract It is well known that the ANOVA estimator of the random effects variance component in one-way random effects model can assume negative values. It is also well known that nonnegative quadratic unbiased estimators do not exist for estimating the random effects variance component (LaMotte, 1973). LaMotte (1985) indicated the possibility that nonnegative invariant quadratic estimator of the random effects variance component uniformly better than the ANOVA estimator may exist for the balanced one-way random effects model. Mathew et al. (1992a) have shown that such estimator exists only when the number of treatments is 9 or less. As noted in Herbach (1959) (see also Thompson, 1962), a simple truncation of the ANOVA estimator at zero yields uniform improvement over the ANOVA estimator. The estimators suggested by Herbach and Thompson are, in fact, restricted maximum likelihood estimators, and they are nonquadratic by nature. In this paper, we discuss a multivariate technique which always yields positive estimate of the random effects variance component in one-way random effects model. The multivariate approach exploits the estimates of the eigenvalues of the covariance matrix of the model in estimating the variance components including the error variance. The resulting estimates are nonquadratic. The success of this multivariate approach depends on the precise estimation of the eigenvalues. Since there does not exist any unbiased estimation procedure in the small sample case for the estimation of the eigenvalues, we use a delete-d jackknife procedure to estimate them. This delete-d jackknife based multivariate approach yields better estimates (in terms of mean squared error) for the random effects variance component than the restricted maximum likelihood estimation as well as Chow and Shao's (1988) nonquadratic estimation approaches, which is shown through a simulation study for the cases with number of treatments up to 20. Keywords: Covariance matrix of the model; Eigenvalues; Positive estimates of variance components; Delete jackknife; Restricted maximum likelihood estimates; Monte-Carlo experiment
1. Introduction Consider the one-way random effect model Yij=g+ri+eij,
i = 1 . . . . . k; j =
1. . . . . n,
(1)
where, in the usual notation, Yij is the j t h observation on the ith treatment, # is the overall mean, zi is the random effect with mean zero and variance cry, and eij are error terms with mean zero and variance a 2. For 0167-7152/97/$17.00 (~) 1997 Elsevier Science B.V. All rights reserved PH S01 67-71 52 (96)001 44-7
B.C. Sutradhar I Statistics & Probability Letters 33 (1997) 333 339
334
i ~ i', i,i' = 1. . . . . k, zi and and j. Hence
Cov(Yij, Yi'j' ) =
/7it
are uncorrelated. Also
0 -2 -[- O"z2
for i = i', j = j ' ,
0-z2
for i = i', j ¢ j ' ,
0
for i ¢ i'.
Ti
and e are uncorrelated, i.e., E('ci~i,j) = 0 for all i,i'
(2)
As is well known, the usual ANOVA estimates of the variance components of the model (1) are ~2 = MSSE,
^2zA = (MSST - MSSE)/n, ~T
(3)
where k MSSE
n
=
(rij
-
L)2/k(n
- 1),
i-j j=l k
MSST = n ~
(Yi. - Y..)2/( k - 1),
,=1
Yi.
=~
YiJ/n"
j=l
Also it is well known that when normality assumptions are made for the treatment effect zi and the errors e~j of the model (1), the maximization of the log-likelihood function of p, cr2 and G2 provides the solutions ^2 GL = {(1 - 1/k)MSST - MSSE}/n
62 = MSSE,
(4)
for the variance components of the model. Both of the ANOVA (3) and maximum likelihood (4) estimators of the random effect variance component, 0-2, can assume negative values. Herbach (1959) suggested to use the maximum likelihood estimator 0-,u^2 = max{0,(1/n)[(1 - 1/k)MSST - MSSE]} (see also Thompson, 1962), ^2 which yields uniform improvement over the ANOVA estimator GA in the sense that GH^2has always smaller mean square error than 0-rA" ^2 This estimator, unlike the ANOVA estimator, is, however, nonquadratic by nature. Chow and Shao (1988, p. 351, Eq. (2.2)) suggested an alternative nonquadratic estimator, 0-rCS, ^2 given by
^2
k-
1
nk(k + 1)
MSST
k+l if MSST ~<~
MSSE,
O'rc S
~ I ~ 2 MSST - MSSE],
otherwise.
^2 ^2 ^2 has uniformaly smaller mean square error than GH as well as 0-~A. It was shown in Chow and Shao that Gcs Turning to the quadratic estimators, it is well known that nonnegative quadratic unbiased estimators do not exist for estimating the random effect variance component (LaMotte, 1973). Later on, LaMotte (1985) indicated the possibility that nonnegative invariant quadratic estimators uniformly better than OzA ^2 may exist for the one-way random effects model. Recently, Mathew et al. (1992a) have, however, shown that such estimators exist only when the number of treatments, k, is 9 or less. In fact, these authors considered a general balanced mixed linear model, and obtained necessary and sufficient conditions under which the ANOVA estimator of any random effect variance component can be uniformly improved (in terms of mean squared error) using a nonnegative invariant quadratic estimator. For similar results on the estimation of the variance components of the unbalanced mixed models, we refer to Mathew et al. (1992b).
B. C Sutradhar I Statistics & Probability Letters 33 (1997) 333-339
335
In the present paper, we propose certain nonquadratic estimators for the random effect variance component and the error variance of the one-way random effect model. The estimators are obtained by exploiting the estimates of the eigenvalues of the covariance matrix of the model. The proposed approach is quite novel as it always yields nonnegative nonquadratic estimates of the variance components. The success of this estimation technique, however, depends on the precise estimates of the eigenvalues. But, as there is no unbiased estimation technique available to estimate such eigenvalues in the small sample case, we use a delete-d jackknife procedure to obtain their estimates. This jackknife procedure may, however, yield negative estimate for the variance components of the random effects. A simple truncation of such negative estimates at zero appears to yield uniform improvement in the sense of smaller mean squared error (with substantially smaller bias and the same or slightly more standard error), over the nonquadratic estimates 6~H due to Herbach (1959) and ^2 a~cs due to Chow and Shao (1988). This we have verified in the paper through a simulation study for the small sample cases with number of treatment groups, k, up to 20.
2. Non-quadratic estimates of variance components: an eigenvalue approach It is important to emphasize on the fact that in the random effects model (1), the observations under a treatment are correlated. Let Yi = (yil . . . . . Yij . . . . . Yin)' denote the n-dimensional vector of observations under the ith treatment. Then it follows from (2) that for i = 1. . . . . k, yi has the n x n covariance matrix 2
(5)
S = a2In + a, Un,
where In is the n x n identity matrix and U, is the n x n unit matrix. It also follows from (2) that the n-dimensional vectors Yl . . . . Yi . . . . . Yk are uncorrelated. Let al/, be the ( j , f ) t h element of Z in Eq. (5) for all j , j ' = 1 . . . . . n. Since a11, = E ( Y i j - / z ) ( y i j , - / ~ ) for all i = 1. . . . . k, we may write a sample covariance matrix
s* = (~j,)n×~,
(6)
* , k corresponding to S, where sjj,, the ( j , j ' ) t h element of S*, is given by s j?., = ~-~i=1 (Yij-
f:.. )(Yij'-
Y.. ) / ( k - 1 ) ,
n
with )5. = ~ ' i k l ~ j = l yij/kn. The eigenvalues of this sample covariance matrix S* will be exploited to estimate the variance components a 2 and a 2. We mainly deal with the case k > n. This case has much appeal in estimating a2. This is because, as k, the number of treatments gets larger, one naturally gets more information about treatment variation. Note that by writing y = (Yl . . . . . Yi . . . . yk), an n x k dimensional observation matrix, one may express S* in (6) as (k - 1)S* = y y ' - ( k n ) - l ( l ' ~ y l k ) [ y l k
® l~n) + ( l ~ y ' ® ln)] + ( n 2 k ) - l ( l ' n y l k ) 2 U n ,
(7)
where In = (1, 1. . . . . 1)' : n x 1, and ® denotes the Kronecker product. Furthermore, one may write S* = S + S~, where (k - 1)S = y [ I k -- k - l lk l~]y,, and (k - 1)$1 = k - l [ y l k - n -1 ln(l~nylk)] [ylk -- n - I ln(l~nylk)] ', S being the usual sample covariance matrix.
(8)
B.C. Sutradhar I Statistics & Probability Letters 33 (1997) 333 339
336
Under the assumption that the joint distribution of the n-dimensional vectors Yl . . . . . Yi. . . . . Yk is absolutely continuous (with respect to nk-dimensional Lebesque measure), it follows from Okamoto (1973, Theorem, p. 764) that the eigenvalues of the sample covariance matrix S, as well as the eigenvalues of the covariance typed matrix & in (8), are positive and distinct with probability 1. Thus, the eigenvalues [1, [2 .... , In of the sample covariance matrix S* = S + Sa are also positive and distinct. Let [(1),[(2) . . . . . [(n) denote the ordered eigenvalues of S*. Since the largest eigenvalue of the matrix S in (5) can be shown as 21 = a z + na 2 and the remaining (n - 1) eigenvalues are 22 = ).3 . . . . . )~n = 2 = a z, we propose to use the estimate of the eigenvalues
21 = [0),
2=~
[(h)/(n - 1),
(9)
h=2
2 given by to obtain the estimates of a 2 and a T d2=
~,
0",2z ~--_ [ ~ l -
2]/n.
(10)
The estimates obtained from (10) are always positive. We remark here that for the k ~[2 ~> ' " ~->[n~>0. Let [(1), [(2) . . . . . [(r) be the ordered distinct nonnegative characteristic roots of S* with multiplicities ql, q2,..., qr respectively, with 1 ~r<~n and Y'~-h=lqh = n. Then, in a manner similar to that for the k > n, the estimates for a 2 and a 2T are given by
62= ~-~ qh[(h)/ ~-~ qh, h=2
(11)
h=2
and a~2 z ~--- [[(1) --a2]/n,
(12)
where d 2 is given by (11). The estimates obtained from (11) and (12) are always nonnegative.
3. Delete-d jackknife estimates of variance components The pseudo-Wishart matrix case k ~
sample. Since the d groups of observations may be deleted in t = ( k ) ways, u ranges from l to t. For the
B.C. Sutradhar / Statistics & Probability Letters 33 (1997) 333 339
uth jackknife sample, compute the eigenvalues f(l)~ > ' " and 2 as Zlu = f(l)u,
,~(u) =
A(h)u/(n
-
337
> f(,)~ of S(*~) and estimate the eigenvalues h 1
1).
(13)
h=2
We now compute the jackknife estimate of the bias of )~1 and ~ as ),lb = ( k - d ) ( , ~ l -
'~1.),
(14)
~(.)),
(15)
and ~b = ( k
- d)(~. -
where }.1. = Ztu=l follows:
).l./t,
and ,~(.) =
d~j = z + ~-b, ^2 J = [/~1 + ~lb O'rD
t ~u=l 2(u)/t,
which yield the delete-d jackknife estimates of a 2 and a~2 as
(16) •2Dj]/n.
(17)
^2
Since OrDJ can assume negative values, we truncate it to zero when it is negative, and obtain the delete-d
jackknife estimates of a 2 and a~2 as ,2 = O~j, ODj
(18)
,2 J = max[0, O'rD ^2 J ]. O'zD
(19)
These estimators are nonquadratic in nature. In the next section, we compare the performance of these estimators with other competitive nonquadratic estimators, for the small sample case, which is of practical interest.
4. Simulation study The parameters to be controlled in the Monte-Carlo study are as follows: (a) k, the number of treatments in the random effect model (1) which in the present set-up is the sample size; (b) n, the number of observations under each treatment which corresponds to the dimension of the observations in the sample of size k; (c) a 2, the error variation; (d) a 2, the random effects variance; (e) /~, the overall mean, which can be set to zero without any loss of generality. Note that irrespective of n~>2, the sample covariance matrix S* in (6) is more informative for a 2 when the sample size k is sufficiently large. We therefore consider k ~ oc as the asymptotic case, and choose k = 5, 10, 15 and 20 in the simulation study as the small sample cases. As the dimension n of E in (5) has little effect on the estimation of a 2, we consider small n, for convenience n = 2. We choose a z = 1.0, and 3 values of a~: 0.05, 0.1 and 0.15, for which A N O V A and MLE approaches might produce negative estimates for a~2 with large probability. For n = 2, for each pair of (a 2, a~), 2 we generate 5000 samples each of size k from the bivariate normal distribution with mean 0 and covariance matrix Z = a212 + a 2 U2. We then compute the 5000 values of azH ^2 , O'zC A2 S and O'.2 A2 is the TDJ using their formulas from Sections 1 and 3. Recall that OrH -2 s is the nonquadratic estimator of restricted maximum likelihood estimate of a~2 due to Herbach (1959), arc
338
B.C. Sutradhar / Statistics & Probability Letters 33 (1997) 333 339 Table 1 Empirical means and standard errors (SE) of ~r~H, O'~C A2 S and O',2 ~DJ for selected d's for the one-way random effects model with ~r2=l, n = 2, based on 5000 simulations Methods of estimation
,2 ffrDJ k
~ 5
0.05 0.10 0.15
10
0.05 0.10 0.15
15
0.05 0.10 0.15
20
0.05 0.10 0.15
a, H^2
a,cs^2
d= 1
d=2
d=3
d=4
Mean SE Mean SE Mean SE
0.1439 0.0676 0.1697 0.0850 0.1968 0.1043
0.1737 0.0604 0.1994 0.0763 0.2262 0.0942
0.3699 0.1686 0.3932 0.1925 0.4185 0.2201
0.1976 0.1147 0.2194 0.1364 0.2427 0.1657
0.0830 0.0623 0.0991 0.0788 0.1167 0.0984
Mean SE Mean SE Mean SE
0.1287 0.0399 0.1588 0.0515 0.1909 0.0644
0.1463 0.0358 0.1755 0.0467 0.2066 0.0588
0.2841 0.0788 0.3048 0.0910 0.3292 0.1051
0.1950 0.0683 0.2151 0.0816 0.2380 0.0972
0.1346 0.0550 0.1541 0.0680 0.1757 0.0833
0.0930 0.0423 0.1099 0.0542 0.1296 0.0681
Mean SE Mean SE Mean SE
0.1160 0.0277 0.1472 0.0368 0.1810 0.0469
0.1285 0.0250 0.1588 0.0336 0.1914 0.0433
0.2312 0.0488 0.2524 0.0568 0.2777 0.0668
0.1646 0.0436 0.1840 0.0529 0.2077 0.0641
0.1201 0.0365 0.1383 0.0460 0.1603 0.0574
0.0894 0.0297 0.1066 0.0387 0.1266 0.0497
Mean SE Mean SE Mean SE
0.1080 0.0218 0.1398 0.0293 0.1745 0.0376
0.1174 0.0199 0.1482 0.0271 0.1819 0.0351
0.2023 0.0355 0.2223 0.0423 0.2476 0.0505
0.1460 0.0324 0.1656 0.0398 0.1897 0.0491
0.1086 0.0275 0.1270 0.0381 0.1500 0.0445
0.0821 0.0228 0.0994 0.0300 0.1206 0.0392
,2 J is the delete-d jackknife estimate o f a~. 2 For k = 5, we computed o.z2 due to Chow and Shao (1988), and O'rD ,2
,2
the values of OzDJ with d = 1,2 and 3, and for each k = 10, 15 and 20, we computed the values of a~o J -2 ^2 s and O'~D ,2 J (for selected d's) with d = 1,2,3 and 4. The empirical means and standard deviations of arH, arc are summarized in Table 1. We also computed the empirical means and standard deviations of 62 under these three approaches, but not reported here. Herbach's, and Chow and Shao's procedures produce the same estimate of a 2 as in ANOVA. It was found that the A N O V A estimator of a 2, which is quadratic in nature, was uniformly better (in the sense of mean squared) than the delete-d jackknife nonquadratic estimator. For further improved estimator than the A N O V A estimator for a 2, we refer to Mathew et al. (1992c). Turning back to the performance of the delete-d jackknife estimator of the random effects variance component, it is .2 clear from Table 1 that for a suitable selection of d, the empirical means of aro J are substantially closer to ^2 ^2 the parameter values o f a~2 than those means o f a~H and arc s. For k = 5, 10, 15 and 20, the values o f d for ,2 ^2 ^2 which CrrDJ performs better than a~n and arc s are 3, 4, 4 and 4, respectively. Note that the selection of d does not satisfy the condition x/~ < d < k always, as mentioned in Efron and Tibsharani (1993). Furthermore, .2 these values of d are pretty close to x/~. For the above choices o f d, the standard deviations o f a~i~j are found to be the same or slightly more than the standard deviations of O'^2.cH and O'zC ^2 s . Thus, in the sense o f
B.C Sutradhar / Statistics & Probability Letters 33 (1997) 333 339
339
minimum square error, the delete-d jackknife procedure appears to perform better than its competitors, OrH ^2 ^2 due to Herbach (1959) and arc s due to Chow and Shao (1988), in estimating the random effects variance component a~.
5. Remarks Although the delete-d jackknife procedure, as compared to its competitors, was found to be superior in estimating the random effects variance component, the method is, however, not problem free. It was found that the delete-d jackknife procedure yields more negative estimates for a~ than the maximum likelihood or A N O V A estimates. The simple truncation o f the negative estimate to zero, yielded this better performance o f O'zD *: J as shown in Table 1 Next, the choice o f d in the delete-d jackknife procedure is arbitrary. For other selections o f a 2, in particular when a~2 is large, the choice o f d as shown in Table 1 m a y not be appropriate. Further investigation for the selection criterion o f d is necessary to make the delete-d jackknife useful procedure for the practitioners.
Acknowledgements The research was partially supported by a grant from the Natural Sciences and Engineering Research Council o f Canada. The author thanks Sean M a y for carrying out the computations. The author also thanks the referee for comments and suggestions.
References Chow, S.C. and J. Shao (1988), A new procedure for the estimation of variance components, Statist. Probab. Lett. 6, 349-355. Efron, B. and R.J. Tibshirani (1993), An Introduction to the Bootstrap (Chapman & Hall, New York). Herbach, L.H. (1959), Properties of type II analysis of variance test, Ann. Math. Statist. 30, 939-959. LaMotte, L.R. (1973), On non-negative quadratic unbiased estimation of variance components, J. Amer. Statist. Assoc. 68, 728-730. LaMotte, L.R. (1985), Admissibility, unbiasedness, and nonnegativity in the balanced random one-way ANOVA model, in: T. Calinski and W. Klonecki, eds., Linear Statistical Inference (Springer, New York) pp. 184-199. Mathew, T., B.K. Sinha and B.C. Sutradhar (1992a), Nonnegative estimation of variance components in general balanced mixed models, in: A.K. Md and E. Saleh, eds., Nonparametric Statistics and Related Topics (Elsevier, North-Holland) pp. 281-295. Mathew, T., B.K. Sinha and B.C. Sutradhar (1992b), Nonnegative estimation of variance components in unbalanced mixed models with two variance components, J. Multivariate Anal 42, 77-101. Mathew, T., B.K. Sinha and B.C. Sutradhar (1992c), Improved estimation of error variance in general balanced mixed models, Statist. Decisions 10, 227-238. Okamoto, M. (1973), Distinctness of the eigenvalues of a quadratic form in a multivariate sample, Ann. Statist. l, 763 765. Shao, J. and C.F.J. Wu (1989), A general theory for jackknife variance estimation, Ann. Statist. 17, 1176-1197. Sutradhar, B.C. and R.F. Bartlett (1990), An empirical study of the estimation of eigenvalues in connection with nonnegative estimation of variance components of a one-way random effect model, Comput. Statist. Data Anal 9, 21-35. Thompson, W.A. Jr. (1962), The problem of negative estimates of variance components, Ann. Math. Statist. 33, 273-289.