importDittoBulk {dittoSeq} | R Documentation |
import bulk sequencing data into a format that dittoSeq functions expect.
importDittoBulk(x, ...) ## S4 method for signature 'SummarizedExperiment' importDittoBulk(x, reductions = NULL, metadata = NULL, combine_metadata = TRUE) ## S4 method for signature 'DGEList' importDittoBulk(x, reductions = NULL, metadata = NULL, combine_metadata = TRUE) ## S4 method for signature 'list' importDittoBulk(x, reductions = NULL, metadata = NULL)
x |
a |
... |
For the generic, additional arguments passed to specific methods. |
reductions |
a named list of dimensionality reduction embeddings matrices.
names will become the names of the dimensionality reductions and how each will be used with the |
metadata |
a data.frame like object containing columns of extra information about the cells/samples (rows). The names of these columns can then be used to tretrieve and plot such data in any dittoSeq visualizations. |
combine_metadata |
Logical which sets whether original colData (DESeqDataSet/SummarizedExperiment) or $samples (DGEList) from x should be retained. |
A SingleCellExperiment
object containing all assays (DESeqDataSet or SummarizeedExperiment) or all common slots (DGEList) of the input x
,
as well as any dimensionality reductions provided to reductions
, and any provided metadata
stored in colData.
When combine_metadata
is set to FALSE, metadata inside x
(colData or $samples) is ignored entirely.
When combine_metadata
is TRUE (the default), metadata inside x
is combined with what is provided to the metadata
input; but names must be unique, so when there are similarly named slots, the values provided to the metadata
input are used.
One recommended assay to create if it is not already present in your dataset, is a log-normalized version of the counts data. The logNormCounts function of the scater package is an easy way to make such a slot. dittoSeq defaults to grabbing expression data from an assay named logcounts > normcounts > counts
SingleCellExperiment
for more information about this storage system.
## Bulk data is stored as a SingleCellExperiment library(SingleCellExperiment) # Generate some random data nsamples <- 60 exp <- matrix(rpois(1000*nsamples, 20), ncol=nsamples) colnames(exp) <- paste0("sample", seq_len(ncol(exp))) rownames(exp) <- paste0("gene", seq_len(nrow(exp))) logexp <- log2(exp + 1) # Dimensionality Reductions pca <- matrix(runif(nsamples*5,-2,2), nsamples) tsne <- matrix(rnorm(nsamples*2), nsamples) # Some Metadata conds <- factor(rep(c("condition1", "condition2"), each=nsamples/2)) timept <- rep(c("d0", "d3", "d6", "d9"), each = 15) genome <- rep(c(rep(TRUE,7),rep(FALSE,8)), 4) grps <- sample(c("A","B","C","D"), nsamples, TRUE) clusts <- as.character(1*(tsne[,1]>0&tsne[,2]>0) + 2*(tsne[,1]<0&tsne[,2]>0) + 3*(tsne[,1]>0&tsne[,2]<0) + 4*(tsne[,1]<0&tsne[,2]<0)) ### We can import the counts directly, or as a SummarizedExperiment myRNA <- importDittoBulk( x = list(counts = exp, logcounts = logexp)) ### Adding metadata & PCA or other dimensionality reductions # We can add these directly during import, or after. myRNA <- importDittoBulk( x = list(counts = exp, logcounts = logexp), metadata = data.frame( conditions = conds, timepoint = timept, SNP = genome, groups = grps), reductions = list( pca = pca)) myRNA$clustering <- clusts myRNA <- addDimReduction( myRNA, embeddings = tsne, name = "tsne") # (other packages SCE manipulations can also be used) ### When we import from SummarizedExperiment, all metadata is retained. # The object is just 'upgraded' to hold extra slots. # The input is the same, aside from a message when metadata are replaced. se <- SummarizedExperiment( list(counts = exp, logcounts = logexp)) myRNA <- importDittoBulk( x = se, metadata = data.frame( conditions = conds, timepoint = timept, SNP = genome, groups = grps, clustering = clusts), reductions = list( pca = pca, tsne = tsne)) myRNA ### For DESeq2, how we might have made this: # DESeqDataSets are SummarizedExperiments, and behave similarly # library(DESeq2) # dds <- DESeqDataSetFromMatrix( # exp, data.frame(conditions), ~ conditions) # dds <- DESeq(dds) # dds_ditto <- importDittoBulk(dds) ### For edgeR, DGELists are a separate beast. # dittoSeq imports what I know to commonly be inside them, but please submit # an issue on the github (dtm2451/dittoSeq) if more should be retained. # library(edgeR) # dgelist <- DGEList(counts=exp, group=conditions) # dge_ditto <- importDittoBulk(dgelist)