... | @@ -73,13 +73,13 @@ Next we merge the reads to produce single sequence variants: |
... | @@ -73,13 +73,13 @@ Next we merge the reads to produce single sequence variants: |
|
|
|
|
|
> mergers <- mergePairs(dadaFs, filtFs, dadaRs, filtRs, verbose=TRUE)
|
|
> mergers <- mergePairs(dadaFs, filtFs, dadaRs, filtRs, verbose=TRUE)
|
|
|
|
|
|
Be careful inspecting this object, it's really big! Have a go at subsetting it using square brackets and use head() to inspect a part of the output.
|
|
Be careful inspecting this object, it's really big! Have a go at subsetting it using square brackets and use head() to inspect just a part of the output.
|
|
|
|
|
|
The next step is to make the ASV table (formatted like a standard OTU table with sequences and counts across samples):
|
|
The next step is to make the ASV table (formatted like a standard OTU table with sequences and counts across samples):
|
|
|
|
|
|
> seqtab <- makeSequenceTable(mergers)
|
|
> seqtab <- makeSequenceTable(mergers)
|
|
|
|
|
|
Then we need to remove chimeric sequences. DADA2 does this by aligning sequences to each other and checking if any have two "parents" - i.e. they match perfectly over some of their length to two different unique sequences.
|
|
Then we need to remove chimeric sequences. DADA2 does this by aligning sequences to each other and checking if any have two "parents" - i.e. they match perfectly over some of their length to two different unique sequences. This step takes an age, so if it's late in the day just leave it running and we can return to it tomorrow.
|
|
|
|
|
|
> seqtab.nochim <- removeBimeraDenovo(seqtab, method="consensus", multithread=8, verbose=TRUE)
|
|
> seqtab.nochim <- removeBimeraDenovo(seqtab, method="consensus", multithread=8, verbose=TRUE)
|
|
|
|
|
... | @@ -93,7 +93,8 @@ Don't worry about what this code is actually doing (clearly the developer knew w |
... | @@ -93,7 +93,8 @@ Don't worry about what this code is actually doing (clearly the developer knew w |
|
> rownames(track) <- sample.names
|
|
> rownames(track) <- sample.names
|
|
> head(track)
|
|
> head(track)
|
|
|
|
|
|
The final step is then to classify the sequences, using the silva database we downloaded before.
|
|
How do the numbers look? It's ok to lose a few more to chimeras than the other steps, so long as there's not anywhere else that we lose more than 50% of our reads.
|
|
|
|
The final step is then to classify the sequences, using the silva database we downloaded before. This also takes a while, so again, we might start it and then let it run while we do something else.
|
|
|
|
|
|
> taxa <- assignTaxonomy(seqtab.nochim, "silva_nr_v132_train_set.fa.gz", multithread=8)
|
|
> taxa <- assignTaxonomy(seqtab.nochim, "silva_nr_v132_train_set.fa.gz", multithread=8)
|
|
|
|
|
... | | ... | |