0
$\begingroup$

I am analyzing tumor DNA sequencing data to perform variant calling. I need to find the minimum number of DNA strands (sample size; independent samples, 'depth of coverage') needed to detect mutations occurring at frequencies of 2-5%, with 95% confidence, given a 1% background mutation rate, at power levels of 0.8, 0.9, and 0.99.

Past collaborators did this, using a "cloglog Binomial distribution"(?), and got some of the following results:

frequency: 0.02
power: 0.8
alpha: 0.05
sample size: 1239

frequency: 0.03
power: 0.8
alpha: 0.05
sample size: 423

frequency: 0.04
power: 0.90
alpha: 0.05
sample size: 299

frequency: 0.05
power: 0.99
alpha: 0.05
sample size: 315

I am trying to replicate their analysis in R, to validate and fill in more values, but it seems like I am doing something wrong because I am not getting the same values. Using the pwr library:

library("pwr")
pwr.p.test(h = 0.02, 
           sig.level = 0.05, 
           power = 0.80, 
           alternative = "greater")

output:

 proportion power calculation for binomial distribution (arcsine transformation) 

          h = 0.02
          n = 15456.39
  sig.level = 0.05
      power = 0.8
alternative = greater

Here, it is giving me an n of 15456, when the value should be 1239

As per the docs for this package, the 'effect size' is important, so I am wondering if that might be the source of the discrepancy? And I am not sure how the 'coglog Binomial distribution' plays into it, especially since pwr says it uses an 'arcsine transformation' instead.

$\endgroup$
  • 1
    $\begingroup$ cloglog = "complementary log log" so that the probability is modeled linearly on this scale $\log(-\log(p)) = \beta_0 + \beta_1 X$. It is good for discrete time survival. There's a nice discussion on these types of models in Applied Survival Analysis by Hosmer, Lemeshow, and May 2nd ed or later. $\endgroup$ – AdamO Mar 5 '18 at 16:28
1
$\begingroup$

Found the solution, using the binom R package instead. For context:

  • Null Hypothesis: No variant present, background error rate only (1%)

  • Alternative Hypothesis: A variant is present at the given frequency (Variant Allele Frequency; 2%)


code:

library("binom")  
VAF <- 0.02
background_seq_error_rate <- 0.01
alpha <- 0.05
conf_level <- 1 - alpha
power <- 0.8

cloglog.sample.size(p.alt = VAF, 
                    p = background_seq_error_rate, 
                    power = power,
                    alpha = alpha)

output:

  p.null p.alt delta alpha power    n phi
1   0.01  0.02  0.01  0.05   0.8 1239   1

Software:

  • R version 3.2.3 (2015-12-10)

  • binom_1.1-1

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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