R design on 2-way interaction not displaying all options
1
0
Entering edit mode
14 months ago
fr ▴ 170

I have a design given by treatment and categorical stages, each with 5 replicates that I want to model in limma, that looks like this:

     treatment      stage
1    ctrl              A
2    ctrl              B
3    ctrl              C
4    ctrl              D
1    drug              A
2    drug              B
3    drug              C
4    drug              D


I want to model the individual effect of treatment, stage, as well as their interaction. To do so I'm considering design = ~0+ treatment*stage. I then specify the contrasts that I'm interested, namely:

drug - ctrl
D - B
(drug_D - drug_B) - (ctrl_D - ctrl_B)


so I specify model_mat = model.matrix( eval( parse ( text = "~0+ treatment + stage" ) ), data=some_df) .

This generates the model matrix that then I want to use to build my contrasts, for the 3 different contrasts above. However, some fail because not all columns are present in model_mat:

#works
makeContrasts(contrasts = 'treatmentdrug - treatmentctrl', levels = model_mat)

#works
makeContrasts(contrasts = 'stageD - stageB', levels = model_mat)

# no columns in model_mat contain treatmentdrg:stage
# only columns with treatmentctrl:stage are present
# (fails because treatmentdrug:stageD and treatmentdrug:stageB are not present)
#
makeContrasts(
contrasts = '(treatmentdrug:stageD - treatmentdrug:stageB) - (treatmentctrl:stageD - treatmentctrl:stageB)',
levels=model_mat)


I understand that some of the columns are not there as they can be given by combinations of the other factors, but this becomes trickier in more complex designs such as the last contrast that I want to make. Is there a way to overcome this issue?

I see that the alternative given the tests I want to make could be to specify design1=~0 + treatment + stage which would let me model the first 2 contrats, and design1=~0 + treatment:stage which would let me model the last. But in the case of design2, then the model does not take into account the individual effect of each treatment and stage, but only their combination, which I don't think it is correct mathematically.

In other words, what I am asking is how to express all contrasts that are not present? (e.g. treatmentdrug:stageB).

From what I understand, makeContrasts gives everything with reference to the first alphabetical group, so no treatmentdrug will show up for the interactions if I take a different design such as those above. In this case how would you compute treatmentdrug:stageD - treatmentdrug:stageB) - (treatmentctrl:stageD - treatmentctrl:stageB)

Any suggestions on how to solve these issues? Thanks in advance

modeling design limma deseq edger • 470 views
0
Entering edit mode
14 months ago
ashish ▴ 640

If you do treatment * stage it will use both main effects and their interaction. It is equivalent to treatment + stage + treatment:stage i.e both main effects and their interaction.

0
Entering edit mode

Yes, I am aware of this, this is not what I'm asking. I'm asking is how to express all contrasts that are not present? makeContrasts gives everything with reference to the first alphabetical group, so no treatmentdrug will show up for the interactions if I take any of your options.

I understand this, everything is given with respect to a reference. But then how do you make the contrast of (treatmentdrug:stageD - treatmentdrug:stageB) - (treatmentctrl:stageD - treatmentctrl:stageB)`

I've made this clearer in my question.