#Programm zur Auswertung von Hybridisierungen auf dem Weizen-Chip, Version 1.0"

library(limma)

# input of the target file
targets <- readTargets("CMDON4h targetSW.txt")


# show the targets
targets

# read all columns of the RG object from the .gpr files, spots with flag get a weight"
RG <- read.maimages(targets$FileName, source="genepix", wt.fun=wtflags(0))
#read gal file
RG$genes <- readGAL("PICME_wheat_3K_1_3_6_gal.txt")

# set print layout information
RG$printer <- getLayout(RG$genes)


# check that data are read in correctly
show (RG)

# scatter plot 
plotMA(RG, array=1)
plotMA(RG, array=2)
plotMA(RG, array=3)
plotMA(RG, array=4)




#show the content of the RG object from spot 1:3
RG$genes[1:3]


#imageplot(log2(RG$Rb[,1]), RG$printer, low="white", high="red", zlim=c(-3,3), main="red background, Array 1")
#imageplot(log2(RG$Rb[,1]), RG$printer, low="white", high="red", zlim=c(-1,1))
#imageplot(log2(RG$Gb[,1]), RG$printer, low="white", high="green", zlim=c(-1,1))
#fit the array to the three replicates, change RG to RG1 and only read rows up to 18
RG1 <- RG
RG1$genes <- RG$genes[RG$genes$Row<19,]
RG1$M <- RG$M[RG$genes$Row<19,]
RG1$A <- RG$A[RG$genes$Row<19,]
RG1$weights <- RG$weights[RG$genes$Row<19,]
RG1$R <- RG$R[RG$genes$Row<19,]
RG1$Rb <- RG$Rb[RG$genes$Row<19,]
RG1$G <- RG$G[RG$genes$Row<19,]
RG1$Gb <- RG$Gb[RG$genes$Row<19,]
RG1$nspot.r <- 18
RG1$printer$nspot.r <-18
#show file structure (RG1$printer)
RG1$printer
#show file structure (RG1)
RG1
objects()
dim(RG1$genes)
dim(RG1$R)
dim(RG1$Rb)
dim(RG1$Gb)
dim(RG1$G)
# normalization within each array
MA <- normalizeWithinArrays(RG1, method="loess")

# show plot with print tip
plotPrintTipLoess(MA)

# showing the normalized data (each array are in one box, if the variation is different proceed with normalizatio between arrays))"
boxplot(MA$M~col(MA$M),names=colnames(MA$M))

# normalization between all arrays
MA <- normalizeBetweenArrays(MA)

# showing the normalized data (box plot should be similar)
boxplot(MA$M~col(MA$M),names=colnames(MA$M))


# show the image of the array with M(logR/logG array1) 
imageplot(MA$M[,1], MA$printer, low="green", high="red", main="jo9")

# show the image of the array with M(logR/logG array 2
imageplot(MA$M[,2], MA$printer, low="green", high="red", main="jo2")


#design <- modelMatrix(targets, ref="CMH20")
design

# linear model with replication in the arrays
fit <- lmFit(MA, design, ndups=3, spacing=1,)



fit2 <-contrasts.fit(fit, contrast.matrix)
fit2 <-eBayes(fit)

#showing the best genes
topTable(fit2, n=35, adjust="fdr")

#writing the results to an external file in the R folder
write.fit(fit2, file="CMDON4h limma.txt", digits=2, adjust="fdr")