View source: R/24.filter_testOut.R
filter_testOut | R Documentation |
filter results of test_dPDUI()
filter_testOut(
res,
gp1,
gp2,
outdir = getInPASOutputDirectory(),
background_coverage_threshold = 2,
P.Value_cutoff = 0.05,
adj.P.Val_cutoff = 0.05,
dPDUI_cutoff = 0.2,
PDUI_logFC_cutoff = log2(1.5)
)
res |
a UTR3eSet object, output of |
gp1 |
tag names involved in group 1. gp1 and gp2 are used for filtering purpose if both are specified; otherwise only other specified thresholds are used for filtering. |
gp2 |
tag names involved in group 2 |
outdir |
A character(1) vector, a path with write permission for storing InPAS analysis results. If it doesn't exist, it will be created. |
background_coverage_threshold |
background coverage cut off value. for each group, more than half of the long form should greater than background_coverage_threshold. for both group, at least in one group, more than half of the short form should greater than background_coverage_threshold. |
P.Value_cutoff |
cutoff of P value |
adj.P.Val_cutoff |
cutoff of adjust P value |
dPDUI_cutoff |
cutoff of dPDUI |
PDUI_logFC_cutoff |
cutoff of PDUI log2 transformed fold change |
A data frame converted from an object of GenomicRanges::GRanges.
Jianhong Ou, Haibo Liu
test_dPDUI()
library(limma)
path <- system.file("extdata", package = "InPAS")
load(file.path(path, "eset.MAQC.rda"))
tags <- colnames(eset@PDUI)
g <- factor(gsub("\\..*$", "", tags))
design <- model.matrix(~ -1 + g)
colnames(design) <- c("Brain", "UHR")
contrast.matrix <- makeContrasts(
contrasts = "Brain-UHR",
levels = design
)
res <- test_dPDUI(
eset = eset,
method = "limma",
normalize = "none",
design = design,
contrast.matrix = contrast.matrix
)
filter_testOut(res,
gp1 = c("Brain.auto", "Brain.phiX"),
gp2 = c("UHR.auto", "UHR.phiX"),
background_coverage_threshold = 2,
P.Value_cutoff = 0.05,
adj.P.Val_cutoff = 0.05,
dPDUI_cutoff = 0.3,
PDUI_logFC_cutoff = .59
)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.