## Beginning of article

Suppose that an estimator X and an estimator V of its variance are available. Is the estimator V unbiased?

This question is frequently one aspect of computer simulation studies. A large number, N, of independent replicates ([X.sub.i], [V.sub.i]) are generated, and the average value of the [V.sub.i]'s, [bar]V, is to be used to test the hypothesis [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

If [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] were known, the test could be based on the statistic

(1) [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

which under the null hypothesis has an asymptotic standard normal distribution. Typically, however, neither [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] nor [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is known. The sample variance [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] can be substituted for [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] in the denominator of (1) without changing the limiting distribution of [Z.sub.1]; it would seem natural to substitute [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] for [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] in the numerator of (1) as well, and to consider the statistic

(2) [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

on the presumption that since [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is a consistent estimator of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], this statistic is a good approximation to [Z.sub.1] and should have the same limiting distribution. This, in fact, is not the case.

In Section 1, I show that Z is asymptotically distributed as a mean-zero normal random variable, but with nonunit variance. In Section 2, two examples illustrate the effect of incorrectly assuming that Z has a standard normal distribution.

Throughout this article, [S.sup.2], [[sigma].sup.2], and [theta] denote the sample variance, true variance, and mean for the random variable that subscripts them. In addition, [[tau].sup.2] denotes the difference between the fourth central moment and the square of the variance. All of these parameters are assumed to be finite.

1. ASYMPTOTIC DISTRIBUTION OF Z

Straightforward algebraic manipulations yield

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII],

where [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Under Ho: [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], the central limit theorem states that [Z.sub.1] is asymptotically distributed as a standard normal variate. [Z.sub.2] and [Z.sub.3] also have asymptotic standard normal distributions (this consequence of the U-statistic theorem is problem 3.2.6, p. 76 of Randles and Wolfe 1979).

Since [[alpha].sub.N] and [[beta].sub.N] converge in probability to 1 and 0, several applications of Slutzky's lemma (Bickel and Doksum 1977, p. 461) allow the conclusion that Z has the same limiting distribution as

[Z.sub.1] - [k.sub.1][Z.sub.2],

where [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

Note that [Z.sub.1] and [Z.sub.2] are not independent; their covariance is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

The asymptotic distribution of ([Z.sub.1], [Z.sub.2]) is bivariate normal with zero means, unit variances, and covariances equal to [k.sub.2] (see problem 73, p. 405 of Lehmann 1975). Thus [Z.sub.1] - [k.sub.1][Z.sub.2] (and consequently Z) is asymptotically normally distributed with mean zero and variance given by

(3) [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Expressed another way, the asymptotic variance of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is

(4) [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

which can be consistently estimated using moment estimators.

2. EXAMPLES

2.1 Estimation of the Mean for a Normal Population

Suppose that X is the mean of k [is less than] 1 independent and identically distributed normal random variables with mean [mu] and variance [[sigma].sup.2]. It is well known that V = [S.sup.2]/k is unbiased for the variance of X; nevertheless, it is instructive to consider the result of testing [H.sub.0]: E(V) = var(X) as described above.

Since the sample mean and variance are independent for normal random variables, the covariance term in (3) vanishes. X is normally distributed, so its fourth central moment is three times the square of its variance. Thus [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Since (k - 1)[S.sup.2]/[[sigma].sup.2] has a chi-squared distribution with k - 1 degrees of freedom,

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

From (3) it follows that var(Z) = k, rather than 1. In this case, treating Z as standard normal leads to a grossly inflated Type I error rate.

2.2 Maximum Likelihood Estimation of Parameters of a Multinomial Distribution

Let Y = ([Y.sub.1], [Y.sub.2], ..., [Y.sub.k])' be a multinomial random variable, describing the frequencies with which M independent trials result in each of k mutually exclusive and exhaustive outcomes. Suppose that f and [theta] are unknown constants, 0 [is less than] f, [theta] [is less than] 1, and let [p.sub.i] (the probability that an individual trial results in outcome i) be given by

[p.sub.i] = f[[theta].sup.i - 1] for i = 1, 2, ..., k - 1,

and [p.sub.k] = 1 - [p.sub.1] - [p.sub.2] - ... - [p.sub.k-1]. Maximum likelihood estimators (MLE's) f and [theta] are not available in closed form. Estimates of the variances of f and [theta] are taken from the estimated information matrix.

I conducted a simulation study to assess the performance of this estimation procedure. One thousand vectors Y were generated with M = 250, k = 7, f = .10, and [theta] = .70. For each, the MLE [[theta].sub.i] of [theta] and its estimated variance [V.sub.i] were obtained. For consistency with the notation of Section 1, let [X.sub.i] = [[theta].sub.i].

The mean estimate of [theta] was .70093, and the sample variance of these estimates was [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. The average estimate of variance was [bar]V = 2.897 x [10.sup.-3], with estimated standard error of [S.sub.V]/[square root of 1,000] = 1.142 x [10.sup.-5]. Computation of the Z statistic (2) yields a value of -9.88, which would be regarded as strong evidence that the information matrix estimate of var([theta]) is too small, if it is compared against standard normal quantiles.

Using method of moments estimates of the quantities on the right side of (3), the standard deviation of Z was estimated as 12.50, so the value Z = -9.88 is not significant evidence against the unbiasedness of V.

3. DISCUSSION

Under [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], the statistic

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

has an asymptotic standard normal distribution. The sample variance [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] can be substituted for [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] denominator without changing the limiting distribution of [Z.sub.1], but substitution of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] for [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] in the numerator of [Z.sub.1] changes the asymptotic variance by a factor of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. Thus, when [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is unknown [H.sub.0] should be tested using the statistic

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII],

where [[psi].sup.2] is a consistent estimator of

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

Sample sizes typically encountered in simulation studies should be adequate to ensure the normal approximation to the distribution of [Z.sub.4]. For example, N = 1,000 pairs (X, V) were generated according to the specifications of the example in Section 2.1, and used to calculate the statistic [Z.sub.4]. This procedure was replicated R = 1,000 times. The observed rejection rates of [H.sub.0] at nominal levels [alpha] = .05 and .10 (in favor of the two-tailed alternative) were .046 and .101. A similar simulation, under the conditions of Section 2.2, with N = 250 and R = 100, yielded observed rejection rates of .05 and .09.

The foregoing discussion has assumed that [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is estimated using paired observations ([X.sub.i], [V.sub.i]), i = 1, 2, ..., N. A referee has suggested that alternative sampling schemes might be most cost effective, in terms of the variance of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] as a function of the cost in time for data generation. Three options for data generation are:

1. Generate [N.sub.1] pairs (X, V) and, independently, [N.sub.2] X's. Use all ([N.sub.1] + [N.sub.2]) X's to compute [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

2. Generate [N.sub.1] V's and, independently, [N.sub.2] X's.

3. Generate N pairs (X, V).

The variance of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], under Option (1),

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], under Option (2),

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], under Option (3),

where [rho] denotes the correlation between [(X - [[theta].sup.x]).sup.2] and V.

If the costs in computing time for generating X, V, and (X, V) are a, b, and c, then a + b [is greater than or equal to] c, with a, b [is less than] c. The costs in computing time are [N.sub.1]c + [N.sub.2]a for Option (1), [N.sub.1]b + [N.sub.2]a for Option (2), and Nc for Option (3). The relative performance of the three options depends on a/c, b/c, [rho], and [theta] = [[sigma].sub.v]/[[tau].sub.x].

Options (1) and (2) can be compared to Option (3) by minimizing the associated variances subject to the restriction that the total cost be Nc. Thus the minimum variance under Option (2), divided by the variance under Option (3), is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII].

Calculation of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is more complicated; I used numerical minimization for the calculations that follow.

To get an idea of the relative values of a, b, and c, I generated a million pairs (X, V), a million independent X's, and a million independent V's, according to the example in Section 2.1, with k = 10. These three tasks required 302.97 seconds, 261.17 seconds, and 288.58 seconds using the GAUSS programming language on a 486/33 IBM compatible PC. Thus a/c = .862 and b/c = .953. Under the conditions of the example in Section 2.1, [rho] = 0 and [theta] = 1/3. In this case, the efficiency of Option (3) relative to Option (1) is 1.433; relative to Option (2) it is 1.415.

Under the conditions of the example in Section 2.2, a = b = c; the parameter estimate and its associated variance estimator are obtained simultaneously. Bootstrap estimates of [theta] and [rho] are .0793 (.0043) and .2610 (.0496). Using these estimates, the efficiency of Option (3) relative to Option (1) is 1.164; relative to Option (2) it is 1.207.

None of the three options is uniformly the best (for example, with a/c = b/c = .60, and [theta] = .20, the most efficient option is #2, #1, #3, as [rho] [is less than] .33, [rho] [epsilon] (.33, .65), and [rho] [is greater than] .65). However, a/c and b/c will usually be close to 1, so the most relevant comparison is when a = b = c. In this case

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII],

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII],

each of which is greater than one for all [rho][epsilon](-1, 1). Thus when a = b = c, Option (3) is uniformly more efficient than either of the other two options.

If a/c or b/c is less than one, the best choice of [N.sub.1] and [N.sub.2] under Options (1) and (2) will depend on the unknown parameters [theta] and [rho]. The effort necessary to obtain planning values of [theta] and [rho] so as to obtain the best values of [N.sub.1] and [N.sub.2] will probably offset whatever gains could be obtained using Options (1) or (2) rather than (3). Option (3) seems to be the natural and best procedure to follow.

[Received April 1991. Revised May 1992.]

REFERENCES

Bickel, P. J., and Doksum, K. A. (1977), Mathematical Statistics: Basic Ideas and Selected Topics, San Francisco: Holden-Day.

Lehmann, E. L. (1975), Nonparametrics: Statistical Methods Based on Ranks, New York: John Wiley.

Randles, R. H., and Wolfe, D. A. (1979), Introduction to the Theory of Nonparametric Statistics, New York: John Wiley.

* William A. Link is Mathematical Statistician, U.S. Fish and Wildlife Service, Patuxent Wildlife Research Center, Laurel, MD 20708. The author thanks R. M. Korwar, Kongming Wang, J. D. Nichols, and J. R. Sauer for comments on an earlier version of this article. Two anonymous reviewers and an associate editor are also thanked for helpful comments and suggestions which shaped its final form.