Dada2: Filtering on length post mergePairs

Created on 8 Nov 2017  路  3Comments  路  Source: benjjneb/dada2

Hello everyone,
I am very happy with dada2 results but I am scratching my head for a couple of hours thinking how I could filter uniques obtained with paired end reads merging (hence the resulting sequence table) based on MinLength and / or maxLength. I could not find such function in dada2 and my skills are not good enough to imagine doing that from the sequence table itself. I think that many people would need this as even if paired end reads passed through filtering and they merged fine, you can get in particular shorter (to what you expect) uniques due to the PCR system you are using and its imperfections. Of course you have some grasp on the raw reads minLength already but none that I know on the merged Uniques.
Any idea anyone?
Thanks,
Jef

Most helpful comment

Hi Jef,
Pretty much any type of filtering you could imagine is available through basic R functions. The sequence table is just a matrix, and so rows (samples) and columns (ASVs) can be selected or removed at will as for any matrix. For you particular questions...

Filtering based on sequence length:

# Assumes st is your sequence table of merged sequences
MINLEN <- 240
MAXLEN <- 260
seqlens <- nchar(getSequences(st))
st.filt <- st[,seqlens >= MINLEN & seqlens <= MAXLEN]

FIltering based on total number of reads across samples:

MINABUND <- 10
abundances <- colSums(st)
st.filt <- st[,abundances >= MINABUND]

Admittedly this is less obvious than a bunch of named specific-purpose functions for people unfamiliar with R. But R's facilities for manipulating/filtering/etc on matrices are so flexible and powerful that it doesn't make sense to try to reimplement them.

All 3 comments

just a little add, it would be useful also to be able to filter the merged uniques on their minNumber of reads. I would also get any idea on how to do that easyly :)
Thanks,
Jef

Hi Jef,
Pretty much any type of filtering you could imagine is available through basic R functions. The sequence table is just a matrix, and so rows (samples) and columns (ASVs) can be selected or removed at will as for any matrix. For you particular questions...

Filtering based on sequence length:

# Assumes st is your sequence table of merged sequences
MINLEN <- 240
MAXLEN <- 260
seqlens <- nchar(getSequences(st))
st.filt <- st[,seqlens >= MINLEN & seqlens <= MAXLEN]

FIltering based on total number of reads across samples:

MINABUND <- 10
abundances <- colSums(st)
st.filt <- st[,abundances >= MINABUND]

Admittedly this is less obvious than a bunch of named specific-purpose functions for people unfamiliar with R. But R's facilities for manipulating/filtering/etc on matrices are so flexible and powerful that it doesn't make sense to try to reimplement them.

Thanks a lot for this fast and helpful answer, I fully realise this is basic R commands, I'm on my ath to learn as quick as I can :)

Was this page helpful?
0 / 5 - 0 ratings