Dear Fishpond team.

I've been taking a look at the Swish (fishpond) package for my DTE/DGE analysis and it seem to fit my needs perfectly, however I have a small question about the scaled counts output.

As I understand, it is possible to plot the results of each transcript/gene with plotInfReps via idx but is it possible to also extract a matrix containing the scaled counts for all transcripts/genes? I would like to make my own plots of the scaled counts result.

Best regards, Fredrik

Thanks for your interest in using Swish.

So here is a typical run of swish:

library(SummarizedExperiment)
library(fishpond)
y <- makeSimSwishData()
y <- scaleInfReps(y)
## Progress: 20 on 20
y <- labelKeep(y)
set.seed(1)
y <- swish(y, x="condition")
## Generating test statistics over permutations
## Progress: 100 on 100

Then we can look at a particular feature (say a gene):

plotInfReps(y, which.min(mcols(y)$qvalue), "condition")

This plots the distribution of the Gibbs samples (or bootstraps) for each sample.

The values that are plotted are stored in a series of matrices: the matrices named "infRep..." in assayNames(y).

If you want the values for a particular feature you can use the code:

infRepMatrix <- unlist( assays(y[idx,])[ grep("infRep", assayNames(y)) ] )

This matrix is nreps x nsamples. You could then compute the column median (the thick lines in the boxplot) or the column mean of these.

If you want to compute this for all features at once you can use the abind function to make an array and then use apply:

infReps <- assays(y)[ grep("infRep", assayNames(y)) ]
infArray <- abind::abind( as.list(infReps), along=3 )

infArray is nfeature x nsamples x nreps. You can compute summaries like so:

infMed <- apply(infArray, 1:2, median)
infMean <- apply(infArray, 1:2, mean)