Sample Size

For the calculations I borrowed the idea from a paper by Zhu and Lakkis (2014) titled Sample size calculation for comparing two negative binomial rates. Which has also been incoorporated as a Package (MKmisc) in R.

1) power.nb.test function from MKmisc Package

Assumptions in the calculation of the sample sizes are as follows;-

I looked at two scenerios:-
- Fixed Power at \(80\%\) but varied the other parameters to calculate the sample sizes.
- Fixed the N’s, and calculated the power.

\(mu0\) - The expected rate of colds per time unit for Control group 0.

\(mu1\) - The expected rate of cold per time unit for Treatment group 1.

Effect size in terms of Rate Ratio(RR) - i.e, Reduction in the rates .
- RR = \(mu1/mu0\).
- I considered 0.33, 0.50, and 0.67.

Duration - This is the average Treatment duration in Years. 2 durations have been considered.
- 0.25 years - This is assuming no drop-outs so everyone was observed for the whole study period of 90 days which \(\approx 0.25yrs\).
- 0.18 years - This is assuming the drop-outs. That is, only 70% complete the 90 days.

\(\theta\) - Parameter of the Negative Binomial distribution - The dispersion parameter

# Libraries used

lib <- c("epiDisplay","tidyverse","dplyr","tidyr","plyr","kableExtra")
lapply(lib, library, character.only = TRUE)

# install.packages("MKmisc")

library(MKmisc)
# _______________________
# Arguments used by the Function
# power.nb.test(n = NULL, 
#               mu0, #expected rate of events per time unit for group 0/control
#               mu1, #expected rate of events per time unit for group 1
#               RR,  #ratio of expected event rates: mu1/mu0
#               duration = , - #(average) treatment duration
#               theta =, -  parameter of negative binomial distribution; 
#               ssize.ratio = Allocation ratio,
#               sig.level = 0.05,
#               power = NULL, Power of 0.8 = if we performed a study 1000 times,
#                       we would see a statistically significant difference 80% 
#                       of the time
#               alternative = c("two.sided", "one.sided"),
#               approach = )
# _______________________

# example of how to use the function
n <-power.nb.test(mu0=0.7, RR=0.67,theta = 1/0.1,power=0.8,
                 duration = 0.18, ssize.ratio = 3,approach = 2)


# _______________________ ____________
# Giving values to the Parameters

# mu1  <- c(0.469,0.201,0.268,0.335,0.402,0.4355)
mu0         <- c(0.7,0.3,0.5,0.65)
RR          <- c(0.33,0.50,0.67) #  mu1/mu0 r1/r0
power       <- 0.80
ssize.ratio <- 3 # 1:1:1:1 Allocation ratio - 1 control:3 trt
theta       <- c(0.1,0.5) # k dispersion perameter
approach    <- 2
duration    <- c(0.25,0.18) # average exposure time yrs

# _______________________ ___________________________________________
#  Get N given power and other parameters

  L1 <- length(mu0)
  L2 <- length(RR)
  L3 <- length(theta)
  L4 <- length(duration)

  
# Data frame to collect results
  par.vec <- data.frame(Mu0      = rep(mu0,times=L2),
                        RR       = rep(RR,each=L1),
                        Theta    = rep(round(1/theta,2),each=(L1*L2)),
                        Duration = rep(duration, each=L1*L2*L3)) %>%
               dplyr::distinct(Theta,RR,Mu0,Duration, .keep_all = TRUE)%>% #drop duplicates
               dplyr::arrange(desc(Duration),RR)

    
  par.vec <- par.vec%>%
                mutate(Power       = power,
                       Row         = 1:(L1*L2*L3*L4),
                       Ssize.ratio = ssize.ratio,
                       Approach    = approach
                      ) %>% 
                select(Row,Duration,RR,Theta,Mu0,everything())
  
  
s.size <- function(x){ 
# Sample size using power.nb.test  
k <- power.nb.test(mu0 = x$Mu0, RR = x$RR, theta = x$Theta, 
                   duration = x$Duration, ssize.ratio = x$Ssize.ratio,
                   power = x$Power, approach = x$Approach,
                   alternative = "two.sided")
  # save the N's in a df
k1 <- data.frame(n=round(k$n,0), n1 = round(k$n1,0))
}


# Call the function
n.n <- par.vec %>% group_by(Row) %>% summarise(Row = Row,
                                                 n = s.size(.)$n,
                                                n1 = s.size(.)$n1)

# Merge in the N's infor
par.vec1 <- par.vec %>% full_join(n.n,by="Row")

# _______________________ ___________________________________________




# _______________________ ___________________________________________

#  Get POWER given N and other parameters

n          <- c(55,70,100,140,175)
sig.level  <- c(0.05,0.1)

L5 <- length(n)
L6 <- length(sig.level)

# Data frame to collect results
 par.vec2 <- data.frame(Mu0      = rep(mu0,times=L5),
                        n        = rep(n, each=L1),
                        RR       = rep(RR,each=L1*L5),
                        Sig.level= rep(sig.level,each=L1*L2*L5),
                        Theta    = rep(round(1/theta,2),each=(L1*L2*L5)),
                        Duration = rep(duration, times=L1*L2*L5)) %>%
               dplyr::distinct(.keep_all = TRUE)%>% #drop duplicates
               dplyr::arrange(n,desc(Duration),Sig.level,RR,Mu0)

    
  par.vec2 <- par.vec2%>%
                mutate(Row      = 1:(L1*L2*L5*L6),
                       Ssize.ratio = ssize.ratio,
                       Approach    = approach
                      ) %>% 
                  select(Row,Duration,n,Sig.level,RR,Mu0,Theta,everything())


s.size2 <- function(x){ 
# Power using power.nb.test  
  k1 <- power.nb.test(mu0 = x$Mu0, RR = x$RR, theta = x$Theta, 
                      duration = x$Duration, ssize.ratio = x$Ssize.ratio,
                      n = x$n, approach = x$Approach,
                      alternative = "two.sided")
  # save the Power in a df
  k2 <- data.frame(Power=round(k1$power,2)*100)
}

# Call the function
n.n2 <- par.vec2 %>% group_by(Row) %>% summarise(Row   = Row,
                                                Power = s.size2(.)$Power)
# Merge in the N's infor
par.vec3 <- par.vec2 %>% full_join(n.n2,by="Row")

N’s given the other parameters

N’s given the other parameters
Duration	RR	Theta	Mu0	Power	Ssize.ratio	n	n1
No Drop-Out
0.25	0.33	10	0.70	0.8	3	74	223
0.25	0.33	2	0.70	0.8	3	78	233
0.25	0.33	10	0.65	0.8	3	80	240
0.25	0.33	2	0.65	0.8	3	83	250
0.25	0.33	10	0.50	0.8	3	104	311
0.25	0.33	2	0.50	0.8	3	107	321
0.25	0.33	10	0.30	0.8	3	172	516
0.25	0.33	2	0.30	0.8	3	175	526
0.25	0.50	10	0.70	0.8	3	158	473
0.25	0.50	2	0.70	0.8	3	166	499
0.25	0.50	10	0.65	0.8	3	170	509
0.25	0.50	2	0.65	0.8	3	178	535
0.25	0.50	10	0.50	0.8	3	220	660
0.25	0.50	2	0.50	0.8	3	229	686
0.25	0.50	10	0.30	0.8	3	365	1096
0.25	0.50	2	0.30	0.8	3	374	1122
0.25	0.67	10	0.70	0.8	3	425	1276
0.25	0.67	2	0.70	0.8	3	451	1354
0.25	0.67	10	0.65	0.8	3	458	1373
0.25	0.67	2	0.65	0.8	3	484	1451
0.25	0.67	10	0.50	0.8	3	593	1778
0.25	0.67	2	0.50	0.8	3	619	1857
0.25	0.67	10	0.30	0.8	3	984	2951
0.25	0.67	2	0.30	0.8	3	1010	3029
Considering Drop-Outs
0.18	0.33	10	0.70	0.8	3	103	308
0.18	0.33	2	0.70	0.8	3	106	318
0.18	0.33	10	0.65	0.8	3	111	332
0.18	0.33	2	0.65	0.8	3	114	342
0.18	0.33	10	0.50	0.8	3	143	430
0.18	0.33	2	0.50	0.8	3	147	441
0.18	0.33	10	0.30	0.8	3	239	716
0.18	0.33	2	0.30	0.8	3	242	726
0.18	0.50	10	0.70	0.8	3	218	655
0.18	0.50	2	0.70	0.8	3	227	681
0.18	0.50	10	0.65	0.8	3	235	705
0.18	0.50	2	0.65	0.8	3	244	731
0.18	0.50	10	0.50	0.8	3	305	914
0.18	0.50	2	0.50	0.8	3	313	940
0.18	0.50	10	0.30	0.8	3	506	1519
0.18	0.50	2	0.30	0.8	3	515	1545
0.18	0.67	10	0.70	0.8	3	588	1764
0.18	0.67	2	0.70	0.8	3	614	1843
0.18	0.67	10	0.65	0.8	3	633	1899
0.18	0.67	2	0.65	0.8	3	659	1977
0.18	0.67	10	0.50	0.8	3	821	2462
0.18	0.67	2	0.50	0.8	3	847	2541
0.18	0.67	10	0.30	0.8	3	1364	4091
0.18	0.67	2	0.30	0.8	3	1390	4169

NOTE: n = sample size of control group, n1 = sample size of treatment group

Calculating Power given N’s and other parameters

I have dropped cases where the resulting power was too low, i.e, <50.

Power given N’s and other parameters
Duration	n	Sig.level	RR	Mu0	Theta	Ssize.ratio	Power
No Drop-Out
0.25	55	0.05	0.33	0.50	10	3	53
0.25	55	0.05	0.33	0.70	10	3	67
0.25	55	0.10	0.33	0.50	2	3	52
0.25	55	0.10	0.33	0.70	2	3	65
0.25	70	0.05	0.33	0.50	10	3	63
0.25	70	0.05	0.33	0.70	10	3	78
0.25	70	0.10	0.33	0.50	2	3	62
0.25	70	0.10	0.33	0.70	2	3	76
0.25	100	0.05	0.33	0.50	10	3	79
0.25	100	0.05	0.33	0.70	10	3	90
0.25	100	0.05	0.50	0.70	10	3	61
0.25	100	0.10	0.33	0.50	2	3	77
0.25	100	0.10	0.33	0.70	2	3	89
0.25	100	0.10	0.50	0.70	2	3	58
0.25	140	0.05	0.33	0.50	10	3	90
0.25	140	0.05	0.33	0.70	10	3	97
0.25	140	0.05	0.50	0.50	10	3	61
0.25	140	0.05	0.50	0.70	10	3	75
0.25	140	0.10	0.33	0.50	2	3	89
0.25	140	0.10	0.33	0.70	2	3	96
0.25	140	0.10	0.50	0.50	2	3	59
0.25	140	0.10	0.50	0.70	2	3	73
0.25	175	0.05	0.33	0.50	10	3	95
0.25	175	0.05	0.33	0.70	10	3	99
0.25	175	0.05	0.50	0.50	10	3	70
0.25	175	0.05	0.50	0.70	10	3	84
0.25	175	0.10	0.33	0.50	2	3	95
0.25	175	0.10	0.33	0.70	2	3	99
0.25	175	0.10	0.50	0.50	2	3	69
0.25	175	0.10	0.50	0.70	2	3	82
Considering Drop-Outs
0.18	55	0.05	0.33	0.65	10	3	51
0.18	70	0.05	0.33	0.65	10	3	61
0.18	70	0.10	0.33	0.65	2	3	59
0.18	100	0.05	0.33	0.65	10	3	76
0.18	100	0.10	0.33	0.65	2	3	75
0.18	140	0.05	0.33	0.30	10	3	57
0.18	140	0.05	0.33	0.65	10	3	88
0.18	140	0.05	0.50	0.65	10	3	58
0.18	140	0.10	0.33	0.30	2	3	57
0.18	140	0.10	0.33	0.65	2	3	87
0.18	140	0.10	0.50	0.65	2	3	57
0.18	175	0.05	0.33	0.30	10	3	67
0.18	175	0.05	0.33	0.65	10	3	94
0.18	175	0.05	0.50	0.65	10	3	68
0.18	175	0.10	0.33	0.30	2	3	66
0.18	175	0.10	0.33	0.65	2	3	93
0.18	175	0.10	0.50	0.65	2	3	66

Sample Size

Tuesday May 29, 2018