Deseq2 Paired Samples Design
1
0
Entering edit mode
Rita • 0
@rita-24908
Last seen 3.7 years ago
United States

Hi everyone,

I am an undergraduate at the University of Virginia and new to RStudio. I am in need of assistance regarding my DESeq2 analysis. I would like to use DESeq2 to process paired RNASeq samples. I have over 100 samples from male and females patients. For each patient, I have two treatments (With-FBS and Without-FBS). I would like to run DESeq2 to identify the differentially expressed genes between males and females irrespective of their treatment. My question is, how do I tell DESeq2 that each pair of With-FBS and Without-FBS treated samples belong to one patient?

I have included my metadata(coldata) and countdata below in addition to the interaction model I think need based on what I described above.

Whenever I include Donor anywhere in the interaction model, I get this error message

                    ddsBA <- DESeqDataSetFromMatrix(countData = GctscountBA, 
                                                      colData = metaDataBA, 
                                                       design = ~ Donor + Treatment:Sex)

Error in checkFullRank(modelMatrix) : the model matrix is not full rank, so the model cannot be fit as specified. One or more variables or interaction terms in the design formula are linear combinations of the others and must be removed.

                             MetaData

                                  Donor  Treatment Sex 
              A79_WithFBS          A79    WithFBS   M
              A86_WithFBS          A86    WithFBS   M
              A88_WithFBS          A88    WithFBS   M
              A95_WithFBS          A95    WithFBS   M
              A96_WithFBS          A96    WithFBS   F
              B78_WithFBS          B78    WithFBS   F
              A79_WithoutFBS       A79  WithoutFBS  M
              A86_WithoutFBS       A86  WithoutFBS  M
              A88_WithoutFBS       A88  WithoutFBS  M
              A95_WithoutFBS       A95  WithoutFBS  M
              A96_WithoutFBS       A96  WithoutFBS  F
              B78_WithoutFBS       B78  WithoutFBS  F

                             CountData

                             A79_WithFBS A86_WithFBS A88_WithFBS A95_WithFBS A96_WithFBS B78_WithFBS A79_WithoutFBS A86_WithoutFBS A88_WithoutFBS A95_WithoutFBS A96_WithoutFBS B78_WithoutFBS
             WASH7P              20          17          53          19          49          49           28             29           17              9                 44          17
             AL627309.5           5          12          21           4          13           6           13             18            5              6                  9          28
             WASH9P              32          57          19          52          53          43           83             27           57             28                178          43

When I run this model below, it runs however, it treats each patient treatment as individual patients for example, A79_WithFBS and A79_WithoutFBS are considered different patients which is not correct.

                    ddsBA <- DESeqDataSetFromMatrix(countData = GctscountBA, 
                                                      colData = metaDataBA, 
                                                       design = ~ Treatment + Treatment:Sex)

What do I have to change to let DESeq2 know that (A79_WithFBS and A79_WithoutFBS) are from the same donor(A79) and subsequently investigate the differentially expressed genes irrespective of treatment?

Thank you for your help!

DESeq2 • 2.1k views
ADD COMMENT
3
Entering edit mode
@mikelove
Last seen 5 days ago
United States

Actually, if you want to compare across sex, and also control for the two samples from each donor, we don't have a mechanism for this in DESeq2. It is essentially a mixed effects model with a random effect for the donor, nested within the sex grouping.

You can alternatively use limma-voom with the duplicateCorrelation() function (see limma vignette).

Or with DESeq2 you could compare males and females within the control samples, or within the treated samples.

ADD COMMENT
0
Entering edit mode

Hello Mr.Love,

Thank you for your assistance, I am now in the process of running limma-voom. I have already compared males and females within the treatment groups.

Best, Rita

ADD REPLY
0
Entering edit mode

Hello Mr.Love, thank you again for the advice. Follow up questions:

Will my model and contrast look like this (below) in that case? If so, how does limma voom know that I am interested in the DEG between males and females (how will it know that Sex is my parameter)?

design <- model.matrix(~ 0 + Treatment ,metaDatamalesMFBatch128)

contrasts <- makeContrasts(TreatmentWithFBS-TreatmentWithoutFBS, levels=colnames(design))

voom_dge <- voom(GctscountMFBatch128, design, plot=TRUE)

block = metaDatamalesMFBatch128$Donor

cor <- duplicateCorrelation(voom_dge, design, block = block )

cor$consensus.correlation

---------------------------------------------------------------------------------------------------------

OR will my parameter also be a factor, in this case, Sex to get DEG between males and females?

design <- model.matrix(~ 0 + Sex+Treatment ,metaDatamalesMFBatch128)

contrasts <- makeContrasts(SexM-SexF, levels=colnames(design))

voom_dge <- voom(GctscountMFBatch128, design, plot=TRUE)

block = metaDatamalesMFBatch128$Donor

cor <- duplicateCorrelation(voom_dge, design, block = block )

cor$consensus.correlation

---------------------------------------------------------------------------------------------------------

ADD REPLY
0
Entering edit mode

hi Rita,

I'm not a developer of this other package. I'd recommend creating a new post with your description of the problem from the start.

ADD REPLY
0
Entering edit mode

Okay, thank you.

ADD REPLY

Login before adding your answer.

Traffic: 383 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6