[BioC] Paired limma analysis

mb3058 at columbia.edu mb3058 at columbia.edu
Fri Aug 3 00:59:13 CEST 2012


I am running a paired analysis with limma.
My platform is affymetrix mouse gene st.1.0 array
I have a paired group (hi, lo) with duplicates per group.
My script so far is this:
exprs <- read.table("data.txt",as.is=T,sep="\t",header=T,row.names=1)
targets<-readTargets("pairedTargetFile.txt")
Pair <- factor(targets$Group)
Treat <- factor(targets$Treatment, levels=c("hi","lo"))
design <- model.matrix(~Pair+Treat)
fit <- lmFit(exprs, design)
fit <- eBayes(fit)
#I want to print out the limma results for all the probes
Results<- topTable(fit,number=241425,adjust.method="BH")
write.table(Results,"paired_Results.txt",sep="\t")

ID	X.Intercept.	Pair2	Treatlo	AveExpr	F	P.Value	adj.P.Val
10469363	10.590048	2.186194	-0.636504	11.364893	1566.532105	0.000000469	0.00000942
10365564	10.12985925	2.3047915	-0.2214415	11.17153425	1528.543391	0.000000495	0.00000942
10469367	9.8800415	2.440832	-0.652667	10.774124	1417.807515	0.000000583	0.00000942
10517519	10.1252575	2.183145	0.212485	11.3230725	1403.521849	0.000000596	0.00000942
10469380	9.61192575	2.9290725	-0.8228465	10.66503875	1396.225476	0.000000603	0.00000942

Q1.I am trying to interpret some of the fields. Could someone help me  
understand what the fields X.Intercept.	Pair2	Treatlo stand for?



Q2. When I run the following line to get results
Results <- topTable(fit, coef="Treatlo"), I get the following results:
ID	logFC	AveExpr	t	P.Value	adj.P.Val	B
10354461	4.390143	5.666826	14.06184	8.62E-05	0.5279716	-0.7055664
10583079	4.249853	6.16608	13.51989	1.02E-04	0.5279716	-0.7296509
10587318	4.180087	4.887135	13.34857	1.08E-04	0.5279716	-0.7378131
10370298	4.0321	5.502014	12.69299	1.33E-04	0.5279716	-0.7717809
10392551	4.259277	4.614051	12.44214	1.45E-04	0.5279716	-0.7860399
10368320	4.387814	5.911313	12.36954	1.49E-04	0.5279716	-0.7903085
10587807	3.805685	4.826062	12.30147	1.52E-04	0.5279716	-0.7943703
10504815	3.793085	5.657937	12.27763	1.54E-04	0.5279716	-0.7958069
10357489	3.968974	4.666258	12.23582	1.56E-04	0.5279716	-0.7983445
10481521	3.757882	7.856306	12.15789	1.60E-04	0.5279716	-0.8031341

However when I run the following command to get results for all the  
probes, I get different values for the same probes. See PValue, Adj. P  
value
ID	X.Intercept.	Pair2	Treatlo	AveExpr	F	P.Value	adj.P.Val
10354461	2.963005	1.0175	4.390143	5.6668265	441.9193584	7.37E-06	1.48E-05
10583079	5.2467745	-2.411241	4.249853	6.1660805	516.5649337	5.25E-06	1.24E-05
10587318	2.88690775	-0.1796325	4.1800875	4.88713525	334.0547737	1.35E-05	2.19E-05
10370298	4.54191225	-2.1118965	4.0321005	5.50201425	408.7688335	8.74E-06	1.64E-05
10392551	2.4626765	0.043472	4.259277	4.614051	261.1260083	2.31E-05	3.27E-05
10368320	3.76314325	-0.0914745	4.3878135	5.91131275	377.3505032	1.04E-05	1.83E-05
10587807	2.69620575	0.4540285	3.8056845	4.82606225	325.4919972	1.43E-05	2.28E-05
10504815	5.27270925	-3.0226305	3.7930855	5.65793675	458.5242305	6.81E-06	1.42E-05
10357489	2.90779925	-0.4520565	3.9689735	4.66625775	286.4197106	1.89E-05	2.80E-05
10481521	6.637718	-1.320705	3.757882	7.8563065	794.2002708	2.06E-06	9.55E-06



More information about the Bioconductor mailing list