Statistical Tests

The MULTTEST Procedure

Statistical Tests

The following section discusses the statistical tests performed in the MULTTEST procedure. For continuous data, a t-test for the mean is available. For discrete variables, available tests are the Cochran-Armitage (CA) linear trend test, the Freeman-Tukey (FT) double arcsine test, the Peto mortality-prevalence test, and the Fisher exact test.

Throughout this section, the discrete and continuous variables are denoted by S_vgsr and X_vgsr, respectively, where v is the variable, g is the treatment group, s is the stratum, and r is the replication. A plus sign (+) subscript denotes summation over an index.

Cochran-Armitage Linear Trend Test

The Cochran-Armitage linear trend test (Cochran 1954; Armitage 1955; Agresti 1990) is implemented using a Z-score approximation, an exact permutation distribution, or a combination of both.

Z-Score Approximation

Let m_vgs denote the sample size for a binary variable v within group g and stratum s. The pooled probability estimate for variable v and stratum s is

p_vs = [(S_v+s+)/(m_v+s)]

The expected value (under constant within-stratum treatment probabilities) for variable v, group g, and stratum s is

E_vgs = m_vgs p_vs

The test statistic for variable v has numerator

$N_v = \sum_{s} \sum_{g} t_g (S_{vgs+} - E_{vgs})$

where t_g denotes a trend coefficient (specified by the CONTRAST statement). The binomial variance estimate for this statistic is

$V_v = \sum_{s} p_{vs} (1 - p_{vs}) \sum_{g} m_{vgs} (t_s - {\bar{t}}_{vs})^2$

where

${\bar{t}}_{vs} = \sum_{g} \frac{m_{vgs} t_g}{m_{v+s}}$

The hypergeometric variance estimate (the default) is

$V_v = \sum_{s} \{m_{v+s}/(m_{v+s}-1)\}p_{vs} (1 - p_{vs}) \sum_{g} m_{vgs} (t_s - {\bar{t}}_{vs})^2$

For any strata s with $m_{v+s} \le 1$ , the contribution to the variance is taken to be zero.

PROC MULTTEST computes the Z-score statistic

$Z_v = \frac{N_v}{\sqrt{V_v}}$

The p-value for this statistic comes from the standard normal distribution. Whenever a 0 is computed for the denominator, the p-value is set to 1. This p-value approximates the probability obtained from the exact permutation distribution, discussed in the following text.

The Z-score statistic can be continuity-corrected to better approximate the permutation distribution. With continuity correction c, the upper-tailed p-value is computed from

$Z_v = \frac{N_v - c}{\sqrt{V_v}}$

For two-tailed, noncontinuity-corrected tests, PROC MULTTEST reports the p-value as 2 min(p, 1 - p), where p is the upper-tailed p-value. The same formula holds for the continuity-corrected test, with the exception that when the noncontinuity-corrected Z and the continuity-corrected Z have opposite signs, the two-tailed p-value is 1.

When the PERMUTATION= option is specified and no STRATA variable is specified, PROC MULTTEST uses a continuity correction selected to optimally approximate the upper-tail probability of permutation distributions with smaller marginal totals (Westfall and Lin 1988). Otherwise, the continuity correction is specified using the CONTINUITY= option in the TEST statement.

The CA Z-score statistic is the Hoel-Walburg (Mantel-Haenszel) statistic reported by Dinse (1985).

Exact Permutation Test

When you use the PERMUTATION= option for CA in the TEST statement, PROC MULTTEST computes the exact permutation distribution of the trend score

$T_v = \sum_{s}\sum_{g} t_g S_{vgs+}$

and then compares the observed value of this trend with the permutation distribution to obtain the p-value

$p_v = \Pr(X \ge { observed } T_v)$

where X is a random variable from the permutation distribution and where upper-tailed tests are requested. This probability can be viewed as a binomial probability, where the within-stratum probabilities are constant and where the probability is conditional with respect to the marginal totals S_v+s+. It also can be considered a rerandomization probability.

Because the computations can be quite time-consuming with large data sets, specifying the PERMUTATION=number option in the TEST statement limits the situations where PROC MULTTEST computes the exact permutation distribution. When marginal total success or total failure frequencies exceed number for a particular stratum, the permutation distribution is approximated using a continuity-corrected normal distribution. You should be cautious in using the PERMUTATION= option in conjunction with bootstrap resampling because the permutation distribution is recomputed for each bootstrap sample. This recomputation is not necessary with permutation resampling.

The permutation distribution is computed in two steps:

The permutation distributions of the trend scores are computed within each stratum.
The distributions are convolved to obtain the distribution of the total trend.

As long as the total success or failure frequency does not exceed number for any stratum, the computed distributions are exact. In other words, if $S_{v+s+} \le$ number or $(m_{v+s} - S_{v+s+}) \le$ number for all s, then the permutation trend distribution for variable v is computed exactly.

In step 1, the distribution of the within-stratum trend

$\sum_{g} t_g S_{vgs+}$

is computed using the multivariate hypergeometric distribution of the S_vgs+, provided number is not exceeded. This distribution can be written as

$\Pr(S_{v1s+}, S_{v2s+}, ... , S_{vGs+}) = \prod^G_{g=1} \frac{ ( m_{vgs} \ S_{vgs+} ) } { ( m_{v+s} \ S_{v+s+} ) }$

The distribution of the within-stratum trend is then computed by summing these probabilities over appropriate configurations. For further information on this technique, refer to Bickis and Krewski (1986) and Westfall and Lin (1988). In step 2, the exact convolution distribution is obtained for the trend statistic summed over all strata having totals that meet the threshold criterion. This distribution is obtained by applying the fast Fourier transform to the exact within-stratum distributions. A description of this general method can be found in Pagano and Tritchler (1983) and Good (1987).

The convolution distribution of the overall trend is then computed by convolving the exact distribution with the distribution of the continuity-corrected standard normal approximation. To be more specific, let S₁ denote the subset of stratum indices that satisfy the threshold criterion, and let S₂ denote the subset of indices that do not satisfy the criterion. Let T_v1 denote the combined trend statistic from the set S₁, which has an exact distribution obtained using Fourier analysis as previously outlined, and let T_v1 denote the combined trend statistic from the set S₂. Then the distribution of the overall trend T_v = T_v1 + T_v2 is obtained by convolving the analytic distribution of T_v1 with the continuity-corrected normal approximation for T_v2. Using the notation from the "Z-Score Approximation" section, this convolution can be written as

$\Pr(T_{v1} + T_{v2} \ge u) & = & \sum_{u1} \Pr(T_{v1} + T_{v2} \ge u | T_{v1} = u1) \Pr(T_{v1} = u1) \ & \approx & \sum_{u1} \Pr(Z \ge z) \Pr(T_{v1} = u1) \$

where Z is a standard normal random variable, and

$z = \frac{1}{\sqrt{V_v}} (u - u1 - \sum_{S_2} p_{vs} \sum_{g} t_g m_{vgs} - c )$

In this expression, the summation of s in V_v is over S₂, and c is the continuity correction discussed under the Z-score approximation.

When a two-tailed test is requested, the expected trend

$E_v = \sum_{s} \sum_{g} t_g E_{vgs}$

is computed, and the two-tailed p-value is reported as the permutation tail probability for the observed trend T_v plus the permutation tail probability for 2E_v - T_v, the reflected trend.

Freeman-Tukey Double Arcsine Test

For this test, the trend scores t₁, ... , t_G are centered to the values c₁, ... , c_G, where $c_g = t_g - \bar{t}$ , $\bar{t} = \sum_{g} t_{g} / G$ , and G is the number of groups. The numerator of this test statistic is

$N_v = \sum_{s} m_{v+s} \sum_{g} c_{g} f(S_{vgs+} , m_{vgs})$

and is weighted by the within-strata sample size (m_v+s) to ensure comparability with the ordinary C-A trend statistic.

The function f(r,n) is the double arcsine transformation:

$f(r,n) = \arcsin ( \sqrt{ \frac{r}{n+1} } ) + \arcsin ( \sqrt{ \frac{r+1}{n+1} } )$

The variance estimate is

$V_v = \sum_{s} m_{v+s}^2 \sum_{g} \frac{c_g^2}{m_{vgs} + \frac{1}2}$

and the test statistic is

$Z_v = \frac{N_v}{\sqrt{V_v}}$

The Freeman-Tukey transformation and its variance are described by Freeman and Tukey (1950) and Miller (1978). Since its variance is not weighted by the pooled probabilities, as is the CA test, the FT test can be more useful than the CA test for tests involving only a subset of the groups.

Peto Mortality-Prevalence Test

The Peto test is a modified Cochran-Armitage procedure incorporating mortality and prevalence information. It represents a special case in PROC MULTTEST because the data structure requirements are different, and the resampling methods used for adjusting p-values are not valid. The TIME= option variable is required to specify "death" times or, more generally, time of occurrence. In addition, the test variables must assume one of the following three values.

0 = no occurrence
1 = incidental occurrence
2 = fatal occurrence

Use the TIME= option variable to define the mortality strata, and use the STRATA statement variable to define the prevalence strata.

The Peto test is computed like two Cochran-Armitage Z-score approximations, one for prevalence and one for mortality.

In the following notation, the subscript v represents the variable, g represents the treatment group, s represents the stratum, and t represents the time. Recall that a plus sign (+) in a subscript location denotes summation over that subscript.

Let S^P_vgs be the number of incidental occurrences, and let m^P_vgs be the total sample size for variable v in group g, stratum s, excluding fatal tumors.

Let S^F_vgt be the number of fatal occurrences in time period t, and let m^F_vgt be the number alive at the end of time t-1.

The pooled probability estimates are

$p^P_{vs} & = & \frac{S^P_{v+s}}{m^P_{v+s}} \ p^F_{vt} & = & \frac{S^F_{v+t}}{m^F_{v+t}} \$

The expected values are

$E^P_{vgs} & = & m^P_{vgs} p^P_{vs} \ E^F_{vgt} & = & m^F_{vgt} p^F_{vt} \$

Define the numerator terms:

$N_v^P & = & \sum_{s}\sum_{g} t_g ( S^P_{vgs} - E^P_{vgs} ) \ N_v^F & = & \sum_{t}\sum_{g} t_g ( S^P_{vgt} - E^P_{vgt} ) \$

where t_g denotes a trend coefficient. Define the denominator variance terms (using the binomial variance) :

$V_v^P & = & \sum_{s} p^P_{vs} (1 - p^P_{vs}) [ ( \sum_{g} m^P_{vgs} {t_g}^2 ) ... ...{g} m^F_{vgt} {t_g}^2 ) - \frac{1}{m^F_{v+t}} ( \sum_{g} m^F_{vgt} t_g )^2 ] \$

The hypergeometric variances (the default) are calculated by weighting the within-strata variances as discussed in the "Z-Score Approximation" section.

The Peto statistic is computed as

$Z_v &=& \frac{N^P_v + N^F_v -c} {\sqrt{V^P_v + V^F_v}}$

where c is a continuity correction. The p-value is determined from the standard normal distribution unless the PERMUTATION=number option is used. When you use the PERMUTATION= option for PETO in the TEST statement, PROC MULTTEST computes the ``discrete approximation'' permutation distribution described by Mantel (1980) and Soper and Tonkonoh (1993). Specifically, the permutation distribution of

$\sum_s \sum_g t_g S^P_{vgs} + \sum_t \sum_g t_g S^F_{vgt}$

is computed, assuming that $\{\sum_g t_g S^P_{vgs}\}$ and $\{\sum_g t_g S^F_{vgt}\}$ are independent over all s and t. The p-values are exact under this independence assumption. However, the independence assumption is valid only asymptotically, which is why these p-values are called "approximate."

An exact permutation distribution is available only under the assumption of equal risk of censoring in all treatment groups; even then, computing this distribution can be cumbersome. Soper and Tonkonoh (1993) describe situations where the discrete approximation distribution closely fits the exact permutation distribution.

Fisher Exact Test

The CONTRAST statement in PROC MULTTEST enables you to compute Fisher exact tests for two-group comparisons. No stratification variable is allowed for this test. Note, however, that the FISHER exact test is a special case of the exact permutation tests performed by PROC MULTTEST and that these permutation tests allow a stratification variable. Recall that contrast coefficients can be -1, 0, or 1 for the Fisher test. The frequencies and sample sizes of the groups scored as -1 are combined, as are the frequencies and sample sizes of the groups scored as 1. Groups scored as 0 are excluded. The -1 group is then compared with the 1 group using the Fisher exact test.

Letting x and m denote the frequency and sample size of the 1 group, and y and n denote those of the -1 group, the p-value is calculated as

$\Pr(X \ge x | X + Y = x + y) = \sum^m_{i=x} \frac{ ( m \ i ) ( n \ x + y - i ) } { ( m + n \ x + y ) }$

where X and Y are independent binomially distributed random variables with sample sizes m and n and common probability parameters. The hypergeometric distribution is used to determine the stated probability; Yates (1984) discusses this technique. PROC MULTTEST computes the two-tailed p-values by adding probabilities from both tails of the hypergeometric distribution. The first tail is from the observed x and y, and the other tail is chosen so that the resulting probability is as large as possible without exceeding the probability from the first tail.

t-Test for the Mean

For continuous variables, PROC MULTTEST automatically centers the trend coefficients, as in the Freeman-Tukey test. These centered coefficients c_g are then used to form a t-statistic contrasting the within-group means. Let n_vgs denote the sample size within group g and stratum s; it depends on variable v only when there are missing values. Define

${\bar{X}}_{vgs+} = \frac{1}{n_{vgs}} \sum_{r}X_{vgsr}$

as the sample mean within a group-and-stratum combination, and define

$s_v^2 = \frac{\displaystyle\sum_{s}\sum_{g}\sum_{r} (X_{vgsr}-{\bar{X}}_{vgs+} )^2 } {\displaystyle\sum_{s}\sum_{g}( n_{vgs} - 1 )}$

as the pooled sample variance. Assume constant variance for all group-and-stratum combinations. Then the t-statistic for the mean is

$M_v = \frac{\displaystyle\sum_{s} n_{v+s} \sum_{g}c_g {\bar{X}}_{vgs+}} {\displaystyle\sqrt{ s_v^2 (\sum_{s} n_{v+s}^2 \sum_{g}\frac{c_g^2}{n_{vgs}} ) } }$

and is weighted by the within-strata sample size (n_v+s) to ensure comparability with the C-A trend and Freeman-Tukey statistics.

Let $\mu_{vgs}$ denote the treatment means. Then under the null hypothesis that

$\sum_{s} n_{v+s} \sum_{g} c_g \mu_{vgs} = 0$

and assuming normality, independence, and homoscedasticity, M_v follows a t-distribution with $\sum_{s}\sum_{g} ( n_{vgs} - 1 )$ degrees of freedom.

Whenever a denominator of 0 is computed, the p-value is set to 1. When missing data force n_vgs = 0, then the contribution to the denominator of the pooled variance is 0 and not -1. This is also true for degrees of freedom.

Chapter Contents
Previous
Next
Top