R/allClasses.R, R/show-methods.R, R/subset-methods.R
MultiAmplicon-class.RdThe MultiAmplicon class is a container that stores at least primer
pairs, read files and progressively processed data in an 'amplicon
x samples' format. The slots in this object are incrementally
filled with by running wrappers functions (mostly around functions
from the dada2 package). The object is treated (subsetted
etc.) like a (pseudo) matrix, colums are samples, rows are
different amplicons.
Subsetting for MultiAmplicon objects should conveniently subset all (potentially) filled slots
MultiAmplicon(
PrimerPairsSet = PrimerPairsSet(),
PairedReadFileSet = PairedReadFileSet(),
.Data = matrix(ncol = 0, nrow = 0),
stratifiedFiles = list(),
sampleData = new("sample_data", data.frame(row.names = names(PairedReadFileSet),
readsF = PairedReadFileSet@readsF, readsR = PairedReadFileSet@readsF)),
derep = list(),
dada = list(),
mergers = list(),
sequenceTable = list(),
sequenceTableNoChime = list(),
taxonTable = list()
)
getPrimerPairsSet(MA)
getPairedReadFileSet(MA)
getRawCounts(MA)
getStratifiedFilesF(MA, simplify = TRUE)
getStratifiedFilesR(MA, simplify = TRUE)
getDerepF(MA, simplify = TRUE)
getDerepR(MA, simplify = TRUE)
getDadaF(MA, simplify = TRUE)
getDadaR(MA, simplify = TRUE)
getMergers(MA, simplify = TRUE)
getSequenceTable(MA, simplify = TRUE)
getSequenceTableNoChime(MA, simplify = TRUE)
getTaxonTable(MA, simplify = TRUE)
getSequencesFromTable(MA)
# S4 method for MultiAmplicon
getSequencesFromTable(MA)
# S4 method for MultiAmplicon
show(object)
# S4 method for MultiAmplicon,index,index,ANY
[(x, i, j, ..., drop = FALSE)
# S4 method for MultiAmplicon,index,missing,ANY
[(x, i, j, ..., drop = FALSE)
# S4 method for MultiAmplicon,missing,index,ANY
[(x, i, j, ..., drop = FALSE)
| PrimerPairsSet | a set of primer pairs specifiying your
amplicons see |
|---|---|
| PairedReadFileSet | a set of paired end sequencing data files
|
| .Data | Users should not supply this parameter, the slot
is created by |
| stratifiedFiles | Users should not supply this parameter, the
slot is created by |
| sampleData | Users should not supply this parameter. It's
filled with a sample_data object from
|
| derep | Users should not supply this parameter, the slot is
created by |
| dada | Users should not supply this parameter, the slot is
created by |
| mergers | Users should not supply this parameter, the slot is
created by |
| sequenceTable | Users should not supply this parameter, the
slot is created by |
| sequenceTableNoChime | Users should not supply this parameter,
the slot is created by |
| taxonTable | Users should not supply this parameter, the slot
is created by |
| MA | MultiAmplicon-class object |
| simplify | Should a list of objects be simplified to only one object if it has length one? |
| object | A |
| x | MultiAmplicon-class object |
| i | numeric, logical or names vector for subsetting rows (== amplicons) |
| j | numeric, logical or names vector for subsetting columns (== read files, corresponding usually to samples) |
| ... | not used |
| drop | should not be used |
MultiAmplicon: Constructor for
MultiAmplicon-class
PrimerPairsSetThe primer pairs used in your experiment to
specify amplicons stored in a
PrimerPairsSet-class object.
PairedReadFileSetThe (quality filtered) fastq files (one file pair for each sample) that store your sequencing data.
.DataA numeric matrix of sequencing read counts per
amplicon and sample. Created by the function
sortAmplicons in the MultiAmplicon pipeline.
sampleDataA sample_data object from
phyloseq. The slot is created from
sample names (names of the PrimerPairsSet, which
have tto be the same as colnames(MA)). More data can be
added by addSampleData.
stratifiedFilestemporary files as a result of stratifying
into amplicons and samples using the MultiAmplicon pipeline
function sortAmplicons. Forward (sometimes
called R1) and reverse (sometimes called R2) files are stored
as a (amplicons x samples) matrix of
PairedReadFileSet-class objects.
derepA list of PairedDerep-class objects
containing pairs of derep-class objects created by
dada2’s derepFastq function or
withing the MultiAmplicon pipeline by
derepMulti.
dadaA list of PairedDada-class object
containing pairs of dada-class objects created by
dada2’s dada function. Within the
MultiAmplicon pipeline this slot is filled by
dadaMulti.
mergersA list of objects containing merged pairs of forward
and reverse reads as created by by dada2’s
mergePairs function. Within the
MultiAmplicon pipeline this slot is filled by
mergeMulti.
sequenceTableA list of matrix objects created by
dada2’s makeSequenceTable.
Samples (in rows) and amplified sequence variants (ASVs) in
columns. Within the MultiAmplicon pipeline this slot is
filled by makeSequenceTableMulti.
sequenceTableNoChimeA list of matrix objects created by
dada2’s removeBimeraDenovo.
Samples (in rows) and ASVs screened for PCR chimeras in
columns. Within the MultiAmplicon pipeline this slot is filled
by removeChimeraMulti.
taxonTableA list of matrix objects created by a function
for taxonomical annotation (for example
blastTaxAnnot. ASVs are in rows and taxnomical
ranks are in columns. MultiAmplicon(PrimerPairsSet, PairedReadFileSet)
Emanuel Heitlinger
primerF <- c("AGAGTTTGATCCTGGCTCAG", "ACTCCTACGGGAGGCAGC", "GAATTGACGGAAGGGCACC", "YGGTGRTGCATGGCCGYT") primerR <- c("CTGCWGCCNCCCGTAGG", "GACTACHVGGGTATCTAATCC", "AAGGGCATCACAGACCTGTTAT", "TCCTTCTGCAGGTTCACCTAC") PPS <- PrimerPairsSet(primerF, primerR) fastq.dir <- system.file("extdata", "fastq", package = "MultiAmplicon") fastq.files <- list.files(fastq.dir, full.names=TRUE) Ffastq.file <- fastq.files[grepl("F_filt", fastq.files)] Rfastq.file <- fastq.files[grepl("R_filt", fastq.files)] PRF <- PairedReadFileSet(Ffastq.file, Rfastq.file) MA <- MultiAmplicon(PPS, PRF) ## sort into amplicons MA1 <- sortAmplicons(MA, filedir=tempfile(pattern = "dir"))#>#> #> #> #>#> #> #> #>#> #> #> #>#> #> #> #>#> #> #> #>#> #> #> #>#> #> #> #>#> #> #> #>## Only after sorting the MultiAmplicon object is really poplated ## with sensible data, now matrix-like access to different ## amplicons (primer pairs) and different sequencing read files ## (usually samples) is implemented ## the number of amplicons (primer pairs) nrow(MA)#> [1] 0#> [1] 0## dereplication is currently not supported ## MA2 <- derepMulti(MA1) ### use dada directly after sorting MA3 <- dadaMulti(MA1, selfConsist = TRUE)#> #> #>#>#> Initializing error rates to maximum possible estimate. #> selfConsist step 1 ...... #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> selfConsist step 6 #> Convergence after 6 rounds. #> Initializing error rates to maximum possible estimate. #> selfConsist step 1 ...... #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> Convergence after 5 rounds.#> #> #>#>#> Initializing error rates to maximum possible estimate. #> selfConsist step 1 ..... #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> selfConsist step 6 #> Convergence after 6 rounds. #> Initializing error rates to maximum possible estimate. #> selfConsist step 1 ..... #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> Convergence after 5 rounds.#> #> #>#>#> Initializing error rates to maximum possible estimate. #> selfConsist step 1 .. #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> Convergence after 5 rounds. #> Initializing error rates to maximum possible estimate. #> selfConsist step 1 .. #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> Convergence after 5 rounds.#> #> #>#>#> Initializing error rates to maximum possible estimate. #> selfConsist step 1 ..... #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> selfConsist step 5 #> Convergence after 5 rounds. #> Initializing error rates to maximum possible estimate. #> selfConsist step 1 ..... #> selfConsist step 2 #> selfConsist step 3 #> selfConsist step 4 #> Convergence after 4 rounds.#> #>#>#> #>#>#> #>#>#> #>#>#>#>#>#>#>#>#>#>