Hello, I have a question about which analysis to do in LIMMA. I did three experiments (three biological replicates) each with one treatment microcosm ("L") and one control microcosm ("D"). For each experiment, I sampled both the treatment and control microcosms at 5 time points (0 hr, 2 hr, 6 hr, 12 hr, 24 hr). With 3 replicate experiments * (1 treatment + 1 control) * 5 time points, I have data from 30 microarrays. Based on section 8.6 of the LIMMA manual (which is incredibly helpful thank you!) here is the way I set up my targets table: FileName Target EF731D0.CEL D.0hr EF731D2.CEL D.2hr EF731D6.CEL D.6hr EF731D12.CEL D.12hr EF731D24.CEL D.24hr EF731L0.CEL L.0hr EF731L2.CEL L.2hr EF731L6.CEL L.6hr EF731L12.CEL L.12hr EF731L24.CEL L.24hr EF813D0.CEL D.0hr EF813D2.CEL D.2hr EF813D6.CEL D.6hr EF813D12.CEL D.12hr EF813D24.CEL D.24hr EF813L0.CEL L.0hr EF813L2.CEL L.2hr EF813L6.CEL L.6hr EF813L12.CEL L.12hr EF813L24.CEL L.24hr EF815D0.CEL D.0hr EF815D2.CEL D.2hr EF815D6.CEL D.6hr EF815D12.CEL D.12hr EF815D24.CEL D.24hr EF815L0.CEL L.0hr EF815L2.CEL L.2hr EF815L6.CEL L.6hr EF815L12.CEL L.12hr EF815L24.CEL L.24hr I would like to determine which genes changed over time in the treatment relative to the control at each time point. Is the following code correct to answer this question?? cont.dif=makeContrasts(Dif2hr=(L.2hr-L.0hr)-(D.2hr-D.0hr), Dif6hr=(L.6hr-L.0hr)-(D.6hr-D.0hr), Dif12hr=(L.12hr-L.0hr)-(D.12hr-D.0hr), levels=design) fitdif=contrasts.fit(fit, cont.dif) fitdif=eBayes(fitdif) topTable(fitdif, genelist=genetext, adjust="BH") I also would like to determine which genes were differentially expressed between the treatment (L) and control (D) at each time point but I'm not sure which of the following two options to use?? option #1: contLD.matrix=makeContrasts(L.0hr-D.0hr, L.2hr-D.2hr, L.6hr-D.6hr, L.12hr-D.12hr, L.24hr-D.24hr, levels=design) fitLDcontmatrix=contrasts.fit(fit, contLD.matrix) fitLDcontmatrixB=eBayes(fitLDcontmatrix) OR option #2: Since each biological replicate experiment starts with the same bacterial culture (half exposed to the treatment while half are kept as a control), I started to think that the treatment and control samples from each experiment are paired. This made me think about doing paired sample tests for each time point. So, following the paired sample section of the LIMMA manual I created a different targets table and did the following code for each time point: targets: FileName Rep Treatment EF731D0.CEL 1 D EF731L0.CEL 1 L EF813D0.CEL 2 D EF813L0.CEL 2 L EF815D0.CEL 3 D EF815L0.CEL 3 L code: Rep = factor(targets$Rep) # this is like SibShip in the LIMMA example Treat = factor(targets$Treatment, levels=c("D","L")) design = model.matrix(~Rep+Treat) fit0 = lmFit(eset_t0, design) fit0 = eBayes(fit0) topTable(fit0, coef="TreatL", genelist=genetext, adjust="BH") Thanks in advance! Lauren

Hi Lauren, On Thursday, October 24, 2013 11:53:12 PM, Lauren Sassoubre wrote: > Hello, > I have a question about which analysis to do in LIMMA. I did three experiments (three biological replicates) each with one treatment microcosm ("L") and one control microcosm ("D"). For each experiment, I sampled both the treatment and control microcosms at 5 time points (0 hr, 2 hr, 6 hr, 12 hr, 24 hr). With 3 replicate experiments * (1 treatment + 1 control) * 5 time points, I have data from 30 microarrays. > > Based on section 8.6 of the LIMMA manual (which is incredibly helpful thank you!) here is the way I set up my targets table: > FileName Target > EF731D0.CEL D.0hr > EF731D2.CEL D.2hr > EF731D6.CEL D.6hr > EF731D12.CEL D.12hr > EF731D24.CEL D.24hr > EF731L0.CEL L.0hr > EF731L2.CEL L.2hr > EF731L6.CEL L.6hr > EF731L12.CEL L.12hr > EF731L24.CEL L.24hr > EF813D0.CEL D.0hr > EF813D2.CEL D.2hr > EF813D6.CEL D.6hr > EF813D12.CEL D.12hr > EF813D24.CEL D.24hr > EF813L0.CEL L.0hr > EF813L2.CEL L.2hr > EF813L6.CEL L.6hr > EF813L12.CEL L.12hr > EF813L24.CEL L.24hr > EF815D0.CEL D.0hr > EF815D2.CEL D.2hr > EF815D6.CEL D.6hr > EF815D12.CEL D.12hr > EF815D24.CEL D.24hr > EF815L0.CEL L.0hr > EF815L2.CEL L.2hr > EF815L6.CEL L.6hr > EF815L12.CEL L.12hr > EF815L24.CEL L.24hr > > I would like to determine which genes changed over time in the treatment relative to the control at each time point. Is the following code correct to answer this question?? Yes, if I understand your question correctly. These are all interaction terms, where you are looking for genes that respond differently over time to treatment as compared to control. > > cont.dif=makeContrasts(Dif2hr=(L.2hr-L.0hr)-(D.2hr-D.0hr), Dif6hr=(L.6hr-L.0hr)-(D.6hr-D.0hr), Dif12hr=(L.12hr-L.0hr)-(D.12hr-D.0hr), levels=design) > fitdif=contrasts.fit(fit, cont.dif) > fitdif=eBayes(fitdif) > topTable(fitdif, genelist=genetext, adjust="BH") > > I also would like to determine which genes were differentially expressed between the treatment (L) and control (D) at each time point but I'm not sure which of the following two options to use?? > > option #1: > contLD.matrix=makeContrasts(L.0hr-D.0hr, L.2hr-D.2hr, L.6hr-D.6hr, L.12hr-D.12hr, L.24hr-D.24hr, levels=design) > fitLDcontmatrix=contrasts.fit(fit, contLD.matrix) > fitLDcontmatrixB=eBayes(fitLDcontmatrix) > > OR > > option #2: Since each biological replicate experiment starts with the same bacterial culture (half exposed to the treatment while half are kept as a control), I started to think that the treatment and control samples from each experiment are paired. This made me think about doing paired sample tests for each time point. So, following the paired sample section of the LIMMA manual I created a different targets table and did the following code for each time point: > With cell culture, you are somewhere in between biological and technical replicates, regardless. In other words each cell culture is in its own milieu, so hypothetically you are introducing some sort of biological variability, but that is far different from one would usually consider a true biological replicate. As with all statistical analyses the first question you have to answer is what assumptions you are willing to make, and why? You could blindly assume that the samples are highly correlated, and pair them, or you could assume that they are not that correlated, and assume independence. Or you could use the duplicateCorrelation function in limma to estimate the correlation structure, and determine from that what assumptions you think are reasonable. Best, Jim > targets: > FileName Rep Treatment > EF731D0.CEL 1 D > EF731L0.CEL 1 L > EF813D0.CEL 2 D > EF813L0.CEL 2 L > EF815D0.CEL 3 D > EF815L0.CEL 3 L > > code: > Rep = factor(targets$Rep) # this is like SibShip in the LIMMA example > Treat = factor(targets$Treatment, levels=c("D","L")) > design = model.matrix(~Rep+Treat) > fit0 = lmFit(eset_t0, design) > fit0 = eBayes(fit0) > topTable(fit0, coef="TreatL", genelist=genetext, adjust="BH") > > Thanks in advance! > Lauren > _______________________________________________ > Bioconductor mailing list > Bioconductor at r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor -- James W. MacDonald, M.S. Biostatistician University of Washington Environmental and Occupational Health Sciences 4225 Roosevelt Way NE, # 100 Seattle WA 98105-6099