###################################################
### chunk number 1: lkc
###################################################
#line 61 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
library(encoDnaseI)
data(rawCD4)
rawCD4


###################################################
### chunk number 2: lkd
###################################################
#line 69 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
c19g = rawCD4[ chrnum(19) ]
c19g


###################################################
### chunk number 3: lklk
###################################################
#line 75 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
c19gxy = getTrkXY(c19g)
plot(c19gxy)


###################################################
### chunk number 4: doinf
###################################################
#line 86 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
clipSnps = function( sms, chrn, lo, hi ) {
 allp = getSnpLocs(sms)
 allp = allp-allp[1] # relative
 ok = allp >= lo & allp <= hi
 thesm = smList(sms)[[1]]
 rsn = colnames(thesm)
 rid = rsn[which(ok)]
 thesm = thesm[, rid, drop=FALSE]
 nn = new.env()
 tmp = list(thesm)
 names(tmp) = as.character(chrn)
 assign("smList", tmp, nn)
 sms@smlEnv = nn
 sms@activeSnpInds=which(ok)
 sms
}

rangeX = function(htrk) {
 range(getTrkXY(htrk)$x)
}


###################################################
### chunk number 5: dogg
###################################################
#line 111 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
library(GGtools)
library(GGdata)
h19 = getSS("GGdata", "19")
rs19g = rangeX(c19g)
library(SNPlocs.Hsapiens.dbSNP.20090506)
c19l = getSNPlocs("chr19")
h19locs = rbind(rsid=as.numeric(c19l[,"RefSNP_id"]), loc=as.numeric(c19l[,"loc"]))
#h19locs = getSnpLocs(hmceuB36[chrnum(19),])[[1]]
goodlocs = which(h19locs[2,] >= rs19g[1] & h19locs[2,] <= rs19g[2])
h19rsn = paste("rs", h19locs[1,goodlocs], sep="")
h19trim = h19[rsid(h19rsn),]
#ok
#c19gf = clipSnps( hmceuB36[chrnum(19),], chrnum(19), rs19g[1], rs19g[2] )
#c19gf


###################################################
### chunk number 6: lkmxi1
###################################################
#line 129 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
oo = options()  # don't take warnings on multiple probes... caveat emptor
options(warn=0)
library(GGtools)
showMethods("gwSnpTests")
smxi1 = gwSnpTests(genesym("MXI1")~1-1, h19trim, chrnum(19))
smxi1
#plot(smxi1)
options(oo)


###################################################
### chunk number 7: doj
###################################################
#line 144 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
print(juxtaPlot( c19g, smxi1, h19locs ))


###################################################
### chunk number 8: donex
###################################################
#line 150 "vignettes/encoDnaseI/inst/doc/dnaseUse.Rnw"
oo = options()
options(warn=0)
sOSR2 = gwSnpTests(genesym("OSR2")~1-1, h19trim, chrnum(19))
print(juxtaPlot( c19g, sOSR2, h19locs ))
options(oo)