[BioC] M values

[Sean Davis]

On Tue, Sep 30, 2008 at 9:42 AM, Jose Francisco Rodriguez
<francquimico at gmail.com> wrote:
> Hi:
> I performed a microarray (5 biological replicates) and analyzed using
> Limma.  The software gives me the M value as a combination from the 5
> biological replicates microarray.  In order to submit the data to GEO
> database we need the M-value for each array ( 5 M-values).  How can I do
> that? I'll appreciate your help.

MA$M contains the M values for each array.

Sean

>  Below there is a copy of the sequence for the analysis:
>
> library(limma)
> library(geneplotter)
>
> ###Lee targets y image files
> targets=readTargets("maTargetFKS12.txt")
> RG=read.maimages(files=cbind(targets$FileNameCy3,targets$FileNameCy5),
>  source="imagene")
>
> ##lear gal
> gal=readGAL("013384_D_20050601.gal")
> gal=gal[1:nrow(RG$R),]
> ##asignar el gal
> RG$genes=gal
>
> ###definir los Spot Types y colores
> spottypes <- readSpotTypes() ##usamos el default que es SpotTypes
> RG$genes$Status <- controlStatus(spottypes, RG)
>
> ###definir el layour
> printer= list(ngrid.r=1,ngrid.c=1,nspot.r=107,nspot.c=101)
> RG$printer=structure(printer, class = "PrintLayout")
>
> ###normalize: background y normlize within
> RGoriginal = RG
> RG <- backgroundCorrect(RGoriginal, method="normexp", offset=50)
> MA = normalizeWithinArrays(RG,method="loess",bc.method="none")
>
> ##escoger los genes nadamas
> geneIndex = which(MA$genes$Status=="gene")
>
> ##preparar el disenyo para fit
> design = modelMatrix(targets,ref="WT")
>
> ###ver todos los MA plots como QC: todos OK!
> for(i in 1:6){
>  plotMA(MA[geneIndex,i],main=targets$Name[i],ylim=range(MA$M[geneIndex,]))
>  scan()
> }
>
> ###estimar los paramteros usando eBayes (moderated t-stat) y calcular FDR
> fit1 = lmFit(MA,design)
> fit1 = eBayes(fit1)
> tt=topTable(fit1,adjust="fdr",number=3000)
> tt=tt[tt$P.Value<=0.01,]
>
> ###Escoger genes para incluir en las graficas
> Index=as.numeric(rownames(tt))
> interestGenes <-
> c("CHS3","SLT2","PIR3","ECM4","ECM13","SPI1","ORF:YHR097C","MYO1","PIR3","HOG1","SWE1")
> interestIndex <- lapply(interestGenes,function(x)
> which(MA$genes$GeneName==x))
> names(interestIndex) <- interestGenes
>
> ##MA-PLOT
> plot(fit1$Amean[geneIndex],fit1$coef[geneIndex],pch=".",col="grey",ylab="M",xlab="A")
> points(fit1$Amean[Index],fit1$coef[Index],pch=16,col="blue",cex=.5)
> for(i in seq(along=interestGenes)){
>
> text(fit1$Amean[interestIndex[[i]]],fit1$coef[interestIndex[[i]]],interestGenes[i],cex=.5)
> }
>
> ##VOLCANO
> plot(fit1$coef[geneIndex],-log10(fit1$p.value[geneIndex]),pch=".",col="grey",xlab="M",ylab="-log_10
> p-value")
> points(fit1$coef[Index],-log10(fit1$p.value[Index]),pch=16,col="blue",cex=.5)
> for(i in seq(along=interestGenes)){
>
> text(fit1$coef[interestIndex[[i]]],-log10(fit1$p.value[interestIndex[[i]]]),interestGenes[i],cex=.5)
> }
>
>
> ###hacer xls file
> todos = topTable(fit1,adjust="fdr",number=10807)
> write.csv(todos,file="myo1.csv")
> --
> José F. Rodríguez Quiñones
> Graduate student
> School of Medicine
> Medical Sciences Campus
> University of Puerto Rico
>
>        [[alternative HTML version deleted]]
>
>
> _______________________________________________
> Bioconductor mailing list
> Bioconductor at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor
>