edgeR testing interactions between two factors while controlling for blocking factor
0
0
Entering edit mode
Mafalda • 0
@89d39968
Last seen 12 weeks ago
Sweden

Hi,

I am sorry if this is a repeated question, but I couldn't find any example that fits my exact case.

I have a dataset of individuals with alternative genotypes (North or Heterozygous). Each genotype was reared under different temperatures (Seven or Ten). I am interested in understanding the effect of temperature by genotype and of genotype by temperature. Initially I constructed the following model matrix in edgeR following section 3.3.1:

design$group<-factor(paste0(design$TEMP, ".", design$GENO))
design_model <- model.matrix(~ 0 + group)
colnames(design_model) <- levels(group)

  SEVEN.Heterozygous SEVEN.North TEN.Heterozygous TEN.North
1                   1           0                0         0
2                   1           0                0         0
3                   1           0                0         0
4                   1           0                0         0
5                   1           0                0         0
6                   1           0                0         0
7                   1           0                0         0
8                   1           0                0         0
9                   1           0                0         0
10                  1           0                0         0
11                  1           0                0         0
12                  1           0                0         0
13                  0           1                0         0
14                  0           1                0         0
15                  0           1                0         0
16                  0           1                0         0
17                  0           1                0         0
18                  0           1                0         0
19                  0           0                1         0
20                  0           0                1         0
21                  0           0                1         0
22                  0           0                1         0
23                  0           0                1         0
24                  0           0                1         0
25                  0           0                1         0
26                  0           0                1         0
27                  0           0                1         0
28                  0           0                1         0
29                  0           0                1         0
30                  0           0                0         1
31                  0           0                0         1
32                  0           0                0         1
33                  0           0                0         1
34                  0           0                0         1
35                  0           0                0         1
36                  0           0                0         1
37                  0           0                0         1
38                  0           0                0         1
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$group
[1] "contr.treatment"

And testing differential expression for the effect of, let's say, genotype for a given temperature:

# Calculate dispersion
df <- estimateDisp(df, model_design, robust=TRUE)
df <- estimateGLMCommonDisp(df, model_design)
# Fit the model
fit <- glmQLFit(df, model_design, robust=TRUE)
# Make a contrast
con.SEVEN.North_vs_Het <- makeContrasts(groupSEVEN.North - groupSEVEN.Heterozygous, levels=design_model)
qlf.SEVEN.North_vs_Het <- glmQLFTest(fit, contrast=con.SEVEN.North_vs_Het)

Now I wanted to make this a bit more complicated:

1. Testing for overall effect of temperature or genotype

With the model matrix above, is it correct to test for the overall effect of temperature like this?

qlf.Ten_vs_Seven <- glmQLFTest(fit, contrast=c(-1,-1,1,1)

Or should I instead construct a model matrix that includes temperature and genotype explicitly (see code bellow)? But then how do I test for the effect of genotype within temperature, or temperature within genotype, like above; doo I need to have a separate model matrix for this interaction?

TEMP<-factor(paste(design$TEMP))
GENO<-factor(paste(design$GENO))
design_model <- model.matrix(~ GENO+TEMP)

2. How to add a blocking factor?

Sex doesn't seem to have a strong effect on gene expression in this dataset (all groups of GENO x TEMP have equivalent numbers of males and females), but I thought I should still control for this effect. What would be the best way to do this in the model matrix? From the first model matrix, should I just add SEX as a blocking factor in the design matrix above like so:

design$group<-factor(paste0(design$TEMP, ".", design$GENO))
SEX<-factor(paste(design$SEX))
design_model <- model.matrix(~ 0 + group + SEX)

And then proceeding as above to make the contrasts that I want? I am able to fit the model with this model matrix, but I am not sure it is the correct way of adding the blocking factor.

So I also tried this code following the section 3.3.2 of the manual and inspired by this question, but I get a warning that my Design matrix is not of full rank

SEX<-factor(paste(design$SEX))
TEMP<-factor(paste(design$TEMP))
GENO<-factor(paste(design$GENO))
design_model <- model.matrix(~ 0 + SEX + GENO:TEMP)

  SEXF SEXM GENOHeterozygous:TEMPSEVEN GENONorth:TEMPSEVEN GENOHeterozygous:TEMPTEN GENONorth:TEMPTEN
1     0    1                          1                   0                        0                 0
2     1    0                          1                   0                        0                 0
3     0    1                          1                   0                        0                 0
4     1    0                          1                   0                        0                 0
5     0    1                          1                   0                        0                 0
6     0    1                          1                   0                        0                 0
7     1    0                          1                   0                        0                 0
8     1    0                          1                   0                        0                 0
9     0    1                          1                   0                        0                 0
10    1    0                          1                   0                        0                 0
11    0    1                          1                   0                        0                 0
12    0    1                          1                   0                        0                 0
13    1    0                          0                   1                        0                 0
14    0    1                          0                   1                        0                 0
15    1    0                          0                   1                        0                 0
16    1    0                          0                   1                        0                 0
17    1    0                          0                   1                        0                 0
18    0    1                          0                   1                        0                 0
19    1    0                          0                   0                        1                 0
20    0    1                          0                   0                        1                 0
21    0    1                          0                   0                        1                 0
22    0    1                          0                   0                        1                 0
23    1    0                          0                   0                        1                 0
24    1    0                          0                   0                        1                 0
25    0    1                          0                   0                        1                 0
26    1    0                          0                   0                        1                 0
27    1    0                          0                   0                        1                 0
28    1    0                          0                   0                        1                 0
29    0    1                          0                   0                        1                 0
30    1    0                          0                   0                        0                 1
31    0    1                          0                   0                        0                 1
32    1    0                          0                   0                        0                 1
33    0    1                          0                   0                        0                 1
34    1    0                          0                   0                        0                 1
35    0    1                          0                   0                        0                 1
36    0    1                          0                   0                        0                 1
37    0    1                          0                   0                        0                 1
38    0    1                          0                   0                        0                 1
attr(,"assign")
[1] 1 1 2 2 2 2
attr(,"contrasts")
attr(,"contrasts")$SEX
[1] "contr.treatment"

attr(,"contrasts")$GENO
[1] "contr.treatment"

attr(,"contrasts")$TEMP
[1] "contr.treatment"

Then I tried explicitly adding GENO and TEMP:

SEX<-factor(paste(design$SEX))
TEMP<-factor(paste(design$TEMP))
GENO<-factor(paste(design$GENO))
design_model <- model.matrix(~ SEX + GENO + TEMP + GENO:TEMP)

  (Intercept) GENONorth TEMPTEN SEXM GENONorth:TEMPTEN
1            1         0       0    1                 0
2            1         0       0    0                 0
3            1         0       0    1                 0
4            1         0       0    0                 0
5            1         0       0    1                 0
6            1         0       0    1                 0
7            1         0       0    0                 0
8            1         0       0    0                 0
9            1         0       0    1                 0
10           1         0       0    0                 0
11           1         0       0    1                 0
12           1         0       0    1                 0
13           1         1       0    0                 0
14           1         1       0    1                 0
15           1         1       0    0                 0
16           1         1       0    0                 0
17           1         1       0    0                 0
18           1         1       0    1                 0
19           1         0       1    0                 0
20           1         0       1    1                 0
21           1         0       1    1                 0
22           1         0       1    1                 0
23           1         0       1    0                 0
24           1         0       1    0                 0
25           1         0       1    1                 0
26           1         0       1    0                 0
27           1         0       1    0                 0
28           1         0       1    0                 0
29           1         0       1    1                 0
30           1         1       1    0                 1
31           1         1       1    1                 1
32           1         1       1    0                 1
33           1         1       1    1                 1
34           1         1       1    0                 1
35           1         1       1    1                 1
36           1         1       1    1                 1
37           1         1       1    1                 1
38           1         1       1    1                 1
attr(,"assign")
[1] 0 1 2 3 4
attr(,"contrasts")
attr(,"contrasts")$GENO
[1] "contr.treatment"

attr(,"contrasts")$TEMP
[1] "contr.treatment"

attr(,"contrasts")$SEX
[1] "contr.treatment"

And this works with fitting the model, but then I am not sure how to test for differences between temperatures within genotypes, or differences between genotypes within temperatures because I only have the GENONorth:TEMPTEN, and not GENONorth:TEMPSEVEN, GENOHeterozygous:TEMPTEN, or GENOHeterozygous:TEMPSEVEN...

Any advice?

Thank you so much!

edgeR ExpressionData RNASeqData • 107 views
ADD COMMENT

Login before adding your answer.

Traffic: 436 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