Power / Sample size calculation for comparing two non-normally distributed continuous variables using rank-sum test

[Adam Badenoch]

Dear Statalist users,

I would like to perform a power/sample size calculation for a study comparing two non-normally distributed numerical variables (a numerical [integer] score ranging from 0-28).

I can't see an easy solution within stata's inbuilt commands including the "power" command.

Searching further I found a suggestion to do this using simulation in the Stata journal (The Stata Journal (2002) 2, Number 2, pp. 107–124). Running the code from page 122-3 in Stata unedited resulted in an error (‘medrat’ invalid name). I have tried altering any references to arguments and locals in the code from ‘medrat’ format to `mediate' as I thought they should be but this doesn't seem to help.

Does anyone have any tips on how to run the program presented in the Stata journal? Once I can do this I should be able to adapt it to my situation (equal group sizes and 2-tailed test). Alternatively are there any simpler solutions to my problem?

Many thanks,

Adam

[William Lisowski]

First off, you have been done in by Statacorp's presentation of the keyboard characters used to dereference local macros

` and '

by the typographical characters

‘ and ’

which Stata does not recognize, as you saw and tried to correct. Below I have placed a copy of rksimsum.do with these corrections made.

For future reference, please review the Statalist FAQ linked to from the top of the page, as well as from the Advice on Posting link on the page you used to create your post. Note especially sections 9-12 on how to best pose your question.

Section 12.1 is particularly pertinent

12.1 What to say about your commands and your problem

Say exactly what you typed and exactly what Stata typed (or did) in response. N.B. exactly!
...
Never say just that something "doesn't work" or "didn't work", but explain precisely in what sense you didn't get what you wanted.

It would have been helpful for you to have shown us your copy of rksimsum.do, and it would have been helpful for you to have told us what precisely went wrong when you tried to run it.

This is old code (version 7) and my attempts at running it failed on various syntax issues that I was not able to overcome with the amount of effort I care to devote to this. To start with, running the command in the paper that runs the do-file

. for n in num 50(10)160: run rksumsim 1.15 0.50 n .05 1000

->  ru50 rksumsim 1.15 0.50 50 .05 1000
command ru50 is unrecognized
r(199);

shows us that the "n" in run, as well as the "n" in the arguments list, is being replaced with the value of the local macro n created by the for loop. Running it under version control did not help.

Perhaps someone with a better grasp of running old code, will be able to take this further. I'll limit my contribution to supplying the code with the correct quotation convention.

/* rksumsim.do
       Rank sum 2-sample rank test power simulation  */
      version 7.0
      args medrat sig n alpha NIT
      /*
         medrat is the assumed ration of medians (men to women)
         sig is the standard deviation of log salaries
         n is sample size for women
         alpha is the level of the test
         NIT is number of simulations
      */
      local dif=log(`medrat') /* difference in means of log salaries */
      tempname buff
      postfile `buff' pv1 zrat using temppow,replace
forv it = 1/`NIT' { drop _all
              set obs `n'
              gen byte group=1
              local N=3*`n'
              set obs `N'
              replace group=2 if group==.
      /* without loss of generality - take mean log salary for women to be zero */
              gen y= 0 + `sig'*invnorm(uniform()) if group==1
      /* mean log salary for men is `dif' */
              replace y= `dif'+`sig'*invnorm(uniform()) if group==2
              replace y=exp(y)
              ranksum y,by(group) /* test medians */
              local zrat=r(z)
scalar pv1=normprob(`zrat') post `buff' (pv1) (zrat) if `it'==`NIT' {
                save temp,replace
} }
postclose `buff'
use temppow,clear
/* 1-sided pvalue, H1: Group2 > Group1 */
/*The dataset in memory now contains NIT simulated p-values. To get an
estimate of the power, say for alpha=.05, just count the proportion of pv's
that are less than 0.05:  */
local n2 = 2*`n'
count if pv1<`alpha'
scalar power=r(N)/`NIT'
noi di "Power is ",%8.4f scalar(power),/*
*/ " for `n' women, `n2' men, ratio = `medrat', alpha = `alpha'"

[Dave Airey]

I think the harder problem is simulating the data generating process for which limited information was provided. How was it non-normal? Is your data (expected to be) distributed in the same non-normal way between your groups?

[Adam Badenoch]

Thank you for your help William and Dave,

I response to Dave:
The data is non-normal in the sense that it is skewed. The distributions are expected to be the same shape in both groups. The exact distribution is unknown as it is a novel, recently proposed clinical scoring system.

In response to William:
The code you present is one version of multiple edited versions of the published code I referenced. I did not supply specific errors for each version because I thought I should keep the original referenced code as the source of reference, for which an error was presented.

As far as explaining precisely in what sense I did not get what I wanted I can rephrase my problem as follows:

I would like to perform a power/sample size calculation for a study comparing two non-normally distributed numerical variables (a skewed numerical [integer] score ranging from 0-28). How can I do this? (The proposed simulation strategy is but one possible option).

Kind Regards,

Adam

[Dave Airey]

The difficulty I mentioned continues in your answer. You have to decide how skewed the distribution is going to be, or perhaps you will investigate multiple skewness amounts (and right or left?). The entry for user made power commands is here (v15):

help power usermethod

I believe SAS has a solution for this problem if you have access. SAS has many nonnormal distributions to choose from.

Don't forget that if resources allow, samples larger than 30 may allow parametric approaches despite the skewness.

Here's another solution in R: https://cran.r-project.org/web/packa...samplesize.pdf, and the reference:

https://www.ncbi.nlm.nih.gov/pubmed/17487941

[Adam Badenoch]

Thanks Dave,

Because I don't have a good estimate of the degree of skewness I will perform sample size/power calculations for a variety of possible distributions.

I'll try the SAS and R options.

Adam

[Joseph Coveney]

The data is non-normal in the sense that it is skewed. The distributions are expected to be the same shape in both groups.

How can that be? If the values are bounded at zero and 28, any medically interesting hypothesized difference in a measure of central tendency would force a difference in shapes, especially if the distributions are expected to be skewed. Your rank-sum test is going to be affected by features of the data that you're ostensibly uninterested in.

Why not consider sticking with parametrics, and use a Johnson skew t-test or a generalized linear model with a suitable link function?

[Adam Badenoch]

Thanks for your suggestions Joseph,

The expectation of distributions being the same in both groups is an expectation under the null hypothesis.

I'm not sure I follow the second point, I would have thought any rank differences in scores between the two groups is as clinically relevant as differences in the central tendency of any given distribution. So I don't see this as a reason to choose a method based on distribution/central tendency shifts rather than score rankings on this basis. Which 'uninteresting' features of the data are you referring to?

If the distributions can be easily transformed to allow reasonable use of parametric methods I'll do that. The original post was designed to obtain help for cases where this isn't easy to do.

I'd be interested if you have any further thoughts?

Adam

[Joseph Coveney]

The expectation of distributions being the same in both groups is an expectation under the null hypothesis.

I wasn't referring to conditions under the null hypothesis.

Which 'uninteresting' features of the data are you referring to?

Shape, as opposed to location.

If the distributions can be easily transformed to allow reasonable use of parametric methods I'll do that. The original post was designed to obtain help for cases where this isn't easy to do.

The methods in the suggestion obviate the need for transformation, easy or otherwise.

[Adam Badenoch]

Thanks Joseph,

Sorry I must have misunderstood your "how can that be?" Did you mean how can the data be skewed? This is simply a guess based on distributions of similar clinical scoring systems.

I think useful information can be gleamed from ranking, central tendency and shape in this instance. The three seem co-dependent in my mind due to the score boundary constraints. Although I agree shape, as a summary statistic, in isolation isn't relevant.

Regarding the methods, sorry, I agree both of your suggestions are sensible solutions. Consider the idea of a simple transformation followed by linear regression a third possibility not an exclusive one. Certainly the GLM option allows for more flexibility in terms of distribution shape.

Returning to the initial question about power/sample size calculation, I guess irrespective of the method used, simulations similar to the earlier ranksumsim code but using the test of choice (Johnson skew, GLM, etc) would appear to be the way to go. Are there any other methods I should be aware of?

[Joseph Coveney]

I'm not aware of any other methods, but It's been years since I've used anything but simulation for power and sample size estimation, and so I'm not up on what's become available in.the interim.

[Adam Badenoch]

Good to know. Thanks Joseph.