traverseDown {DIAlignR} | R Documentation |
Features of the root node are propagated to all leaves node. Aligned features are set/added in the multipeptide environment.
traverseDown( tree, dataPath, fileInfo, multipeptide, prec2chromIndex, mzPntrs, precursors, adaptiveRTs, refRuns, params, applyFun = lapply )
tree |
(phylo) a phylogenetic tree. |
dataPath |
(string) path to xics and osw directory. |
fileInfo |
(data-frame) output of |
multipeptide |
(environment) contains multiple data-frames that are collection of features
associated with analytes. This is an output of |
prec2chromIndex |
(list) a list of dataframes having following columns: |
mzPntrs |
(list) a list of mzRpwiz. |
precursors |
(data-frame) atleast two columns transition_group_id and transition_ids are required. |
adaptiveRTs |
(environment) For each descendant-pair, it contains the window around the aligned retention time, within which features with m-score below aligned FDR are considered for quantification. |
refRuns |
(environment) For each descendant-pair, the reference run is indicated by 1 or 2 for all the peptides. |
params |
(list) parameters are entered as list. Output of the |
applyFun |
(function) value must be either lapply or BiocParallel::bplapply. |
analytes |
(integer) this vector contains transition_group_id from precursors. It must be of the same length as of multipeptide. |
(None)
Shubham Gupta, shubh.gupta@mail.utoronto.ca
ORCID: 0000-0003-3500-8152
License: (c) Author (2020) + GPL-3 Date: 2020-07-01
dataPath <- system.file("extdata", package = "DIAlignR") params <- paramsDIAlignR() fileInfo <- getRunNames(dataPath = dataPath) mzPntrs <- list2env(getMZMLpointers(fileInfo)) features <- list2env(getFeatures(fileInfo, maxFdrQuery = params[["maxFdrQuery"]], runType = params[["runType"]])) precursors <- getPrecursors(fileInfo, oswMerged = TRUE, runType = params[["runType"]], context = "experiment-wide", maxPeptideFdr = params[["maxPeptideFdr"]]) precursors <- dplyr::arrange(precursors, .data$peptide_id, .data$transition_group_id) peptideIDs <- unique(precursors$peptide_id) peptideScores <- getPeptideScores(fileInfo, peptideIDs, oswMerged = TRUE, params[["runType"]], params[["context"]]) peptideScores <- lapply(peptideIDs, function(pep) dplyr::filter(peptideScores, .data$peptide_id == pep)) names(peptideScores) <- as.character(peptideIDs) prec2chromIndex <- list2env(getChromatogramIndices(fileInfo, precursors, mzPntrs)) multipeptide <- getMultipeptide(precursors, features) adaptiveRTs <- new.env() refRuns <- new.env() tree <- ape::read.tree(text = "(run1:9,(run2:7,run0:2)master2:5)master1;") tree <- ape::reorder.phylo(tree, "postorder") ## Not run: ropenms <- get_ropenms(condaEnv = "envName", useConda=TRUE) multipeptide <- traverseUp(tree, dataPath, fileInfo, features, mzPntrs, prec2chromIndex, precursors, params, adaptiveRTs, refRuns, multipeptide, peptideScores, ropenms) multipeptide <- getMultipeptide(precursors, features) multipeptide <- traverseDown(tree, dataPath, fileInfo, multipeptide, prec2chromIndex, mzPntrs, precursors, adaptiveRTs, refRuns, params) # Cleanup for(run in names(mzPntrs)) DBI::dbDisconnect(mzPntrs[[run]]) file.remove(list.files(dataPath, pattern = "*_av.rds", full.names = TRUE)) file.remove(list.files(file.path(dataPath, "xics"), pattern = "^master[0-9]+\\.chrom\\.sqMass$", full.names = TRUE)) ## End(Not run)