Question

DESeq2 - 50 vs 50 analysis

0

Entering edit mode

Karim • 0

@ed4b33d3

Last seen 2.7 years ago

Russia

Hi I need help the analysis of my RNA-seq samples using DESeq2. I uploaded cancer evidence for 100 samples divided into two parts, each section contains 50, the first group for a dead patient and the second for a living patient, all samples are cancerous. Data from the site GDC , HTSeq - Counts. I want to extract the expressed genes.

> -- replacing outliers and refitting for 3827 genes

> -- DESeq argument 'minReplicatesForReplace' = 7 

> -- original counts are preserved in counts(dds)

and running : minReplicatesForReplace = Inf and cooksCutoff = FALSE ,Get only 200

Code :


    #groups <- data.frame(type = c(rep("Dead", 50), rep("Alive", 50)))
    #rownames(groups) <- colnames(counts)

    #dds <- DESeqDataSetFromMatrix(countData = counts, colData = groups, design = ~ type)
    #dds <- DESeq(dds)

DESeq2 • 1.2k views

ADD COMMENT • link updated 3.0 years ago by Michael Love 43k • written 3.0 years ago by Karim • 0

0

Entering edit mode

Did you check the DESeq2 vignette ? http://www.bioconductor.org/packages/release/bioc/vignettes/DESeq2/inst/doc/DESeq2.html

ADD REPLY • link 3.0 years ago Basti ▴ 780

0

Entering edit mode

Thank Basti Yes, I check out the DESeq2 vignette, my problem is that I get as few as 476 genes, and when I use (minReplicatesForReplace = Inf and cooksCutoff = FALSE) I only get 200 genes.

ADD REPLY • link 3.0 years ago Karim • 0

Michael Love · Answer 1 · 2022-02-15

0

Entering edit mode

Michael Love 43k

@mikelove

Last seen 1 day ago

United States

This is a lot of outliers. It sounds like the data are not a match for the model you are trying to fit. I would suggest trying RUV or sva to model unwanted variation (see the rnaseqGene workflow).

ADD COMMENT • link 3.0 years ago Michael Love 43k

0

Entering edit mode

Hello Michael . and thank you for your reply. I read ( the rnaseqGene workflow) , and used the code (

#genes <- rownames(counts)[grep("^ENSG", rownames(counts))]
#x <- as.factor(rep(c("Dead", "Alive"), each=50))
#set <- newSeqExpressionSet(as.matrix(counts),phenoData = data.frame(x, row.names=colnames(counts)))
#set1 <- RUVg(set, genes, k=2)

#dds <- DESeqDataSetFromMatrix(countData = counts(set1),colData = pData(set1),design = ~ W_1 + x)
#dds <- DESeq(dds)

and got 470 gin. Is writing the code in this way correct or is there another way to write the code.

ADD REPLY • link updated 3.0 years ago by Michael Love 43k • written 3.0 years ago by Karim • 0

0

Entering edit mode

Why didn't you follow the code in the workflow? You've arbitrarily skipped steps.

https://www.bioconductor.org/packages/devel/workflows/vignettes/rnaseqGene/inst/doc/rnaseqGene.html#using-ruv-with-deseq2

ADD REPLY • link 3.0 years ago Michael Love 43k