Search
   >> Home >> Resources & support >> FAQs >> Stratified generalized Wilcoxon test
This FAQ is relevant for users of Stata 9.

Why do Stata and SAS differ in the results they report for the stratified generalized Wilcoxon test for time-to-event data?

Title   Stratified generalized Wilcoxon test
Author Ronald Thisted, Departments of Statistics, Health Studies, and Anesthesia & Critical Care, The University of Chicago
Date January 1999; updated June 2005

This FAQ details the reason we believe Stata and SAS differ in the results that they report for the stratified Wilcoxon test. The context is time-to-event data from a stratified analysis such as a randomized multicenter clinical trial comparing an active treatment with a control, stratifying by (controlling for) site-specific differences in underlying event-time distribution. Multicenter randomized clinical trials comparing treatments often assess the treatment effect using the log-rank or a generalized Wilcoxon test, stratified by study site.

Although SAS does not have a procedure for calculating the standard stratified log-rank or Wilcoxon tests, PROC LIFETEST can be used to construct an approximate statistical test of treatment effect in the presence of stratification.

Stata, on the other hand, does provide the option strata() to its sts test command, which allows for the analysis of stratified data. It combines the stratum-specific estimates according to the procedure described for the standard stratified test by Collett (1994, section 2.7) and by Klein and Moeschberger (1997, section 7.5).

A particular example

Using a small published clinical trial dataset, we demonstrate the differences obtained by SAS and Stata and show that the use of PROC LIFETEST produces a different statistical significance of the treatment effect compared with the standard stratified test. The data, given by Lee (1992) and quoted by Collett (1994, p. 47, table 2.8), are reproduced below:

 _____________________________________________________
      21 - 40 	    40 - 60	         61+ 
 _______________	_______________    _______________
   
 BCG   c. parvum	BCG   c. parvum    BCG   c. parvum
 _____________________________________________________
  19     27*	 34*      8   	    10	    25*
  24*    21*          4      11*	     5	     8
   8     18*         17*     23*		    11*	
  17*    16*                 12*
  17*     7                  15*
  34*    12*                  8*
         24                   8*
         8            
         8*            
 _____________________________________________________
        An Asterisk (*) indicates a censored observation.

These data have been entered into Stata:

 . use http://www.stata.com/support/faqs/dta/parvum

 . list agegrp time c if treat==0, sepby(agegrp)
        
      +-------------------+
      | agegrp   time   c |
      |-------------------|
   1. |      1     19   1 |
   2. |      1     24   0 |
   3. |      1      8   1 |
   4. |      1     17   0 |
   5. |      1     17   0 |
   6. |      1     34   0 |
      |-------------------|
  16. |      2     34   0 |
  17. |      2      4   1 |
  18. |      2     17   0 |
      |-------------------|
  26. |      3     10   1 |
  27. |      3      5   1 |
      +-------------------+

 . list agegrp time c if treat==1, sepby(agegrp)

      +-------------------+
      | agegrp   time   c |
      |-------------------|
   7. |      1     27   0 |
   8. |      1     21   0 |
   9. |      1     18   0 |
  10. |      1     16   0 |
  11. |      1      7   1 |
  12. |      1     12   0 |
  13. |      1     24   1 |
  14. |      1      8   1 |
  15. |      1      8   0 |
      |-------------------|
  19. |      2      8   1 |
  20. |      2     11   0 |
  21. |      2     23   0 |
  22. |      2     12   0 |
  23. |      2     15   0 |
  24. |      2      8   0 |
  25. |      2      8   0 |
      |-------------------|
  28. |      3     25   0 |
  29. |      3      8   1 |
  30. |      3     11   0 |
      +-------------------+

The aim of the clinical trial was to compare the ability of two immunotherapies (BCG and Cryptosporidium parvum) to prolong the life of malignant melanoma patients. Thirty patients were randomized to one of the two treatment groups and then further classified into three age groups. This dataset is instructive because (a) it is small enough that the calculations can be replicated by hand, (b) it includes a small number of tied event times, and (c) the intermediate calculations are documented in Collett’s book.

After st setting the data, the command in Stata for comparing via the generalized Wilcoxon test the survival experience of the two treatment groups without stratification is sts test treat, wilcoxon, and correcting (stratifying) for age is: sts test treat, wilcoxon strata(agegrp) where treat is an indicator variable of treatment type and agegrp is the stratification variable. The results of these two commands are

 . stset time, failure(c)

 . sts test treat, wilcoxon                 

          failure _d:  c
    analysis time _t:  time
 
 
 Wilcoxon (Breslow) test for equality of survivor functions
 
       |   Events         Events        Sum of
 treat |  observed       expected        ranks
 ------+--------------------------------------
 0     |         5           3.71           34
 1     |         5           6.29          -34
 ------+--------------------------------------
 Total |        10          10.00            0
 
             chi2(1) =       0.91
             Pr>chi2 =     0.3397


 . sts test treat, wilcoxon strata(agegrp)  

          failure _d:  c
    analysis time _t:  time
 
 
 Stratified Wilcoxon (Breslow) test for equality of survivor functions
 
       |   Events         Events        Sum of
 treat |  observed    expected(*)      ranks(*)
 ------+--------------------------------------
 0     |         5           3.76            6
 1     |         5           6.24           -6
 ------+--------------------------------------
 Total |        10          10.00            0
 
 (*) sum over calculations within agegrp
 
             chi2(1) =       0.18
             Pr>chi2 =     0.6726

These results match those reported by Collett (1994, 47–49).

Now let’s analyze these data using SAS’s PROC LIFETEST. Using treatment group as the stratum designator and the STRATA statement, SAS produces precisely the calculations that Stata performs without the strata() option. That is, the SAS command

 PROC LIFETEST;TIME*C(0);STRATA AGEGRP; 

produces the same results as Stata’s

 . sts test treat, wilcoxon

When there is a separate stratification variable (such as agegrp in our example), the test statistics cannot be calculated directly using PROC LIFETEST. Because SAS is so widely used for clinical trial data management and analysis, the approach taken by users to test the treatment effect is often to specify the stratification variable using the STRATA command, and then to test the treatment effect with the TEST command, hoping the effect of ties is not critical. For our example, this procedure yields

         Univariate Chi-Squares for the WILCOXON Test

                Test       Standard                     Pr >
 Variable    Statistic    Deviation    Chi-Square    Chi-Square

 TREAT          0.9409       1.1889       0.6263       0.4287

The significance probability of 0.4287 is substantially different from the 0.6726 that Stata reports and Collett computes.

The reason for the discrepancy

The difference in the reported results is not due to programming errors. Rather, the difference hinges on how the stratum-specific estimates are calculated and combined to produce the overall statistic.

In the absence of strata, the log-rank and Wilcoxon statistics for the effect of an active treatment to a control are calculated as

       U = Σ cj(dj- wj)

where dj is the number of failures at the jth time point, wj is the expected number of failures in, say, the active treatment group, and cj = 1 for the log-rank test, and cj=nj, the number at risk of failure at the jth time point, for the Wilcoxon–Gehan–Breslow test (Kalbfleisch and Prentice 1980, 234; Collett 1994, 43). The test statistic U2/V, where V is the variance of the U, is compared with the chi-squared distribution with 1 df.

When there is a separate stratification variable (such as study site in multicenter randomized clinical trials or age group in our previous example), tests based on pooling data across all strata can be misleading because of the possibility of stratum-specific differences in survival unrelated to treatment. These differences often reflect differences in patient population, diagnostic methods, or case-mix across the study sites and would generally affect all treatments to an equal degree.

The rank-based statistics described above are easily generalized to the stratified case. As described in detail by Collett (1994, 47) or Klein and Moeschberger (1997, 205), the statistics Us and Vs are calculated within each stratum s, and the test statistic ((Σ Us)2)/Σ Vs is then compared with the chi-squared distribution with 1 df. This is what Stata does when the strata() option is specified.

Using Collett’s data from the previous section, let’s examine these calculations in detail.

                                                        Stratum-
          Manual                                        specific
          Calculation     STATA         SAS’s LIFETEST  Chi-squared
 Age      ------------   -----------  ----------------  -------------
 Stratum    U    V        U    V        U       VT      STATA    SAS
 --------------------------------------------------------------------
 21-40     -3. 155.615   -3. 155.615  -0.1804  0.7667  0.0578  0.0424
 41-60      5.  35.000    5.  35.000   0.4545  0.2975  0.7143  0.6944
 61+        4.  11.000    4.  11.000   0.6667  0.3492  1.4545  1.2728
 --------------------------------------------------------------------
 Total      6. 201.615    6. 201.615   0.9409  1.4134
                                                                 
 Chi-squared    0.1786        0.1786           0.6263
 p-value                      0.6726           0.4287
_______________________________________________________________

Let’s examine how SAS computes the overall test statistic. The U statistic and variance (VT) presented in the table for each stratum are those reported by PROC LIFETEST when the TEST statement is used to analyze each stratum independently. In other words, this is what SAS reports if we run PROC LIFETEST with the "TEST treatment" statement (no STRATA statement) on each stratum by itself. The overall chi-squared and corresponding p-values are those obtained when both the STRATA and the TEST statements are used. We can see that this statistic can be obtained by summing the stratum-specific values and then calculating the overall test statistic as [Σ(U)]2/Σ(VT). This is exactly what SAS does internally. So, SAS and Stata are constructing the overall test in the same way by summing the U statistics across strata, summing the variance, and then constructing the test statistic. The difference lies in what SAS and Stata are summing.

Within each stratum, the chi-squared value from SAS is slightly smaller than that from Stata; this is due to the slightly larger variance produced by the estimator used by SAS. Consequently, within each stratum, SAS produces a result that is (very slightly) more supportive of the null hypothesis of no difference. However, the overall chi-squared statistic from SAS is more than three times as large as the Gehan–Wilcoxon stratified chi-squared produced by Stata, and is much LESS supportive of the null hypothesis than is the Stata result.

According to the SAS manual, the TEST statement is intended to be used for testing the association between a list of “numeric covariates” and the failure time. It performs better, in terms of power, when there are few or no ties and when the variable to be tested is continuous. For a categorical variable, such as a stratification variable and tied event times, the test may produce undesirable results. The test statistic is obtained as the average of the statistics one would have obtained by breaking the ties in all possible ways. This averaging leads to inflated variance estimates, resulting in a more conservative test. The documentation for PROC LIFETEST suggests that the effect of tied event times on the results of the TEST command is likely to be small: “Unless the proportion of ties is large, it is unlikely that this will be a problem” (SAS, 1048). This effect can be seen directly in the last two columns of the table above. Although this may be true when only one stratification variable is present, the difference can be substantial when two or more strata variables are used. In their book, Klein and Moeschberger (1997, 200) recommend against doing this analysis using the TEST statement.

Unfortunately, SAS’s TEST statement when used with STRATA need not produce a conservative test, as the example above clearly demonstrates. The problem with SAS is that it implicitly uses an algorithm that gives inappropriate weights to the different strata. In the example above, the youngest age group has 15 observations, the oldest only 5. The SAS stratified test, however, gives only one-third the weight to the youngest age group that it gives to the oldest! To a close approximation, the U values obtained by SAS are equal to the U values obtained from Stata divided by N+1, where N is the stratum size, and the V values can be obtained by dividing the Stata variances by (N+1)2. While this factor cancels out when calculating a chi-squared statistics within strata, it has the effect of introducing stratum weights that penalize the information from the largest strata.

We are confident that Stata’s strata() option in the sts test command produces the desired overall test statistic for testing the difference in treatment groups in the presence of stratification.

References

Collett, D. 1994.
Modeling Survival Data in Medical Research. London: Chapman & Hall.
Kalbfleisch, J. D. and R. L. Prentice. 2002.
The Statistical Analysis of Failure Time Data. 2nd Edition. New York: Wiley.
Klein, J. P. and M. L. L. Moeschberger. 1997.
Survival Analysis. Techniques for Censored and Truncated Data. New York: Springer.
SAS Institute Inc. 1989.
SAS/STAT User’s Guide, Version 6, Volume 2. 4th ed. Cary, NC: SAS Institute.
The Stata Blog: Not Elsewhere Classified Find us on Facebook Follow us on Twitter LinkedIn Google+ Watch us on YouTube