DESeq2--compare two log2-transformed fold changes in paired and multi factor sample
Entering edit mode
tntntntntn • 0
Last seen 6 months ago


I have RNA-seq dataset from two population (1and 2) including some patients and two conditions (B and A, before and after treatment) for each person.

My interest is whether there is a difference between population in the effect of treatment. In other words, I would like to compare whether log2(1A/1B) equals log2(2A/2B) or not.

I read DESEQ2--compare two log2-transformed fold changes, and I am generating a model using interaction terms. The difference is that there are some patients in each group and I have to consider paring each of them.

My sample table is like this:

population patient condition

1 p1 B

1 p2 B

1 p3 B

1 p4 B

1 p5 B

1 p1 A

1 p2 A

1 p3 A

1 p4 A

1 p5 A

2 p6 B

2 p7 B

2 p8 B

2 p6 A

2 p7 A

2 p8 A

Since patient is nested in population, I set patient.nested.

population patient condition patient.nested

1 p1 B p1

1 p2 B p2

1 p3 B p3

1 p4 B p4

1 p5 B p5

1 p1 A p1

1 p2 A p2

1 p3 A p3

1 p4 A p4

1 p5 A p5

2 p6 B p1

2 p7 B p2

2 p8 B p3

2 p6 A p1

2 p7 A p2

2 p8 A p3

Then, I would like to create a modelMatrix and remove columns of levels without samples. However, I cannot figure out what model is appropriate for my purpose. Does someone havs suggestions on how to design the model in this situation in DESeq2?

Thanks in advance for your help.



deseq2 • 430 views
Entering edit mode
Last seen 12 months ago

Due to the paired nature of your experiment, I would recommend using limma so that you can use its duplicateCorrelation method to take into account the repeated measures.

For the design matrix, I'd combine condition and population first into condPop:

condPop patient
B1 p1
B1 p2
B1 p3
...and so on

Then use design <- model.matrix(0 + condPop) as your design matrix.

For the contrast matrix, do:

makeContrasts(interaction = (condPopB1 - condPopA1) - (conPopB2 - condPopA2), levels = design)

This will test for the population difference in treatment response.

If you want to use DESeq2, I guess you can do the same but use this model matrix instead:

design <- model.matrix(0 + condPop+ patient)

This way, you are blocking on the patient ID.

Entering edit mode

We have a section of the DESeq2 vignette on how to compare across groups when individuals are nested.

Entering edit mode

Thank you for your instruction! I will read it thoroughly and try it.

Entering edit mode

Thank you for your recommendation! I am new to limma, so I continue analysis with DESeq2 first, but I will try limma later.


Login before adding your answer.

Traffic: 332 users visited in the last hour
Help About
Access RSS

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

Powered by the version 2.3.6