popower {Hmisc}  R Documentation 
Power and Sample Size for Ordinal Response
Description
popower
computes the power for a twotailed two sample comparison
of ordinal outcomes under the proportional odds ordinal logistic
model. The power is the same as that of the Wilcoxon test but with
ties handled properly. posamsize
computes the total sample size
needed to achieve a given power. Both functions compute the efficiency
of the design compared with a design in which the response variable
is continuous. print
methods exist for both functions. Any of the
input arguments may be vectors, in which case a vector of powers or
sample sizes is returned. These functions use the methods of
Whitehead (1993).
pomodm
is a function that assists in translating odds ratios to
differences in mean or median on the original scale.
simPOcuts
simulates simple unadjusted twogroup comparisons under
a PO model to demonstrate the natural sampling variability that causes
estimated odds ratios to vary over cutoffs of Y.
propsPO
uses ggplot2
to plot a stacked bar chart of
proportions stratified by a grouping variable (and optionally a stratification variable), with an optional
additional graph showing what the proportions would be had proportional
odds held and an odds ratio was applied to the proportions in a
reference group. If the result is passed to ggplotly
, customized
tooltip hover text will appear.
propsTrans
uses ggplot2
to plot all successive
transition proportions. formula
has the state variable on the
left hand side, the first righthand variable is time, and the second
righthand variable is a subject ID variable.\
multEventChart
uses ggplot2
to plot event charts
showing state transitions, account for absorbing states/events. It is
based on code written by Lucy D'Agostino McGowan posted at https://livefreeordichotomize.com/posts/20200521survivalmodeldetective1/.
Usage
popower(p, odds.ratio, n, n1, n2, alpha=0.05)
## S3 method for class 'popower'
print(x, ...)
posamsize(p, odds.ratio, fraction=.5, alpha=0.05, power=0.8)
## S3 method for class 'posamsize'
print(x, ...)
pomodm(x=NULL, p, odds.ratio=1)
simPOcuts(n, nsim=10, odds.ratio=1, p)
propsPO(formula, odds.ratio=NULL, ref=NULL, data=NULL, ncol=NULL, nrow=NULL )
propsTrans(formula, data=NULL, labels=NULL, arrow='\u2794',
maxsize=12, ncol=NULL, nrow=NULL)
multEventChart(formula, data=NULL, absorb=NULL, sortbylast=FALSE,
colorTitle=label(y), eventTitle='Event',
palette='OrRd',
eventSymbols=c(15, 5, 1:4, 6:10),
timeInc=min(diff(unique(x))/2))
Arguments
p 
a vector of marginal cell probabilities which must add up to one.
For 
odds.ratio 
the odds ratio to be able to detect. It doesn't
matter which group is in the numerator. For 
n 
total sample size for 
n1 
for 
n2 
for 
nsim 
number of simulated studies to create by 
alpha 
type I error 
x 
an object created by 
fraction 
for 
power 
for 
formula 
an R formula expressure for 
ref 
for 
data 
a data frame or 
labels 
for 
arrow 
character to use as the arrow symbol for transitions in

nrow , ncol 
see 
maxsize 
maximum symbol size 
... 
unused 
absorb 
character vector specifying the subset of levels of the left hand side variable that are absorbing states such as death or hospital discharge 
sortbylast 
set to 
colorTitle 
label for legend for status 
eventTitle 
label for legend for 
palette 
a single character string specifying the

eventSymbols 
vector of symbol codes. Default for first two symbols is a solid square and an open diamond. 
timeInc 
time increment for the xaxis. Default is 1/2 the shortest gap between any two distincttimes in the data. 
Value
a list containing power
, eff
(relative efficiency), and
approx.se
(approximate standard error of log odds ratio) for
popower
, or containing n
and eff
for posamsize
.
Author(s)
Frank Harrell
Department of Biostatistics
Vanderbilt University School of Medicine
fh@fharrell.com
References
Whitehead J (1993): Sample size calculations for ordered categorical data. Stat in Med 12:2257–2271.
Julious SA, Campbell MJ (1996): Letter to the Editor. Stat in Med 15: 1065–1066. Shows accuracy of formula for binary response case.
See Also
simRegOrd
, bpower
, cpower
, impactPO
Examples
# For a study of back pain (none, mild, moderate, severe) here are the
# expected proportions (averaged over 2 treatments) that will be in
# each of the 4 categories:
p < c(.1,.2,.4,.3)
popower(p, 1.2, 1000) # OR=1.2, total n=1000
posamsize(p, 1.2)
popower(p, 1.2, 3148)
# If p was the vector of probabilities for group 1, here's how to
# compute the average over the two groups:
# p2 < pomodm(p=p, odds.ratio=1.2)
# pavg < (p + p2) / 2
# Compare power to test for proportions for binary case,
# proportion of events in control group of 0.1
p < 0.1; or < 0.85; n < 4000
popower(c(1  p, p), or, n) # 0.338
bpower(p, odds.ratio=or, n=n) # 0.320
# Add more categories, starting with 0.1 in middle
p < c(.8, .1, .1)
popower(p, or, n) # 0.543
p < c(.7, .1, .1, .1)
popower(p, or, n) # 0.67
# Continuous scale with final level have prob. 0.1
p < c(rep(1 / n, 0.9 * n), 0.1)
popower(p, or, n) # 0.843
# Compute the mean and median x after shifting the probability
# distribution by an odds ratio under the proportional odds model
x < 1 : 5
p < c(.05, .2, .2, .3, .25)
# For comparison make up a sample that looks like this
X < rep(1 : 5, 20 * p)
c(mean=mean(X), median=median(X))
pomodm(x, p, odds.ratio=1) # still have to figure out the right median
pomodm(x, p, odds.ratio=0.5)
# Show variation of odds ratios over possible cutoffs of Y even when PO
# truly holds. Run 5 simulations for a total sample size of 300.
# The two groups have 150 subjects each.
s < simPOcuts(300, nsim=5, odds.ratio=2, p=p)
round(s, 2)
# An ordinal outcome with levels a, b, c, d, e is measured at 3 times
# Show the proportion of values in each outcome category stratified by
# time. Then compute what the proportions would be had the proportions
# at times 2 and 3 been the proportions at time 1 modified by two odds ratios
set.seed(1)
d < expand.grid(time=1:3, reps=1:30)
d$y < sample(letters[1:5], nrow(d), replace=TRUE)
propsPO(y ~ time, data=d, odds.ratio=function(time) c(1, 2, 4)[time])
# To show with plotly, save previous result as object p and then:
# plotly::ggplotly(p, tooltip='label')
# Add a stratification variable and don't consider an odds ratio
d < expand.grid(time=1:5, sex=c('female', 'male'), reps=1:30)
d$y < sample(letters[1:5], nrow(d), replace=TRUE)
propsPO(y ~ time + sex, data=d) # may add nrow= or ncol=
# Show all successive transition proportion matrices
d < expand.grid(id=1:30, time=1:10)
d$state < sample(LETTERS[1:4], nrow(d), replace=TRUE)
propsTrans(state ~ time + id, data=d)
pt1 < data.frame(pt=1, day=0:3,
status=c('well', 'well', 'sick', 'very sick'))
pt2 < data.frame(pt=2, day=c(1,2,4,6),
status=c('sick', 'very sick', 'coma', 'death'))
pt3 < data.frame(pt=3, day=1:5,
status=c('sick', 'very sick', 'sick', 'very sick', 'discharged'))
pt4 < data.frame(pt=4, day=c(1:4, 10),
status=c('well', 'sick', 'very sick', 'well', 'discharged'))
d < rbind(pt1, pt2, pt3, pt4)
d$status < factor(d$status, c('discharged', 'well', 'sick',
'very sick', 'coma', 'death'))
label(d$day) < 'Day'
require(ggplot2)
multEventChart(status ~ day + pt, data=d,
absorb=c('death', 'discharged'),
colorTitle='Status', sortbylast=TRUE) +
theme_classic() +
theme(legend.position='bottom')