Entering edit mode
3.1 years ago
Pahi
▴
30
I'm working on my RNA-seq pipeline, with deseq2, and I would like to use Glimma2, for interactive visualisation. My problem is that, after I do the DESeq analysation, I can't plot it. On different sites, It is said, that it is as simple as: glimmaMA(dds)
But I get the following error:
Error in buildXYData(table, status, main, display.columns, anno, counts, :
Length of groups must be equal to the number of columns in counts.
So far I tried adding annotation to my dds with:
annot <- mapIds(EnsDb.Hsapiens.v86, gene_id, keytype = "GENEID", column="SYMBOL",multiVals="first")
Where gene_id contains all my id-s. But still, I get the error.
Here is my code:
library(DESeq2)
counts <- mydata.txt
conditions <- factor(c(rep("untreated", 3), rep("treated", 3)))
coldata <- data.frame(row.names=colnames(Counts), conditions)
dds_data <- DESeqDataSetFromMatrix(countData = Counts,
colData = coldata,
rowData = gene_id,
design = ~ conditions)
annot <- mapIds(EnsDb.Hsapiens.v86, gene_id, keytype = "GENEID", column="SYMBOL",multiVals="first")
keep <- rowSums(counts(dds_data)) >= 10
dds_data <- dds_data[keep,]
dds <- DESeq(dds)
mcols(dds) <- DataFrame(annot)
res <- results(dds, contrast=c("conditions","treated","untreated"))
resultsNames(dds)
resLFC <- lfcShrink(dds, coef="conditions_untreated_vs_treated", type="apeglm")
resOrdered <- res[order(res$pvalue),]
plotMA(res, ylim=c(-2,2)) #working
plotMA(resLFC, ylim=c(-2,2)) #working
library(Glimma)
glimmaMDS(dds) #working
glimmaMA(dds)
glimmaVolcano(dds)
glimmaMA(dds)
43467 of 60605 genes were filtered out in DESeq2 tests
Error in buildXYData(table, status, main, display.columns, anno, counts, :
Length of groups must be equal to the number of columns in counts.
What could be the problem? Thank you in advance!
Same problem here, did you manage to solve it somehow?
I changed to edger for my analysis, becouse later on with deseq I had a lot of other issues.
If I remember correcly, in glimma you have to define some parameters manually. For my edger input it looks like:
Glimma::glimmaMA(edger_data, dge = y)
, where edger_data is an lrt list of my data, which is separated from the original edger output.