results {denoiSeq}

Compute the test statistic

Description

Extracts posterior samples of the parameters which are returned by denoiseq function and computes the summary and test statistics.

Usage

results(RDobject, steps, burnin = floor(steps/3), rope_limit = 0.5)

Arguments

RDobject

A readsData object with a filled output slot.

steps

An integer representing the number of iterations.

burnin

An integer for the number of iterations to be considered as burn in values. A default value equivalent to a third of steps is used.

rope_limit

A float that delimits the range of the region of practical equivalence, ROPE. A default value of 0.5 is used.

Details

To calculate the test statistic, this function first log2 transforms the posterior samples of the two relevant parameters i.e N_{iA} and N_{iB}. It then randomly subtracts posterior samples of one of the parameters from the other and determines the proportion of this distribution of differences that lies in the region of practical equivalence (ROPE) (Kruschke, 2011). The genes can then be arranged in an ascending order of the ROPE_propn column and we can select the most differentially expressed genes as those whose ROPE_propn is less than a particular threshold value.

Using both real and simulated data, optimal values between 0.0007 and 0.4 were obtained for the threshold.

Value

A dataframe with 3 columns namely; the log2 fold change (log2FC), the standard error of the log2 fold change (lgfcSE) and the test static (ROPE_propn).

Examples

#pre -filtering to remove lowly expressed genes
ERCC <- ERCC[rowSums(ERCC) > 0, ]
RD <- new('readsData', counts = ERCC)
steps <- 30
#30 steps are just for illustration here. At least 5000 steps are adequate.
BI <- denoiseq(RD, steps)
rez <- results(BI, steps)
head(rez)

#Re-ordering according to most differentially expressed
rez <- rez[with(rez, order( ROPE_propn)), ]
head(rez, 10)

#Determine significant genes using a threshold of 0.38.
sgf <- rez[rez$ROPE_propn<0.38, ]
head(sgf)
dim(sgf)