Question: Seeking confirmation on manually setting contrasts for a 3-factor design in DESeq2 with multiple interaction terms
1
gravatar for piensaglobalmente
2.9 years ago by
Australia
piensaglobalmente10 wrote:

Hello,

I am seeking confirmation on manually setting contrasts for a 3-factor design in DESeq2 with multiple interaction terms

I have a 3-factor experiment, with each factor containing 2 levels:

Type: Juvenile vs Sediment

Shore: Inshore vs Offshore

Sector: North vs Central

I am interested in the following contrasts:

  1. Average main effects of Type, Shore and Sector across all the levels.
  2. Interactive effect of Type:Shore, specifically:
    1. Juvenile: Inshore vs Offshore
    2. Sediment: Inshore vs Offshore
  3. Interactive effect of Type: Sector, specifically:
    1. Juvenile: North vs Central
    2. Sediment: North vs Central
  4. Interactive effect of Type vs Shore:Sector, specifically:
    1. Juvenile: Inshore North vs. Inshore Central
    2. Juvenile: Inshore North vs. Offshore Central
    3. Juvenile: Inshore North vs. Offshore North
    4. Juvenile: Offshore Central vs. Inshore Central
    5. Juvenile: Offshore North vs. offshore Central
    6. Juvenile: Offshore North vs. Inshore Central

       5. 1-6: Sediment: with same 6 comparisons above

 

My experimental design formula is as follows:

~Type*Shore*Sector

Base levels for all factors are: Juvenile, Inshore, Central

And my Model Matrix column names outputted with resultsNames function:  
resultsNames(varstab_geomLove_T.S.S)

[1] "Intercept"                            

[2] "Type_Sediment_vs_Juvenile"            

[3] "Shore_Offshore_vs_Inshore"            

[4] "Sector_North_vs_Central"               

[5] "TypeSediment.ShoreOffshore"            

[6]"TypeSediment.SectorNorth"             

[7] "ShoreOffshore.SectorNorth"  
  
[8]"TypeSediment.ShoreOffshore.SectorNorth"

Effect

Function

Ok?

Main effect of Type (aka Juvenile vs Sediments)

results(varstab2, contrast=c("Type", "Juvenile", "Sediment"))

yes

Main effect of Shore (aka Inshore vs Offshore)

results(varstab2, contrast=c("Shore", "Inshore", "Offshore"))

yes

Main effect of Sector (aka North or Central)

results(varstab2, contrast=c("Sector", "North", "Central"))

yes

Juvenile: Inshore vs Offshore

 

Subtract column with all inshore vs offshore comparison from sediment: inshore vs offshore

 

[3] "Shore_Offshore_vs_Inshore"     1,0,1,0,0,0,0,0

[5] "TypeSediment.ShoreOffshore"  1,0,0,0,1,0,0,0

Subtract the two to get the following:

results(varstab2, contrast=c(0,0,1,0,-1,0,0,0))      

 

 

Not sure

Sediment: Inshore vs Offshore

 

results(varstab2, name="TypeSediment.ShoreOffshore")

yes

Juvenile: North vs Central

Subtract column with all north vs central comparison from sediment: north vs central

 

[4] "Sector_North_vs_Central"     1,0,0,1,0,0,0,0

[6]"TypeSediment.SectorNorth"  1,0,0,0,0,1,0,0

Subtract the two to get the following:

results(varstab2, contrast=c(0,0,0,1,0,-1,0,0))      

 

Not sure

Sediment: North vs Central

results(varstab2, name=”TypeSediment.SectorNorth”)

yes

Juvenile: Inshore North vs. Inshore Central

[1] "Intercept"   1,0,0,0,0,0,0,0  which is juvenile, inshore, central

[4] "Sector_North_vs_Central"   0,0,0,1,0,0,0,0

 

Subtract the two to get the following:

results(varstab2, contrast=c(1,0,0,-1,0,0,0,0))

 

Not sure

Juvenile: Inshore North vs. Offshore Central

 

?

Juvenile: Inshore North vs. Offshore North

 

?

Juvenile: Offshore Central vs. Inshore Central

 

?

Juvenile: Offshore North vs. offshore Central

[7] "ShoreOffshore.SectorNorth"   

results(varstab2, name="ShoreOffshore.SectorNorth")

Ok

Juvenile: Offshore North vs. Inshore Central

 

?

Sediments: Inshore North vs. Inshore Central

 

?

Sediments: Inshore North vs. Offshore Central

 

?

Sediments: Inshore North vs. Offshore North

 

?

Sediments: Offshore Central vs. Inshore Central

 

?

Sediments: Offshore North vs. offshore Central

[8]"TypeSediment.ShoreOffshore.SectorNorth"

results(varstab2, name="TypeSediment.ShoreOffshore.SectorNorth")

Ok

Sediments: Offshore North vs. Inshore Central

 

?

 

I feel confident in the comparisons with the a "Yes" in the last column. However, could someone please check that I have specified the contrasts in the rows with "Not sure?" Especially the contrasts "Juvenile: Inshore North vs. Inshore Central." If I can get a confirmation that that contrast is set correctly, then I think I can manage the rest.

Thank you for any insights!

cheers,

k

ADD COMMENTlink modified 2.9 years ago by Michael Love22k • written 2.9 years ago by piensaglobalmente10
Answer: Seeking confirmation on manually setting contrasts for a 3-factor design in DESe
1
gravatar for Michael Love
2.9 years ago by
Michael Love22k
United States
Michael Love22k wrote:

hi,

This is a good question for a statistical consultant or a local statistician who you could collaborate with.

There is nothing specific about DESeq2 here, you would use the same contrasts for a normal linear model. So it's not really a DESeq2 software question.

ADD COMMENTlink written 2.9 years ago by Michael Love22k

Thank you Michael. I can keep looking through the glmm literature now that you confirm that contrasts in DESeq2 and other negative binomial models are specified the same. I just wasn't sure if there was some unique DESeq2 functionality I wasn't finding that would draw out these contrasts.  

cheers,

 

ADD REPLYlink written 2.9 years ago by piensaglobalmente10

Hello Michael:

I would like to propose a possible solution to such contrast questions and run it by you.

To estimate higher order contrasts:

1) Construct a design matrix in GLM coding. An example is given here:

 https://support.sas.com/documentation/cdl/en/statug/63347/HTML/default/viewer.htm#statug_glm_sect030.htm

2) Each categorical column (“A1”, “AB21”, etc) in the design can be seen as a separate binary factor. I can also pretend that every column is a quantitative factor (please let me know if it makes a difference for DESeq2). Either way, an additive model that includes all of such factors is fed to DESeq2. Since interactions are actually present, shrinkage of regression coefficients is disabled (betaPrior = FALSE).

3) Obtain the contrast vector as a difference between the two corresponding lsmeans vectors. E.g. for “(A1, B1) vs (A3, B1)” contrast the vector is LSM(AB11) – LSM(AB31). How to obtain lsmeans is described here:

https://support.sas.com/documentation/cdl/en/statug/63962/HTML/default/viewer.htm#statug_glm_a0000000871.htm

4) Feed the contrast vector to results() function from DESeq2. Here I hope to rely on the option: “a numeric contrast vector with one element for each element in resultsNames(object) (most general case)”.

 

Please tell me if DESeq2 can digest an overparameterized design “as is” and whether the results() will work the way I expect.

ADD REPLYlink written 19 months ago by Nik Tuzov70

You can pass in any full rank design matrix to the 'full' argument, and use numeric contrasts to compare coefficients using results(). The current release of DESeq2 has betaPrior=FALSE by default (see NEWS).

ADD REPLYlink written 19 months ago by Michael Love22k

If the user can only feed full rank designs to DESeq2, do you know if there is an R package that could help figure out how to specify contrast coefficients for non-trivial cases? If not, it's very hard to use the design produced by model.matrix() if the goal is to estimate non-trivial contrasts.

ADD REPLYlink written 19 months ago by Nik Tuzov70

The contrast package works with 'lm' like objects. So it doesn't work with DESeq2 but you could extract the linear contrast from the output of the contrast() function from that package, for example running on fake data for 'y' but using the same model.matrix. See the example of the contrast package here.

For non-trivial contrasts I find it just as hard to specify the desired contrast using that package, as working it out myself.

ADD REPLYlink written 19 months ago by Michael Love22k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 326 users visited in the last hour