swish {fishpond}R Documentation

swish: SAMseq With Inferential Samples Helps

Description

swish: SAMseq With Inferential Samples Helps

Usage

swish(y, x, cov = NULL, pair = NULL, interaction = FALSE,
  nperms = 30, estPi0 = FALSE, qvaluePkg = "qvalue", pc = 5,
  nRandomPairs = 30, fast = 1, quiet = FALSE)

Arguments

y

a SummarizedExperiment containing the inferential replicate matrices of median-ratio-scaled TPM as assays 'infRep1', 'infRep2', etc.

x

the name of the condition variable. A factor with two levels for a two group analysis (possible to adjust for covariate or matched samples, see next two arguments)

cov

the name of the covariate for adjustment. If provided a stratified Wilcoxon in performed. Cannot be used with pair

pair

the name of the pair variable, which should be the number of the pair. Can be an integer or factor. If specified, a signed rank test is used to build the statistic. All samples across x must be pairs if this is specified. Cannot be used with cov.

interaction

logical, whether to perform a test of an interaction between x and cov. These are different than the other tests produced by the software, in that they focus on a difference in the log2 fold change across levels of x when comparing the two levels in cov. If pair is specified, this will perform a Wilcoxon rank sum test on the two groups of matched sample LFCs. If pair is not included, multiple random pairs of samples within the two groups are chosen, and again a Wilcoxon rank sum test compared the LFCs across groups.

nperms

the number of permutations

estPi0

logical, whether to estimate pi0

qvaluePkg

character, which package to use for q-value estimation, samr or qvalue

pc

pseudocount for finite estimation of log2FC, not used in calculation of test statistics, locfdr or qvalue

nRandomPairs

the number of random pseudo-pairs (only used with interaction=TRUE and un-matched samples) to use to calculate the test statistic

fast

an integer, toggles different methods based on speed (fast=1 is default). '0' involves recomputing ranks of the inferential replicates for each permutation, '1' is roughly 10x faster by avoiding re-computing ranks for each permutation. The fast argument is only used/relevant for the following three experimental designs: (1) two group Wilcoxon, (2) stratified Wilcoxon, e.g. cov is specified, and (3) the paired interaction test, e.g. pair and cov are specified. For paired design and general interaction test, there are not fast/slow alternatives.

quiet

display no messages

Value

a SummarizedExperiment with metadata columns added: the statistic (either a centered Wilcoxon Mann-Whitney or a signed rank statistic, aggregated over inferential replicates), a log2 fold change (the median over inferential replicates, and averaged over pairs or groups (if groups, weighted by sample size), the local FDR and q-value, as estimated by the samr package.

References

The citation for swish method is:

Anqi Zhu, Avi Srivastava, Joseph G Ibrahim, Rob Patro, Michael I Love "Nonparametric expression analysis using inferential replicate counts" Nucleic Acids Research (2019). https://doi.org/10.1093/nar/gkz622

The swish method builds upon the SAMseq method, and extends it by incorporating inferential uncertainty, as well as providing methods for additional experimental designs (see vignette).

For reference, the publication describing the SAMseq method is:

Jun Li and Robert Tibshirani "Finding consistent patterns: A nonparametric approach for identifying differential expression in RNA-Seq data" Stat Methods Med Res (2013). https://doi.org/10.1177/0962280211428386

Examples


library(SummarizedExperiment)
set.seed(1)
y <- makeSimSwishData()
y <- scaleInfReps(y)
y <- labelKeep(y)
y <- swish(y, x="condition")

# histogram of the swish statistics
hist(mcols(y)$stat, breaks=40, col="grey")
cols = rep(c("blue","purple","red"),each=2)
for (i in 1:6) {
  arrows(mcols(y)$stat[i], 20,
         mcols(y)$stat[i], 10,
         col=cols[i], length=.1, lwd=2)
}

# plot inferential replicates
plotInfReps(y, 1, "condition")
plotInfReps(y, 3, "condition")
plotInfReps(y, 5, "condition")


[Package fishpond version 1.2.0 Index]