Question: Limma time-course: how to apply differential contrast if time zero is considered both as treatment and control?
0
21 months ago by
harelarik10
harelarik10 wrote:

# Limma time-course: how to apply differential contrast if time zero is considered both as treatment and control?

Goal: We are using limma to analyze time course experiment with 3 time points, and one treatment vs. control. We would like to test this combination of time course and one cold stress as follows: Amb_0hr (no stress time zero), Amb_1hr, Amb_6hr, Cld_1hr (cold stress after 1 hour), Cld_6hr.

Experiment description: we started with 5 flasks of bacteria; extracted RNA for time zero from one of them; applied cold stress only on two of them; extracted RNA after 1hr from one stressed and one control flasks; extracted RNA after 6 hours form one stressed and one control flasks (everything was done in triplicates).

Issue":  In the contrasts we used (below) time zero is considered both as "time zero none stress", and "time zero cold stress". PROBLEM: If we use the differential contrast suggested below [Dif1_hr =( Cld _1hr-Amb_0hr)-(Amb_1hr-Amb_0hr)] we get more upregulated genes then expected in "standard" contrast" (Cld _1hr-Amb_0hr). In specific, we get 204 upregulated genes for the differential contrast Dif1_hr. However using "standard" contrast we get only 79 upregulated genes (see results table below).

#### A.** This is the differential contrast:**

cont.wt <- makeContrasts(

Dif1_hr =(Cld_1hr-Amb_0hr)-(Amb_1hr-Amb_0hr),

Dif_6hr=(Cld _6hr- Cld _1hr)-(Amb_6hr-Amb_1hr),

levels=design)

#### B.** This is the "standard" contrast":**

cont.wt <- makeContrasts( " Cld _1hr-Amb_0hr"," Cld _6hr- Cld _1hr", "Amb_1hr-Amb_0hr","Amb_6hr-Amb_1hr",levels=design)

#### C. Results Table:

Problem: We get 204 upregulated genes for the differential contrast Dif1_hr. However using "standard" contrast (Cld _1hr-Amb_0hr) we get only 79 upregulated genes:

Dif1_hr          Dif6_hr            Cld_1hr-Amb_0hr       Cld_6hr- Cld_1hr   Amb_1hr-Amb_0hr    Amb_6hr-Amb_1hr

54D/204U      212D/138U      367D/79U                   311D/368U            334D/61U                  42D/55U

*D=Downregulated, U=upregulated.

## What is the best contrats to ask "Which genes respond at either 1 hour or 6 hour of cold stress?

### Detailed ANALYSIS:

1. Data taken from RNA-seq experiment.

2. Design matrix (two time points :1,6 hr; treatments: control, coldStress, triplicate for each treatment):

Amb_0hr Amb_1hr Amb_6hr Cld_1hr Cld_6hr

1        1       0       0      0      0

2        1       0       0      0      0

3        1       0       0      0      0

4        0       1       0      0      0

5        0       1       0      0      0

6        0       1       0      0      0

7        0       0       1      0      0

8        0       0       1      0      0

9        0       0       1      0      0

10       0       0       0      1      0

11       0       0       0      1      0

12       0       0       0      1      0

13       0       0       0      0      1

14       0       0       0      0      1

15       0       0       0      0      1

3. Commands used:

lev<-c("Amb_0hr","Amb_1hr","Amb_6hr","Cld_1hr","Cld_6hr")

all_groups_repeats<-c("Amb_0hr", "Amb_0hr", "Amb_0hr", "Amb_1hr", "Amb_1hr", "Amb_1hr", "Amb_6hr", "Amb_6hr", "Amb_6hr", "Cld_1hr",  " Cld _1hr",  " Cld _1hr",  " Cld _6hr",  " Cld _6hr",  " Cld _6hr")

f <- factor(all_groups_repeats, levels=lev)

design <- model.matrix(~0+f)

colnames(design) <- lev

fit <- lmFit(eset, design)

v <- voom(xCurrent, design, plot=FALSE)#Voom fit

fit <- lmFit(v , design)

#### Contrast used:

cont.wt <- makeContrasts(

Dif1_hr =(Cld _1hr-Amb_0hr)-(Amb_1hr-Amb_0hr),

Dif_6hr=(Cld _6hr-Cld _1hr)-(Amb_6hr-Amb_1hr),

levels=design)

fit2 <- contrasts.fit(fit, cont.wt)

fit2 <- eBayes(fit2)

#### #Final DE analysis

wt<-decideTests(fit2,p.value = 0.05, adjust.method = "BH",  method = "global")

summary(wt)#Produce summary table

limma time-course • 384 views
modified 20 months ago • written 21 months ago by harelarik10
Answer: Limma time-course: how to apply differential contrast if time zero is considered
0
21 months ago by
United States
James W. MacDonald50k wrote:

This is just simple algebra. Your Dif1hr contrast is the same thing as

(x - y) - (z - y)

Which reduces to

x - z

Which is not the same thing as any other contrast you have made, so it's not surprising that you get different results when comparing different comparisons.

Thank you very much!!

#### Unfortunately, our problem is not solved yet.  Could you Please suggest how should we find which genes respond to 1hr of cold stress?

Reminder: Experiment description: we started with 5 flasks of bacteria; extracted RNA for time zero from one of them; applied cold stress only on two of them; extracted RNA after 1hr from one stressed and one control flasks; extracted RNA after 6 hours form one stressed and one control flasks (everything was done in triplicates).

Issue": If we want to find which genes respond to 1 hr cold stress we would like to use the "differential contrast" suggested in page 47 of the user guide. While it is straightforward for 6hr cold stress (i.e.,  makeContrasts( Dif_6hr=(Cld _6hr- Cld _1hr)-(Amb_6hr-Amb_1hr), levels=design) it is not trivial for 1hr cold stress. The reason for that is that,  in our time zero is considered both as "time zero none stress", and "time zero cold stress".

#### Considering this. Coluld you please suggest what is the best way to build differential contrast to 1hr hour cold stress?

I was thinking about the two following options:

1. For 1 hour cold stress:

Use this "standard contrat":

cont.wt <- makeContrasts( "Cld _1hr-Amb_0hr"," Cld _6hr- Cld _1hr", "Amb_1hr-Amb_0hr","Amb_6hr-Amb_1hr",levels=design)

And then to substract the list of over expressed genes found in "Amb_1hr-Amb_0hr" from the list of genes found by  "Cld _1hr-Amb_0hr".

#### Is this a correct approach?

2. For 6 hour cold stress:

I have two options:

a. Using the same appraoch used for 1hr cold stress based on the "standard contrat".

OR

b. Using the differential contrast based on page 47 of the user guide:

cont.wt <- makeContrasts( Dif_6hr=(Cld _6hr- Cld _1hr)-(Amb_6hr-Amb_1hr), levels=design)

Which one is better for 6 hours?