7
$\begingroup$

I have around 1200 samples as columns and 60,000 genes with Htseq-Counts data. Before normalization with voom function I want to do filtering step.

I want to remove genes whose expression is == 0 in at least 10 samples.

Can I do this with read counts itself or should I transform counts to cpm?

If filtering with counts data is fine May I know how to do this.

Improve this question
$\endgroup$

2 Answers 2

Reset to default
8
$\begingroup$

If you have your counts in a data.frame called counts, something like this might work:

filtered.counts <- counts[rowSums(counts==0)<10, ]

For example lets assume the following data frame.

> A <- c(0,0,0,0,0)
> B <- c(0,1,0,0,0)
> C <- c(0,2,0,2,0)
> D <- c(0,5,1,1,2)
> 
> counts <- data.frame(A=A, B=B, C=C, D=D)
> counts
  A B C D
1 0 0 0 0
2 0 1 2 5
3 0 0 0 1
4 0 0 2 1
5 0 0 0 2

And now you would want to remove the genes with at least 3 zeroes or more:

> filtered.counts <- counts[rowSums(counts==0)<3, ]
> 
> filtered.counts
  A B C D
2 0 1 2 5
4 0 0 2 1

Every gene with 3 or more zeroes is now removed.

Improve this answer
$\endgroup$
10
  • $\begingroup$ In at least 10 samples. So I guess it should be this filtered.counts <- counts[rowSums(counts==0)>=10, ] $\endgroup$
    – stack_learner
    Commented Sep 1, 2017 at 12:29
  • $\begingroup$ No, you keep all rows with less than 10 zeroes with my example. $\endgroup$
    – benn
    Commented Sep 1, 2017 at 12:32
  • $\begingroup$ In your way I will be having only 1/4 th of genes remaining. I want to remove genes whose expression is == 0 in at least 10 samples. With my command now I have more than half genes remaining after filtering. $\endgroup$
    – stack_learner
    Commented Sep 1, 2017 at 12:38
  • $\begingroup$ If you use it your way, you keep all the genes with more than 10 zeroes, check it for your self if you don't believe me. $\endgroup$
    – benn
    Commented Sep 1, 2017 at 12:50
  • 1
    $\begingroup$ 10 out of 1200 samples is a very strict cut off. I'm not surprised you don't have many genes remaining. $\endgroup$
    – Ian Sudbery
    Commented Sep 1, 2017 at 13:34
3
$\begingroup$

While this answers explains how to do it I want to address when and why and which thresholds to do it.

Filtering the genes with low counts is usually done because the counts are not reliable it would be noise, specially when there are low number of samples these genes disturb power of the analysis. However with 1200 samples having 10 or more samples without any count is valuable information because you have enough samples to estimate if this is due to a latent or known variable or if it is due to noise.

The CPM are used because the size of the library per sample affects the power of the experiment. The higher the library size the easier is to find if a gene is expressed or not. CPM are calculated to take into account that effect. Filtering by counts or CPM will result in different genes filtered if you use other thresholds than 0 (for instance I have filtered those genes that the average CPM where below 1, which would be different to filtering by those that the average counts is below 1.)

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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