> Date: Fri, 27 Apr 2007 02:31:05 +0100 (BST) > From: Lev Soinov <lev_embl1 at="" yahoo.co.uk=""> > Subject: Re: [BioC] paired design, LIMMA > To: BioC Mailing List <bioconductor at="" stat.math.ethz.ch=""> > Message-ID: <802561.11912.qm at web27902.mail.ukl.yahoo.com> > Content-Type: text/plain > > Dear List, > > Some time ago I posted a question about barley arrays, below. > Gordon Smyth kindly suggested to look at a fitted model MA plot in order to find something > abnormal in the data. Unsure how to create a fitted model MA plot, plotMA(fit_pair, array=4) > I just created MA plots for > RMA normalized data using mva.pairs. I cannot find anything extraordinary in the plots except > they indeed show more significantly upregulated than downregulated genes, but the majority of > the genes is concentrated along the A axe as usual (medians are straight and ~0 for all plots, > IQRs of M values are in the range of 208 - 359). > > I would be grateful for any further suggestions. > Sorry to bother everyone. Well, you've looked more closely at your data, and have seen that the genes are indeed mostly upregulated, just as the statistical analysis told you. It's not clear to me what the problem is that you are seeking further suggestions on. Best wishes Gordon > Thank you, > Lev. > > Gordon Smyth <smyth at="" wehi.edu.au=""> wrote: > >>[BioC] paired design, LIMMA >>Lev Soinov lev_embl1 at yahoo.co.uk >>Fri Apr 20 12:37:40 CEST 2007 >> >>Dear List, >> >> I am learning about a simple paired design in LIMMA and am >> playing with a small dataset of 6 Affymetrix barley chips, 3 >> treated and 3 untreated. I have some problems with interpreting the >> results and would be grateful for any comments/suggestions. >> Experiment: sample pairs (treated & untreated) were prepared in >> three biological replicates, using the same protocol (same >> treatments, etc.) but separately from each other (in different times). >> For all genes with negative fold changes, adj. p values for >> moderated t statistics appear to be higher than 0.1 (the smallest >> adj. p value among down-regulated genes is 0.139). Besides, only >> two down-regulated probes have abs(logFC)>log2(1.5). >> >> Questions: >> 1. From your experience, is the fact that among significantly >> regulated genes I only get up-regulated ones an indication of a >> problem with the data (log_intensity plots and boxplots did not >> flag up any significant problems)? > > Well, if this is biologically infeasible, then it would seem to > indicate a problem. > >> 2. With moderated t-statistics I am getting no significantly >> down-regulated genes but with ordinary t-statistics I get more than >> 4000 down-regulated probes with adj. p <0.05. Is this common? > > Ordinary t-tests typically throw up a lot of spuriously DE genes, > which have very small standard deviations, low fold changes and low > expression levels. > > The difference here between moderated and ordinary t-test suggests to > me that all the apparently down-regulated probes are in the lower > expression range. This does suggest to me a problem with the data. A > fitted model MA-plot might throw some light on the situation. > > Best wishes > Gordon > >> I also wonder why the difference between adj. p values for >> moderated and ordinary t statistics is so huge, i.e. moderated adj. >> p values for down-regulated genes are all higher than 0.1, while >> ~4000 down-regulated probes have ordinary adj.p<0.05. >> >> With kind regards, >> Lev. >> Bioinformatician, UK. >> >> I am using the following code (as described in the LIMMA user >> guide, p.40, 8.3 Paired Samples): >> memory.limit(size = 2048) >> data<-ReadAffy(widget=TRUE) >> sampleNames(data) >> temp<-rma(data) >> targets <- readTargets("Targets.txt") >> Pair <- factor(targets$Pair) >> Treat <- factor(targets$Treatment, levels=c("A","B")) >> design <- model.matrix(~Pair+Treat) >> fit_pair <- lmFit(temp, design) >> fit_pair <- eBayes(fit_pair) >> topTable(fit_pair, coef="TreatB") >> >> My targets file is: >> FileName Pair Treatment >> Bar18 1 A >> Bar19 1 B >> Bar20 2 A >> Bar21 2 B >> Bar22 3 A >> Bar23 3 B >> >> >>Ordinary t statistics for paired test were calculated using: >> >> >tstat.ord <- fit_pair$coef / fit_pair$stdev.unscaled / fit_pair$sigma >> >> >p.value.ord <- 2 * pt( abs(tstat.ord), df=fit_pair$df.residual, >> lower.tail=FALSE) >> >> >pvalue.ord.adj <- p.adjust(p.value.ord, method="fdr") >> >> My data and session info: >> >data >> AffyBatch object >> size of arrays=712x712 features (23770 kb) >> cdf=Barley1 (22840 affyids) >> number of samples=6 >> number of genes=22840 >> annotation=barley1 >> > sessionInfo() >> R version 2.4.1 (2006-12-18) >> i386-pc-mingw32 >> attached base packages: >> [1] "tcltk" "tools" "stats" "graphics" "grDevices" >> "utils" "datasets" "methods" "base" >> other attached packages: >> barley1cdf tkWidgets DynDoc >> widgetTools limma affy affyio Biobase >> "1.14.0" "1.12.1" "1.12.0" "1.10.0" "2.9.8" >> "1.12.2" "1.2.0" "1.12.2" >> >> There are NO significantly down-regulated genes, see the TreatB >> column below: >> > results<-decideTests(fit_pair) >> > summary(results) >> (Intercept) Pair2 Pair3 TreatB >> -1 0 407 161 0 >> 0 4 21977 22453 22819 >> 1 22836 456 226 21 >> >> topTable (probe IDs are removed for brevity): >> ID logFC AveExpr t P.Value >> adj.P.Val B >> 1 1.4 5.9 15.2 0.0000002 0.003 5.3 >> 2 6.1 6.7 14.6 0.0000003 0.003 5.1 >> 3 3.1 6.4 14.3 0.0000003 0.003 5.0 >> 4 1.2 7.6 12.6 0.0000010 0.005 4.6 >> 5 3.0 7.7 12.4 0.0000011 0.005 4.5 >> 6 1.5 4.6 11.7 0.0000017 0.007 4.3 >> 7 1.9 6.7 11.1 0.0000026 0.007 4.0 >> 8 1.0 4.6 11.1 0.0000026 0.007 4.0 >> 9 3.3 5.2 11.0 0.0000029 0.007 4.0 >> 10 1.8 7.5 10.1 0.0000055 0.013 3.6 >> 11 1.5 6.4 9.5 0.0000089 0.018 3.3 >> 12 1.1 8.2 9.4 0.0000097 0.018 3.2 >> 13 3.5 6.4 8.4 0.0000223 0.039 2.7 >> 14 1.0 4.4 8.3 0.0000256 0.040 2.6 >> 15 1.1 6.4 8.3 0.0000260 0.040 2.6 >> 16 1.8 4.6 8.1 0.0000290 0.041 2.5 >> 17 0.9 9.4 8.0 0.0000345 0.044 2.4 >> 18 0.9 6.9 7.9 0.0000349 0.044 2.3 >> 19 4.5 7.2 7.9 0.0000372 0.045 2.3 >> 20 1.0 7.4 7.8 0.0000397 0.045 2.3 >> 21 0.8 7.0 7.7 0.0000437 0.048 2.2 >> 22 2.1 3.4 7.6 0.0000484 0.050 2.1 >> 23 0.9 6.5 7.2 0.0000701 0.070 1.8 >> 24 2.6 6.0 7.2 0.0000734 0.070 1.8 >> 25 0.6 9.0 7.1 0.0000781 0.071 1.7 >> 26 1.2 5.1 7.1 0.0000808 0.071 1.7 >> 27 1.3 6.5 7.0 0.0000867 0.072 1.7 >> 28 1.5 4.8 7.0 0.0000891 0.072 1.6 >> 29 3.8 5.2 6.9 0.0000945 0.072 1.6 >> 30 2.3 6.0 6.9 0.0000959 0.072 1.6

An embedded and charset-unspecified text was scrubbed... Name: not available Url: https://stat.ethz.ch/pipermail/bioconductor/attachments/20070428/ a9110903/attachment.pl