<p>If your null hypothesis is that the event of being selected by your program is statistically independent from the event of being non-coding, you can use <a href="https://en.wikipedia.org/wiki/Fisher%27s_exact_test" rel="nofollow noreferrer">Fisher's exact test</a> to reject the null hypothesis.</p>

<p>Suppose you have 10,000 (I'm making that up but presumably you have the real number) total sequences, 8,700 coding and 1,300 noncoding. Given that your program selects 728 sequences, you would expect around 13%*728 ~= 95 noncoding sequences. But your program picked 131 noncoding sequences and you want to know if 131 is far enough away from 95 to conclude that the difference is not just due to chance.</p>

<p>The probability of selecting 131 or more of these noncoding sequences (given that we are choosing 728 total sequences) is given by the so-called <a href="https://en.wikipedia.org/wiki/Hypergeometric_distribution" rel="nofollow noreferrer">hypergeometric distribution</a>. The two-sided test will also incorporate the probably of getting an equally unlikely low number (in this case, somewhere around 57) or less.</p>

<p>In R for example, you can calculate this using the 2-by-2 contingency table as follows:</p>

<pre><code>fisher.test(rbind(c(597, 131), c(8103, 1169)))
</code></pre>

<p>which gives:</p>

<pre><code>p-value = 7.399e-05
</code></pre>

<p>So in this situation (and assuming the probability model described), there's around a 0.0074% chance of observing a ratio as "extreme" as the one you observed, which indicates it's unlikely that the two events are independent.</p>

<p>Note: As others have mentioned, a binomial test could be appropriate too. If your number of total sequences is very large, then the two results will converge. But if it isn't, then Fisher's exact test will be more appropriate (essentially since you're sampling without replacement).</p>

CC BY-SA 4.0
2019-01-12T01:10:09.217Z
2019-01-12T23:23:43.450Z
2019-01-12T23:23:43.450Z
leekaiinthesky
<p>A Poisson or binomial would suffice. </p>

<p>I'm not at all sure the difference will significantly deviate from randomness. </p>

<p>Perl: 95% confidence interval of a non-parametric bootstrap for e.g. 1000 replicates is an alternative approach. I'd be more comfortable coding that solution in Perl than I would in Python, so is unlikely to be of much use here. Anyway, I would use a 2-dimensional array for your total data-set, the second dimension being coding/non-coding, I would then use the inbuilt random-number generator and multiply by the size of the array to randomly pluck out each element. 95% CI of the mean works by ranking the data and removing the top and bottom 2.5% and taking the interval for a single replicate (which is repeated 999 times more).</p>

<p>There will be a deep learning approach to this data, which would be easier to implement in Python, but I'd need to think about this because it is not simple classification.</p>

<p>Again my warning is that given the sample size, the result at best might be borderline. There are likely to be lots of interesting in your data, but coding/non-coding may not be one of them.</p>

CC BY-SA 4.0
2019-01-08T17:28:53.590Z
Hi, Michael. Thank you. Yes, I am trying binomial analysis now, don't worry about whether this results in a null hypothesis. I am simply trying to rule out a possible explanation for some results I obtained in the wet-lab. Whatever result I get is fine by me I am just happy to show I have properly investigated possible explanations.
syn_bio_delta
CC BY-SA 4.0
2019-01-08T17:33:49.447Z
Good to know. I'm definitely advocating DL for molecular solutions, they are very, very cool.
M__
CC BY-SA 4.0
2019-01-08T18:05:37.413Z
Haha, I'm interested in DL however my data science, coding and statistics skills need a lot of work before I can use such a tool. But I read a lot about the applications and I'm keen to get to the point where I can use such tools.
syn_bio_delta
CC BY-SA 4.0
2019-01-08T17:18:57.327Z
2019-01-08T17:18:57.327Z
M__
<p>I have used a python script to identify target sequences in a DNA sequence file. There are two classes of sequence: coding and non-coding. </p>

<p>I have identified 728 sequences of interest. 597 of these fall into the coding regions and 131 of these fall into the non-coding regions. This is the equivalent of 18% non-coding. The total non-coding region in the sequence file is 13%. </p>

<p>Is there a statistical tool to demonstrate the python script identified target sequences in a non-random fashion?<br>
If the script identified sequences that were randomly distributed then 13% of them would have been found in the non-coding region, from a total of 728 this seems like it should be reliable.</p>

CC BY-SA 4.0
2019-01-08T06:10:55.430Z
It sounds like you're fishing for a statistic to prove a point. I wouldn't try to over-analyze your data. You might just consider stating that 18% of the sequences are in non-coding regions which is greater than the expected null hypothesis of 13% of the sequences being in non-coding regions.
conchoecia
CC BY-SA 4.0
2019-01-08T07:55:30.327Z
One could ask how big are those sequences of interest (SOI) and if there are some of them that span between coding and non-coding regions. Also the distribution of SOI could depend on TADs or other higher structures of your sequence (if it is a crhomosome size sequence). To model the distribution of SOI would require knowing a lot more about the SOI and the sequence you are working with and without a null model it cannot be tested if it is randomly following this model or not
llrs
CC BY-SA 4.0
2019-01-08T15:58:16.670Z
Hi again @conchoecia. Thanks for the feedback. With this point, however, we find ourselves in disagreement. There is a chance that if my results are due to random chance there may well be an uneven distribution. Say if my non-coding hits were 14% it would be obvious this would not be reliable evidence in support of my hypothesis. Furthermore, I have no problem with obtaining a negative result but however I interpret my results it needs to be statistically reliable, there is a chance that 18% is not significantly different to 13% and therefore cannot be used in favour of my hypothesis.
syn_bio_delta
CC BY-SA 4.0
2019-01-08T16:30:52.320Z
@llrs Hi. The sequences of interest are 29 bp (promotors) and I am searching chromosome FASTA. There are a some overlaps, but I have not counted them as non-coding for simplicity. The non-coding regions were defined by a different study (https://doi.org/10.1186/s13068-017-0742-z) Table 1.
syn_bio_delta
CC BY-SA 4.0
2019-01-08T16:59:22.477Z
Thanks everyone, question answered on cross-validated.
https://stats.stackexchange.com/questions/386173/genomics-stats-problem/386177#386177
syn_bio_delta
<phylogenetics><phylogeny><r><genomics>
Genomics Statistics Problem
# Genomics Statistics Problem

(Asked by: username_0 on 2019-01-07)

I have used a python script to identify target sequences in a DNA sequence file. There are two classes of sequence: coding and non-coding.

I have identified 728 sequences of interest. 597 of these fall into the coding regions and 131 of these fall into the non-coding regions. This is the equivalent of 18% non-coding. The total non-coding region in the sequence file is 13%.

Is there a statistical tool to demonstrate the python script identified target sequences in a non-random fashion?  
If the script identified sequences that were randomly distributed then 13% of them would have been found in the non-coding region, from a total of 728 this seems like it should be reliable.

**Tags:** phylogenetics, phylogeny, r, genomics
> **Comments:**
> **username_1**: It sounds like you're fishing for a statistic to prove a point. I wouldn't try to over-analyze your data. You might just consider stating that 18% of the sequences are in non-coding regions which is greater than the expected null hypothesis of 13% of the sequences being in non-coding regions.  
> **username_2**: One could ask how big are those sequences of interest (SOI) and if there are some of them that span between coding and non-coding regions. Also the distribution of SOI could depend on TADs or other higher structures of your sequence (if it is a crhomosome size sequence). To model the distribution of SOI would require knowing a lot more about the SOI and the sequence you are working with and without a null model it cannot be tested if it is randomly following this model or not  
> **username_0**: Hi again @conchoecia. Thanks for the feedback. With this point, however, we find ourselves in disagreement. There is a chance that if my results are due to random chance there may well be an uneven distribution. Say if my non-coding hits were 14% it would be obvious this would not be reliable evidence in support of my hypothesis. Furthermore, I have no problem with obtaining a negative result but however I interpret my results it needs to be statistically reliable, there is a chance that 18% is not significantly different to 13% and therefore cannot be used in favour of my hypothesis.  
> **username_0**: @llrs Hi. The sequences of interest are 29 bp (promotors) and I am searching chromosome FASTA. There are a some overlaps, but I have not counted them as non-coding for simplicity. The non-coding regions were defined by a different study (https://doi.org/10.1186/s13068-017-0742-z) Table 1.  
> **username_0**: Thanks everyone, question answered on cross-validated. https://stats.stackexchange.com/questions/386173/genomics-stats-problem/386177#386177  

---

## 2 answers

### Answer by username_3 on 2019-01-12. Score: 3 | Accepted answer

If your null hypothesis is that the event of being selected by your program is statistically independent from the event of being non-coding, you can use Fisher's exact test to reject the null hypothesis.
Suppose you have 10,000 (I'm making that up but presumably you have the real number) total sequences, 8,700 coding and 1,300 noncoding. Given that your program selects 728 sequences, you would expect around 13%*728 ~= 95 noncoding sequences. But your program picked 131 noncoding sequences and you want to know if 131 is far enough away from 95 to conclude that the difference is not just due to chance.
The probability of selecting 131 or more of these noncoding sequences (given that we are choosing 728 total sequences) is given by the so-called hypergeometric distribution. The two-sided test will also incorporate the probably of getting an equally unlikely low number (in this case, somewhere around 57) or less.
In R for example, you can calculate this using the 2-by-2 contingency table as follows:
fisher.test(rbind(c(597, 131), c(8103, 1169)))

which gives:
p-value = 7.399e-05

So in this situation (and assuming the probability model described), there's around a 0.0074% chance of observing a ratio as "extreme" as the one you observed, which indicates it's unlikely that the two events are independent.
Note: As others have mentioned, a binomial test could be appropriate too. If your number of total sequences is very large, then the two results will converge. But if it isn't, then Fisher's exact test will be more appropriate (essentially since you're sampling without replacement).


---

### Answer by username_4 on 2019-01-08. Score: 2

A Poisson or binomial would suffice.

I'm not at all sure the difference will significantly deviate from randomness.

Perl: 95% confidence interval of a non-parametric bootstrap for e.g. 1000 replicates is an alternative approach. I'd be more comfortable coding that solution in Perl than I would in Python, so is unlikely to be of much use here. Anyway, I would use a 2-dimensional array for your total data-set, the second dimension being coding/non-coding, I would then use the inbuilt random-number generator and multiply by the size of the array to randomly pluck out each element. 95% CI of the mean works by ranking the data and removing the top and bottom 2.5% and taking the interval for a single replicate (which is repeated 999 times more).

There will be a deep learning approach to this data, which would be easier to implement in Python, but I'd need to think about this because it is not simple classification.

Again my warning is that given the sample size, the result at best might be borderline. There are likely to be lots of interesting in your data, but coding/non-coding may not be one of them.

> **Comments:**
> **username_0**: Hi, Michael. Thank you. Yes, I am trying binomial analysis now, don't worry about whether this results in a null hypothesis. I am simply trying to rule out a possible explanation for some results I obtained in the wet-lab. Whatever result I get is fine by me I am just happy to show I have properly investigated possible explanations.  
> **username_4**: Good to know. I'm definitely advocating DL for molecular solutions, they are very, very cool.  
> **username_0**: Haha, I'm interested in DL however my data science, coding and statistics skills need a lot of work before I can use such a tool. But I read a lot about the applications and I'm keen to get to the point where I can use such tools.  

---

