Tell me more ×
Cross Validated is a question and answer site for statisticians, data analysts, data miners and data visualization experts. It's 100% free, no registration required.

I am looking for a program (in R or SAS or standalone, if free or low cost) that will do power analysis for ordinal logistic regression.

share|improve this question

1 Answer

up vote 9 down vote accepted

I prefer to do power analyses beyond the basics by simulation. With precanned packages, I am never quite sure what assumptions are being made.

Simulating for power is quite straight forward (and affordable) using R.

  1. decide what you think your data should look like and how you will analyze it
  2. write a function or set of expressions that will simulate the data for a given relationship and sample size and do the analysis (a function is preferable in that you can make the sample size and parameters into arguments to make it easier to try different values). The function or code should return the p-value or other test statistic.
  3. use the replicate function to run the code from above a bunch of times (I usually start at about 100 times to get a feel for how long it takes and to get the right general area, then up it to 1,000 and sometimes 10,000 or 100,000 for the final values that I will use). The proportion of times that you rejected the null hypothesis is the power.
  4. redo the above for another set of conditions.

Here is a simple example with ordinal regression:

library(rms)

tmpfun <- function(n, beta0, beta1, beta2) {
    x <- runif(n, 0, 10)
    eta1 <- beta0 + beta1*x
    eta2 <- eta1 + beta2
    p1 <- exp(eta1)/(1+exp(eta1))
    p2 <- exp(eta2)/(1+exp(eta2))
    tmp <- runif(n)
    y <- (tmp < p1) + (tmp < p2)
    fit <- lrm(y~x)
    fit$stats[5]
}

out <- replicate(1000, tmpfun(100, -1/2, 1/4, 1/4))
mean( out < 0.05 )
share|improve this answer
Thanks @gregSnow. This will be very helpful – Peter Flom Feb 7 '12 at 23:20
2  
+1, this is a very robust, universal approach. I have often used it. I'd like to suggest another feature: You can generate data for the max $N$ you would consider, then fit the model for proportions of those data by successively fitting the 1st n of them at regular intervals up to $N$ (eg, n=100, 120, 140, 160, 180, & 200). Instead of saving a p-value from each generated dataset, you can save a row of p-values. Averaging over each column gives you a quick and dirty sense of how power is changing w/ $N$, and helps you hone in on an appropriate value quickly. – gung Feb 8 '12 at 2:57
1  
@gung: yours comment make sense, would you mind adding your codes so that less experience people in R could also benefit from it? thanks – user12349 Jul 2 '12 at 14:06

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.