###################################################
### chunk: loadandExt
###################################################
library("hgu95av2")
chrVec = unlist(as.list(hgu95av2CHR))
table(chrVec)
class(chrVec)
#names(chrVec)[1:10]


###################################################
### chunk: splitCHR
###################################################
byChr = split(names(chrVec), chrVec)
sapply(byChr, length)


###################################################
### chunk: chrY
###################################################
byChr[["Y"]]


###################################################
### chunk: map
###################################################
library("hgu95av2")
hgu95av2MAP$"1001_at"


###################################################
### chunk: eApplyex
###################################################
myPos = eapply(hgu95av2MAP, function(x) grep("^17p", x, value=TRUE))
myPos = unlist(myPos)
length(myPos)


###################################################
### chunk: ans2
###################################################
myFindMap = function(mapEnv, which) {
  myg = ppc(which)
  a1 = eapply(mapEnv, function(x) grep(myg, x, value=TRUE))
  unlist(a1)
}


###################################################
### chunk: simpleCbind
###################################################
x = matrix(1:6, nc=2, dimnames=list(letters[1:3], LETTERS[1:2]))
y = matrix(21:26, nc=2, dimnames=list(letters[6:8], LETTERS[3:4]))
cbind(x,y)
rbind(x,y)


###################################################
### chunk: Stackex
###################################################
s1 = list(a=1:3,b= 11:12,c= letters[1:6])
ss = stack(s1)
ss
unsplit(s1, ss[,2])


###################################################
### chunk: chrMapEx
###################################################
mapP = as.list(hgu95av2MAP)
mLens = unlist(eapply(hgu95av2MAP, length))


###################################################
### chunk: summaryChr
###################################################
mlt = table(mLens)
mlt


###################################################
### chunk: 
###################################################
len3 = mLens[mLens==3]
hgu95av2SYMBOL[[names(len3)[1]]]
hgu95av2MAP[[names(len3)[1]]]


###################################################
### chunk: 
###################################################
len2 = names(mLens[mLens==2])
len2EG = unlist(mget(len2, hgu95av2ENTREZID))
len2EG = len2EG[!duplicated(len2EG)]
len2 = len2[!duplicated(len2EG)]
mapP = mget(len2, hgu95av2MAP)
hasX = sapply(mapP, function(x) if( length(grep("^X", x)) == 1) TRUE
else FALSE)

hasY = sapply(mapP, function(x) if( length(grep("^Y", x)) == 1) TRUE
else FALSE)
table(hasX & hasY)


###################################################
### chunk: missingMap
###################################################
missingMap = unlist(eapply(hgu95av2MAP, 
    function(x) any(is.na(x))))
table(missingMap)


###################################################
### chunk: chrMapEx2
###################################################
mapPs = sapply(mapP, function(x) x[1])
mapPs = mapPs[!is.na(mapPs)]

mapByPos = split(names(mapPs), mapPs)
table(sapply(mapByPos, length))


###################################################
### chunk: DBIEx
###################################################
library("RSQLite")
m = dbDriver("SQLite")
con = dbConnect(m, dbname="test")
data(USArrests)
dbWriteTable(con, "USArrests", USArrests, overwrite = TRUE)
dbListTables(con)


###################################################
### chunk: DBIresultSets
###################################################
rs = dbSendQuery(con, "select * from USArrests")
d1 = fetch(rs, n = 5)
d1
dbHasCompleted(rs)
dbListResults(con)
d2 = fetch(rs, n = -1)
dbHasCompleted(rs)
dbClearResult(rs)


###################################################
### chunk: DBIsimple
###################################################
dbListTables(con)
dbListFields(con, "USArrests")


###################################################
### chunk: SQLitegetalltables
###################################################
query = paste("SELECT name FROM sqlite_master WHERE",
    "type='table' ORDER BY name;")
rs = dbSendQuery(con, query)
fetch(rs, n= -1)


###################################################
### chunk: conditional selection
###################################################
rs = dbSendQuery(con, 
    "SELECT * FROM USArrests WHERE Murder > 10")


###################################################
### chunk: DBIcleanup
###################################################
unlink("test")


###################################################
### chunk: setupSQLite
###################################################
 library("RSQLite")
 m = dbDriver("SQLite")
 ##open up our test db    
 testDB = system.file("extdata/hgu95av2-sqlite.db", package="RBioinf")
 con = dbConnect(m, dbname = testDB)

 tabs = dbListTables(con)
 tabs 
 dbListFields(con, tabs[2])


###################################################
### chunk: SFex
###################################################
 query = paste("SELECT go_ont.go_id, go_ont.ont,", 
     "go_ont_name.ont_name FROM go_ont,",
     "go_ont_name WHERE (go_ont.ont = go_ont_name.ont)")
 rs = dbSendQuery(con, query)
 f3 = fetch(rs, n=3)
 f3
 dbClearResult(rs)


###################################################
### chunk: ijsol
###################################################
 query = paste("SELECT acc_num, go_id", 
     "FROM acc, go_probe", 
     "WHERE (acc.affy_id = go_probe.affy_id)")
 rs = dbSendQuery(con, query)
 f3 = fetch(rs, n=3)
 f3
 dbClearResult(rs)


###################################################
### chunk: SFex2
###################################################
query = paste("SELECT g1.*, g2.evi FROM go_probe g1,",
    "go_probe g2 WHERE  (g1.go_id = 'GO:0005737' ",
    "AND g2.go_id = 'GO:0005737') ",
    "AND (g1.affy_id = g2.affy_id) ",
    "AND (g1.evi = 'IDA' AND g2.evi = 'ISS')")
 rs = dbSendQuery(con, query)
 fetch(rs)


###################################################
### chunk: loadDBpkg
###################################################
library("hgu95av2.db")
mycon = hgu95av2_dbconn()


###################################################
### chunk: GOEX
###################################################
colnames(hgu95av2GO)
toTable(hgu95av2GO)[1:10,]
Lkeys(hgu95av2GO)[1:10]
Rkeys(hgu95av2GO)[1:10]


###################################################
### chunk: showLinks
###################################################
links(hgu95av2GO)[1:10,]


###################################################
### chunk: revMapEx
###################################################
is(hgu95av2SYMBOL, "Bimap")
rmMAP = revmap(hgu95av2SYMBOL)
rmMAP$"ABL1"


###################################################
### chunk: revmapList
###################################################
 myl=list(a="w", b="x", c="y")
 revmap(myl)


###################################################
### chunk: SYM2EG
###################################################
queryAlias = function(x) {
  it = paste("('", paste(x, collapse="', '"), "'", sep="")
  paste("select _id, alias_symbol from alias",
        "where alias_symbol in", it, ");")
}

queryGeneinfo = function(x) {
  it = paste("('", paste(x, collapse="', '"), "'", sep="")
  paste("select _id, symbol from gene_info where",
        "symbol in", it, ");")
}

queryGenes = function(x) {
  it = paste("('", paste(x, collapse="', '"), "'", sep="")
  paste("select * from genes where _id in", it,  ");")
}

findEGs = function(dbcon, symbols) {
  rs = dbSendQuery(dbcon, queryGeneinfo(symbols))
  a1 = fetch(rs, n=-1)
  stillLeft = setdiff(symbols, a1[,2])

  if( length(stillLeft)>0 ) {
    rs = dbSendQuery(dbcon, queryAlias(stillLeft))
    a2 = fetch(rs, n=-1)
    names(a2) = names(a1)
    a1 = rbind(a1, a2)
  } 
  
  rs = dbSendQuery(dbcon, queryGenes(a1[,1]))
  ans = merge(a1, fetch(rs, n=-1))
  dbClearResult(rs)
  ans
}


###################################################
### chunk: findEGs
###################################################
findEGs(mycon, c("ALL1", "AF4", "BCR", "ABL"))


###################################################
### chunk: combineGO_hgu95av2
###################################################
GOdbloc = system.file("extdata", "GO.sqlite", package="GO.db")
attachSql = paste("ATTACH '", GOdbloc, "' as go;", sep = "")
dbGetQuery(mycon, attachSql)


###################################################
### chunk: makeQuery
###################################################
sql = paste("SELECT DISTINCT a.go_id AS 'hgu95av2.go_id',",
            "a._id AS 'hgu95av2._id',",
            "g.go_id AS 'GO.go_id', g._id AS 'GO._id',",
            "g.ontology",
            "FROM go_bp_all AS a, go.go_term AS g", 
            "WHERE a.go_id = g.go_id LIMIT 10;")
dataOut = dbGetQuery(mycon, sql)
dataOut


###################################################
### chunk: dbschemaShow eval=FALSE
###################################################
## schema = capture.output(hgu95av2_dbschema())
## head(schema, 18)


###################################################
### chunk: QAlisting
###################################################
qcdata = capture.output(hgu95av2())
head(qcdata, 20)


###################################################
### chunk: mapcounts eval=FALSE
###################################################
## hgu95av2MAPCOUNTS
## hgu95av2_dbInfo()


###################################################
### chunk: getIntermediateDB
###################################################
tryCatch(library("human.db0"), error=function(e) {
      source("http://bioconductor.org/biocLite.R")
      biocLite("human.db0")
      library("human.db0") } )


###################################################
### chunk: gidb eval=FALSE
###################################################
##   source("http://bioconductor.org/biocLite.R")
##   biocLite("human.db0")


###################################################
### chunk: makeSqlite
###################################################
hgu95av2_IDs = system.file("extdata", 
                           "hgu95av2_ID", 
                           package="AnnotationDbi")

#Then specify some of the metadata for my database
myMeta = c("DBSCHEMA" = "HUMANCHIP_DB",
    "ORGANISM" = "Homo sapiens",
    "SPECIES" = "Human",
    "MANUFACTURER" = "Affymetrix",
    "CHIPNAME" = "Affymetrix Human Genome U95 Set Version 2",
    "MANUFACTURERURL" = "http:www.affymetrix.com")



###################################################
### chunk: td
###################################################
tmpout = tempdir()
popHUMANCHIPDB(affy = FALSE, prefix = "hgu95av2Test",
    fileName = hgu95av2_IDs, metaDataSrc = myMeta,
    baseMapType = "gb", outputDir = tmpout,
    printSchema = TRUE)


###################################################
### chunk: makeAnnDbPkg
###################################################
seed <- new("AnnDbPkgSeed",
            Package = "hgu95av2Test.db",
            Version = "1.0.0",
            PkgTemplate = "HUMANCHIP.DB",
            AnnObjPrefix = "hgu95av2Test")

makeAnnDbPkg(seed, 
             file.path(tmpout, "hgu95av2Test.sqlite"),
             dest_dir = tmpout)


###################################################
### chunk: SQLForge
###################################################
makeHUMANCHIP_DB(affy=FALSE,
    prefix="hgu95av2",
    fileName=hgu95av2_IDs,
    baseMapType="gb",
    outputDir = tmpout,
    version="2.1.0",
    manufacturer = "Affymetrix",
    chipName = "Affymetrix Human Genome U95 Set Version 2",
    manufacturerUrl = "http://www.affymetrix.com")


###################################################
### chunk: setupFileName
###################################################
Yeastfn = system.file("extdata", "yeast_small-01.xml", package="RBioinf")


###################################################
### chunk: checkNamespace
###################################################
yeastIntAct = xmlTreeParse(Yeastfn)
nsY = xmlNamespaceDefinitions(xmlRoot(yeastIntAct))
ns = getDefaultNamespace(xmlRoot(yeastIntAct))
namespaces = c(ns = ns)


###################################################
### chunk: numNspc
###################################################
length(nsY)
sapply(nsY, function(x) x[[2]])


###################################################
### chunk: DOM2
###################################################
nullf = function(x, ...) NULL
yeast2 = xmlTreeParse(Yeastfn, 
    handlers = list(sequence = nullf,
    organism = nullf, primaryRef = nullf, 
    secondaryRef = nullf,
    names = nullf), asTree=TRUE)


###################################################
### chunk: DOM1
###################################################
object.size(yeastIntAct)
object.size(yeast2)


###################################################
### chunk: DOM3
###################################################
yeast3 = xmlTreeParse(Yeastfn, useInternalNodes=TRUE)
f1 = getNodeSet(yeast3, "//ns:attributeList", namespaces)
length(f1)


###################################################
### chunk: getNodeSet
###################################################
iaM = getNodeSet(yeast3, 
    "//ns:interactionDetectionMethod//ns:fullName", 
    namespaces)
sapply(iaM, xmlValue)

f4 = getNodeSet(yeast3, "//ns:hostOrganism//ns:fullName", 
    namespaces)
sapply(f4, xmlValue)


###################################################
### chunk: interactors
###################################################
interactors = getNodeSet(yeast3, 
   "//ns:interactorList//ns:interactor",
   namespaces)
length(interactors)
interactions = getNodeSet(yeast3, 
   "//ns:interactionList/ns:interaction",
   namespaces)
length(interactions)


###################################################
### chunk: xpathApply
###################################################
 interactors = xpathApply(yeast3,
     "//ns:interactorList//ns:interactor", 
     xmlValue, namespaces = namespaces)


###################################################
### chunk: simpleEvent
###################################################
entSH = function(name, attrs, ...) {
          cat("Starting", name, "\n")
          level <<- attrs["level"]
          minorVersion <<- attrs["minorVersion"]
    }
e2 = new.env()
e2$level = NULL
e2$minorVersion = NULL
environment(entSH) = e2


###################################################
### chunk: twomoreHandlers
###################################################
hOrg = function(name, attrs, ...) {
         taxid <<- c(attrs["ncbiTaxId"], taxid)
		      }
e3 = new.env()
e3$taxid = NULL
environment(hOrg) = e3

hInt = function(name, attrs, ...)
		     numInt <<- numInt + 1

e3$numInt = 0
environment(hInt) = e3


###################################################
### chunk: xmlEventP
###################################################
s1 = xmlEventParse(Yeastfn,
	handlers = list(entrySet = entSH, hostOrganism = hOrg,
		         interactor = hInt))

environment(s1$entrySet)$level
environment(s1$hostOrganism)$taxid
environment(s1$interactor)$numInt


###################################################
### chunk: HTMLparsing
###################################################
url = paste("http://www.bioconductor.org/checkResults/", 
    "2.1/bioc-LATEST/", sep="")
s1 = htmlTreeParse(url, useInternalNodes=TRUE)
class(s1)


###################################################
### chunk: FindNodes
###################################################
f1 = getNodeSet(s1, "//a[@href]")
length(f1)


###################################################
### chunk: RefinedGetNodes
###################################################
f2 = getNodeSet(s1, "//b/a[@href]")
p2 = sapply(f2, xmlValue)
length(p2)
p2[1:10]


###################################################
### chunk: 
###################################################
pkgs = sapply(f1, xmlGetAttr, "href")
pkg = grep("/packages/2.1/bioc/html/", pkgs, fixed=TRUE)



###################################################
### chunk: Entrezex
###################################################
 ezURL = "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/"
 t1 = url(ezURL, open="r")
 if( isOpen(t1) ) {
    z = xmlTreeParse(paste(ezURL, "einfo.fcgi", sep=""),
     isURL=TRUE, handlers=NULL, asTree=TRUE) 

    dbL = xmlChildren(z[[1]]$children$eInfoResult)$DbList

    dbNames = xmlSApply(dbL, xmlValue)
    
    length(dbNames)

    dbNames[1:5]
} 


###################################################
### chunk: listMarts
###################################################
library("biomaRt")
head(listMarts())


###################################################
### chunk: selectMart
###################################################
ensM = useMart("ensembl")

ensData = head(listDatasets(ensM))
dim(ensData)

ensMH = useDataset("hsapiens_gene_ensembl", mart=ensM)


###################################################
### chunk: selM eval=FALSE
###################################################
## ensMH = useMart("ensembl", 
##     dataset = "hsapiens_gene_ensembl")


###################################################
### chunk: filters
###################################################
filterSummary(ensMH)
lfilt = listFilters(ensMH, group="GENE:")
nrow(lfilt)
head(lfilt)


###################################################
### chunk: attributes
###################################################
head(attributeSummary(ensMH))
lattr = listAttributes(ensMH, group="PROTEIN:")
lattr


###################################################
### chunk: bioMIntro
###################################################
entrezID = c("983", "3581", "1017") 
rval = getGene(id=entrezID, type="entrezgene", mart = ensMH) 
unique(rval$hgnc_symbol) 


###################################################
### chunk: INTERPRO
###################################################
ensembl = useMart("ensembl", 
    dataset = "hsapiens_gene_ensembl")
    ipro = getBM(attributes=c("entrezgene","interpro",
    "interpro_description"), 
filters = "entrezgene", values = entrezID, 
    mart = ensembl) 
ipro


###################################################
### chunk: GEOqueryEx
###################################################

 library(GEOquery) 

 gds = getGEO("GDS10") 

 eset = GDS2eSet(gds, do.log2 = TRUE) 



###################################################
### chunk: GEOqueryPrintExSetShow eval=FALSE
###################################################
## s1 = experimentData(eset)
## abstract(s1)
## s1@pubMedIds


###################################################
### chunk: KEGGEx
###################################################
library("KEGG")
library("KEGGSOAP")
KEGGPATHID2NAME$"00740"
SoapAns = get.genes.by.pathway("path:sce00740")
SoapAns


###################################################
### chunk: LocalKEGG
###################################################
SA = gsub("^sce:", "", SoapAns)
localAns = KEGGPATHID2EXTID$"sce00740"
setdiff(SA, localAns)