My question arise from the field of genomics:

Given that every chromosome exists in two copies, a specific feature can be present in

  • 0 chromosomes (p=0)
  • 1 chromosome (p=0.5)
  • both chromosomes (p=1)

However, in real-world (noisy) data, there usually exists a set of DNA pieces (reads) per feature (i.e. genomic variant), e.g. given 5 reads, 2 of them show a certain feature and 3 don't.

Assuming a variant exists in only one chromosome, one can do a binomial test with p=0.5 (e.g. binom.test(2, 5, 0.5) in R), which yields p-value 1. Hence, it is assumed that the feature exists on one chromosome copy.

Now consider there exist 500 reads and 200 show the certain feature. The test binom.test(200, 500, 0.5) yields p-value 8.94e-06, of course.

However, given the choice between 0 and 0.5 (and 1) exclusively one would rather "trust" the second example to originate from a single chromosome feature, not at least due to the larger sample size.

Hence, I wonder if simple binomial testing is really the proper approach for this use case. On the one hand the test does not scale with sample size, as described above. On the other hand, there is no possibility to do the test for p=0 and p=1, as this will always result in p-value=0 for noisy data.

So, my question is, if there are any suggestions out there how binomial testing can be adapted to this special use case. Is there a possibility to test for a certain interval of p ($0\leq p \leq 25$, $0.25\leq p \leq 0.75$, ...), or should confidence intervals be involved somehow to scale with larger sample sizes?

Thanks a lot in advance for your suggestions (...hoping this was not too much biology).

EDIT example data

index   count   reference   variant
1   7   1   6
2   15  6   9
3   43  31  12
4   20  11  9
5   9   0   9

EDIT count refers to the overall count of observed reads covering a certain position in the genome. Column reference shows the counts of reads showing no variant, i.e. being identical to a reference genome. Column variant refers to the count of reads showing a variant (i.e. the "specific feature").

  • $\begingroup$ If the feature is really not present on any chromosome (i.e. p=0), do you ever get positive reads? In other words, is there some false positive detection that needs to be accounted for, or is the detection problem all one sided. If you could provide a snippet of data that demonstrate the distribution of reads in some real samples that would be a huge help. $\endgroup$ – atiretoo Oct 14 '15 at 13:57
  • $\begingroup$ Unfortunately, due to sequencing errors, DNA contamination, ... there will be positive reads even in case $p=0$. $\endgroup$ – corinna Oct 14 '15 at 14:52
  • $\begingroup$ OK, that's fine, just makes things interesting. Your example data don't really clarify things however. Can you relate the quantities in the data to the terms you used in your question? For example, if I pool all the rows with variant == 9, are those the denominator in the binomial test for variant 9? What is the relationship between count and reference, and does that tell me whether a read "shows a feature" or not? $\endgroup$ – atiretoo Oct 15 '15 at 16:34

Your Answer

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

Browse other questions tagged or ask your own question.