File: epi.ssequc.R

package info (click to toggle)
r-cran-epir 2.0.80%2Bdfsg-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 2,332 kB
  • sloc: makefile: 5
file content (109 lines) | stat: -rw-r--r-- 4,134 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
epi.ssequc <- function(treat, control, sigma, delta, n, power, r = 1, type = "equivalence", nfractional = FALSE, alpha){

  if(type == "equality"){
    
    # Sample size:
    if (!is.na(treat) & !is.na(control) & !is.na(power) & is.na(n)){
      z.alpha <- qnorm(1 - alpha / 2, mean = 0, sd = 1)
      beta <- (1 - power)
      z.beta <- qnorm(1 - beta, mean = 0, sd = 1)
      
      # http://powerandsamplesize.com/Calculators/Compare-2-Means/2-Sample-Equality:
      n <- ((sigma * (z.alpha + z.beta) / (abs(treat - control)))^2)
      
      if(nfractional == TRUE){
        n.control <- (1  + 1 / r) * n
        n.treat <- n.control * r
        n.total <- n.treat + n.control
      }
      
      if(nfractional == FALSE){
        n.control <- (1  + 1 / r) * (ceiling(n))
        n.treat <- n.control * r
        n.total <- n.treat + n.control
      }

    }
    
    # Power:
    if (!is.na(treat) & !is.na(control) & !is.na(n) & is.na(power) & !is.na(r) & !is.na(alpha)) {
      z.alpha <- qnorm(1 - alpha / 2, mean = 0, sd = 1)
      beta <- (1 - power)
      z.beta <- qnorm(1 - beta, mean = 0, sd = 1)
      
      z <- (treat - control) / (sigma * sqrt((1 / n.treat) + (1 / n.control))) 
      power <- pnorm(z - z.alpha) + pnorm(-z - z.alpha)  
    }

    rval <- list(n.total = n.total, n.treat = n.treat, n.control = n.control, power = power)
        
  }
  
  if(type == "equivalence"){
    
    # Stop if a negative value for delta entered:
    if (delta <= 0){
      stop("For an equivalence trial delta must be greater than zero.")
    }
    
    # Sample size:
    if (!is.na(treat) & !is.na(control) & !is.na(power) & is.na(n)) {
      z.alpha <- qnorm(1 - alpha, mean = 0, sd = 1)
      beta <- (1 - power)
      z.beta <- qnorm(1 - beta / 2, mean = 0, sd = 1)

      # http://powerandsamplesize.com/Calculators/Compare-2-Means/2-Sample-Equivalence:
      n <- ((sigma * (z.alpha + z.beta) / (abs(treat - control) - delta))^2)
      
      if(nfractional == TRUE){
        n.control <- (1  + 1 / r) * n
        n.treat <- n.control * r
        n.total <- n.treat + n.control
      }
      
      if(nfractional == FALSE){
        n.control <- (1  + 1 / r) * (ceiling(n))
        n.treat <- n.control * r
        n.total <- n.treat + n.control
      }
    }
    
    # Power:
    if (!is.na(treat) & !is.na(control) & !is.na(n) & is.na(power) & !is.na(r) & !is.na(alpha)) {
      
      z.alpha <- qnorm(1 - alpha, mean = 0, sd = 1)
      beta <- (1 - power)
      z.beta <- qnorm(1 - beta / 2, mean = 0, sd = 1)
      
      # Work out the number of subjects in the control group. r equals the number in the treatment group divided by the number in the control group.
      
      if(nfractional == TRUE){
        n.control <- 1 / (r + 1) * n
        n.treat <- n - n.control
        n.total <- n.treat + n.control
      }
      
      if(nfractional == FALSE){
        n.control <- ceiling(1 / (r + 1) * n)
        n.treat <- n - n.control
        n.total <- n.treat + n.control
      }
      
      z <- (abs(treat - control) - delta) / (sigma * sqrt((1 / n.treat) + (1 / n.control))) 
      power <- 2 * (pnorm(z - z.alpha, mean = 0, sd = 1) + pnorm(-z - z.alpha, mean = 0, sd = 1)) - 1    
    }
    rval <- list(n.total = n.total, n.treat = n.treat, n.control = n.control, delta = delta, power = power)
  }
  
  return(rval)
}
    
# Chow S, Shao J, Wang H. 2008. Sample Size Calculations in Clinical Research. 2nd Ed. Chapman & Hall/CRC Biostatistics Series. page 62

# epi.ssequc(treat = 5, control = 4, sigma = 10, delta = 5, n = NA, power = 0.80, r = 1, nfractional = FALSE, alpha = 0.05)
# n.treat = 108, n.control = 108, n.total = 216
# Agrees with http://powerandsamplesize.com/Calculators/Compare-2-Means/2-Sample-Equivalence

# epi.ssequc(treat = 5, control = 4, sigma = 10, delta = 5, n = NA, power = 0.80, r = 2, nfractional = TRUE, alpha = 0.05)
# n.treat = 162, n.control = 81, n.total = 243
# Agrees with http://powerandsamplesize.com/Calculators/Compare-2-Means/2-Sample-Equivalence