### R code from vignette source 'vignettes/motifStack/inst/doc/motifStack.Rnw'

###################################################
### code chunk number 1: DNAseqLogo
###################################################
library(motifStack)
pcm <- read.table(file.path(find.package("motifStack"), "extdata", "bin_SOLEXA.pcm"))
pcm <- pcm[,3:ncol(pcm)]
rownames(pcm) <- c("A","C","G","T")
motif <- new("pcm", mat=as.matrix(pcm), name="bin_SOLEXA")
##pfm object
#motif <- pcm2pfm(pcm)
#motif <- new("pfm", mat=motif, name="bin_SOLEXA")
plot(motif)
#try a different font
plot(motif, font="mono,Courier")
#try a different font and a different color group
motif@color <- colorset(colorScheme='basepairing')
plot(motif,font="Times")


###################################################
### code chunk number 2: fig1
###################################################
opar<-par(mfrow=c(3,1))
motif@color<-colorset(colorScheme='auto')
motif@name="bin_SOLEXA, font='Helvetica', color='auto'"
plot(motif)
motif@name="bin_SOLEXA, font='mono,Courier', color='auto'"
plot(motif, font="mono,Courier")
motif@color <- colorset(colorScheme='basepairing')
motif@name="bin_SOLEXA, font='mono,Courier', color='basepairing'"
plot(motif,font="Times")
par<-opar


###################################################
### code chunk number 3: AAseqLogo
###################################################
library(motifStack)
protein<-read.table(file.path(find.package("motifStack"),"extdata","cap.txt"))
protein<-t(protein[,1:20])
motif<-pcm2pfm(protein)
motif<-new("pfm", mat=motif, name="CAP", 
            color=colorset(alphabet="AA",colorScheme="chemistry"))
plot(motif)


###################################################
### code chunk number 4: fig2
###################################################
plot(motif)


###################################################
### code chunk number 5: seqLogoStack
###################################################
library(motifStack)
#####Input#####
pcms<-readPCM(file.path(find.package("motifStack"), "extdata"),"pcm$")
pcms<-lapply(pcms,function(.ele){.ele<-.ele[,3:ncol(.ele)];rownames(.ele)<-c("A","C","G","T");.ele})
motifs<-lapply(pcms,pcm2pfm)
motifs<-lapply(names(motifs), function(.ele, motifs){new("pfm",mat=motifs[[.ele]], name=.ele)},motifs)

##plot stacks
motifStack(motifs, layout="stack", ncex=1.0)
motifStack(motifs, layout="tree")

###When the number of motifs is too much to be shown in a vertical stack, 
###motifStack can draw them in a radial style.
library("MotifDb")
matrix.fly <- query(MotifDb, "Dmelanogaster")
motifs2 <- as.list(matrix.fly)
motifs2 <- motifs2[grepl("Dmelanogaster\\-FlyFactorSurvey\\-", names(motifs2))]
names(motifs2) <- gsub("Dmelanogaster_FlyFactorSurvey_", "", 
                      gsub("_FBgn\\d+$", "", 
                            gsub("[^a-zA-Z0-9]","_", 
                                 gsub("(_\\d+)+$", "", names(motifs2)))))
motifs2 <- motifs2[unique(names(motifs2))]
pfms <- sample(motifs2, 50)
motifs2 <- lapply(names(pfms), function(.ele, pfms){new("pfm",mat=pfms[[.ele]], name=.ele)},pfms)
library(RColorBrewer)
color <- brewer.pal(12, "Set3")


###################################################
### code chunk number 6: fig3
###################################################
motifStack(motifs, layout="stack", ncex=1.0)


###################################################
### code chunk number 7: fig4
###################################################
motifStack(motifs, layout="tree")


###################################################
### code chunk number 8: fig5
###################################################
motifStack(motifs2, layout="radialPhylog", circle=0.4, cleaves = 0.2, clabel.leaves = 0.5, 
                                 col.bg=rep(color, each=5), col.bg.alpha=0.3, 
                                 col.leaves=rep(color, each=5),
                                 col.inner.label.circle=rep(color,5),
                                 col.outer.label.circle=rep(color,5), outer.label.circle.width=0.1, 
                                 angle=350)


###################################################
### code chunk number 9: motifCloud
###################################################
groups <- rep(paste("group",1:5,sep=""), each=10)
names(groups) <- names(pfms)
group.col <- brewer.pal(5, "Set3")
names(group.col)<-paste("group",1:5,sep="")
jaspar.scores <- MotIV::readDBScores(file.path(find.package("MotIV"), "extdata", "jaspar2010_PCC_SWU.scores"))
d <- MotIV::motifDistances(pfms)
hc <- MotIV::motifHclust(d)
phylog <- hclust2phylog(hc)
leaves <- names(phylog$leaves)
pfms <- pfms[leaves]
pfms <- lapply(names(pfms), function(.ele, pfms){new("pfm",mat=pfms[[.ele]], name=.ele)},pfms)
motifSig <- motifSignature(pfms, phylog, groupDistance=0.1)


###################################################
### code chunk number 10: fig6
###################################################
motifCloud(motifSig, scale=c(6, .5), layout="rectangles", group.col=group.col, groups=groups, draw.legend=T)


###################################################
### code chunk number 11: motifStack.Rnw:227-228
###################################################
sessionInfo()


