### R code from vignette source 'vignettes/mmnet/inst/doc/mmnet.Rnw'
### Encoding: UTF-8

###################################################
### code chunk number 1: mmnet.Rnw:22-31
###################################################
library(knitr)
options(width=64,digits=2)
opts_chunk$set(size="small")
opts_chunk$set(tidy=TRUE,tidy.opts=list(width.cutoff=50,keep.blank.line=TRUE))
opts_knit$set(eval.after='fig.cap')
# for a package vignette, you do want to echo.
# opts_chunk$set(echo=FALSE,warning=FALSE,message=FALSE)
opts_chunk$set(warning=FALSE,message=FALSE)
opts_chunk$set(cache=TRUE,cache.path="cache/mmnet")


###################################################
### code chunk number 2: install-pkg (eval = FALSE)
###################################################
## ## install release version of mmnet
## source("http://bioconductor.org/biocLite.R")
## biocLite("mmnet")
## 
## ##install the latest development version 
## useDevel()
## biocLite("mmnet")


###################################################
### code chunk number 3: load-pkg
###################################################
library(mmnet)


###################################################
### code chunk number 4: load-pkg2 (eval = FALSE)
###################################################
## library(mmnet)


###################################################
### code chunk number 5: sample-download (eval = FALSE)
###################################################
## download.file("ftp://ftp.metagenomics.anl.gov/projects/10/4440616.3/raw/507.fna.gz", 
##               destfile = "obesesample.fna.gz")
## download.file("ftp://ftp.metagenomics.anl.gov/projects/10/4440823.3/raw/687.fna.gz", 
##               destfile = "leansample.fna.gz")


###################################################
### code chunk number 6: login-mgrast (eval = FALSE)
###################################################
## ## login on MG-RAST 
##  login.info <- loginMgrast(user="mmnet",userpwd="mmnet")


###################################################
### code chunk number 7: upload-mgrast (eval = FALSE)
###################################################
## ## select the sample data to upload
## seq.file <- c("obesesample.fna.gz", "leansample.fna.gz")
## ## upload sequence data to MG-RAST
## metagenome.id <- lapply(seq.file, uploadMgrast, login.info = login.info)


###################################################
### code chunk number 8: submit-mgrast (eval = FALSE)
###################################################
## ## submit MGRAST job
##   metagenome.id <- lapply(seq.file, submitMgrastJob, login.info = login.info)
##   show(metagenome.id)


###################################################
### code chunk number 9: check-metagenome (eval = FALSE)
###################################################
## ## check MGRAST project status
##   metagenome.status <- lapply(metagenome.id, checkMgrastMetagenome, login.info = login.info)
## ## apparently, status of completed annotation metagenome is TRUE 
##   show(metagenome.status)


###################################################
### code chunk number 10: getannotation (eval = FALSE)
###################################################
## ## private data
##   private.annotation <- lapply(metagenome.id, getMgrastAnnotation, login.info=login.info)
## ## public annotation data, does not require login.info 
##   public.annotation <- lapply(metagenome.id, getMgrastAnnotation)


###################################################
### code chunk number 11: data-load
###################################################
data(anno)
summary(anno)


###################################################
### code chunk number 12: MG-RAST-mmnet (eval = FALSE)
###################################################
## ## first login on MG-RAST
## login.info <- loginMgrast("mmnet", "mmnet")
## ## prepare the metagenomic sequence for annotation
## seq <- "obesesample.fna.gz"
## ## mgrast annotation
## uploadMgrast(login.info, seq)
## metagenome.id2 <- submitMgrastJob(login.info, seqfile = basename(seq)) 
## while (TRUE) { 
##   status <- checkMgrastMetagenome(metagenome.id = metagenome.id2)
##   if (status)
##     break
##   Sys.sleep(5)
##   cat("In annotation, please waiting...")
##   flush.console()
## }
## ## if annotation profile is public,take login.info as NULL
## anno2 <- getMgrastAnnotation(metagenome.id2, login.info=login.info)


###################################################
### code chunk number 13: estimate
###################################################
 mmnet.abund <- estimateAbundance(anno)
 show(mmnet.abund)


###################################################
### code chunk number 14: MG-RAST-BIOM
###################################################
## download BIOM functional abundance profile of the two sample metagenome from MG-RAST
if(require(RCurl)){
  function.api <- "http://api.metagenomics.anl.gov/1/matrix/function"
  mgrast.abund <- read_biom(getForm(function.api,
                             .params=list(id="mgm4440616.3",id="mgm4440823.3",
                                          source="KO",result_type="abundance")))
}
## obtain the intersect ko abundance of MG-RAST and esimatiAbundance
intersect.ko <- intersect(rownames(mgrast.abund),rownames(mmnet.abund))
## compare the two by taking one metagenome
mgrast.abund1 <- biom_data(mgrast.abund)[,1][intersect.ko]
mmnet.abund1 <- biom_data(mmnet.abund)[,1][intersect.ko]
if(require(ggplot2)){
  p <- qplot(mgrast.abund1,mmnet.abund1) + 
            geom_abline(slope=1, intercep=0,color="blue") +
            ylim(0, 400) + xlim(0, 400)
  print(p)
}


###################################################
### code chunk number 15: initial-data
###################################################
loadMetabolicData()
summary(RefDbcache)


###################################################
### code chunk number 16: construct-network
###################################################
  ssn <- constructSSN(mmnet.abund)
  g <- ssn[[1]]
  summary(g)
  abund <- get.vertex.attribute(g,"abundance",index=V(g))
  summary(abund)