Entering edit mode
6.9 years ago
Lila M
★
1.2k
Hi everybody. I'm working with RNASeq data and I am relatively new in that field. I would like to know if my approach is correct or not (or how I can improve it). Assuming that I have 3 different condition or treatment (high, low and mut) what I did is:
dds <- DESeqDataSetFromMatrix(
countData = count.table,
colData = data.frame(treatment=treatment),
design = ~ treatment)
vsd <- varianceStabilizingTransformation(dds, blind=TRUE)
vst <- assay(vsd)
colnames(vst) <- colnames(count.table)
# To avoid lengthy unnecessary explanations, we simply shift all values so that non expressed genes have a vst
#value of 0
vst <- vst - min(vst)
dds1 <- DESeq(dds)
res.highCNA.low_CNA <- results(dds1, contrast=c("treatment","High_CNA", "low_CNA"))
res.highCNA.mut <- results(dds1, contrast=c("treatment","High_CNA", "mut"))
res.lowCNA.mut <- results(dds1, contrast=c("treatment", "low_CNA", "mut"))
#to select the 500 most expressed genes
mat.high.low <- assay(vsd)[head(order(res.highCNA.low_CNA$padj), 500), ]
pheatmap(mat.high.low)
mat.high.mut <- assay(vsd)[head(order(res.highCNA.mut$padj), 500), ]
pheatmap(mat.high.mut)
mat.low.mut <- assay(vsd)[head(order(res.lowCNA.mut$padj), 500), ]
pheatmap(mat.low.mut)
Thanks in advance!!!