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

###################################################
### code chunk number 1: wavCluster_vignette.Rnw:68-72
###################################################
library(wavClusteR)
filename <- system.file( "extdata", "example.bam", package = "wavClusteR" )
Bam <- readSortedBam(filename = filename)
Bam


###################################################
### code chunk number 2: wavCluster_vignette.Rnw:82-84
###################################################
countTable <- getAllSub( Bam, minCov = 10 )
head( countTable )


###################################################
### code chunk number 3: wavCluster_vignette.Rnw:91-92
###################################################
plotSubstitutions( countTable, highlight = "TC" )


###################################################
### code chunk number 4: wavCluster_vignette.Rnw:107-108 (eval = FALSE)
###################################################
## model <- fitMixtureModel(countTable, substitution = "TC")	#not run


###################################################
### code chunk number 5: wavCluster_vignette.Rnw:112-114
###################################################
data(model)
str(model)


###################################################
### code chunk number 6: wavCluster_vignette.Rnw:124-125
###################################################
(support <- getExpInterval( model, bayes = TRUE ) )


###################################################
### code chunk number 7: wavCluster_vignette.Rnw:132-133
###################################################
(support <- getExpInterval( model, bayes = FALSE, leftProb = 0.9, rightProb = 0.9 ) )


###################################################
### code chunk number 8: wavCluster_vignette.Rnw:141-142
###################################################
plotSubstitutions( countTable, highlight = "TC", model )


###################################################
### code chunk number 9: wavCluster_vignette.Rnw:156-159
###################################################
highConfSub <- getHighConfSub( countTable, 
                               support = support, 
                               substitution = "TC" )


###################################################
### code chunk number 10: wavCluster_vignette.Rnw:162-166
###################################################
highConfSub <- getHighConfSub( countTable, 
                               supportStart = 0.2, 
                               supportEnd = 0.7, 
                               substitution = "TC" )


###################################################
### code chunk number 11: wavCluster_vignette.Rnw:169-170
###################################################
head( highConfSub )


###################################################
### code chunk number 12: wavCluster_vignette.Rnw:179-181
###################################################
coverage <- coverage( Bam )
coverage$chrX


###################################################
### code chunk number 13: wavCluster_vignette.Rnw:186-193
###################################################
clusters <- getClusters( highConfSub = highConfSub,
                         coverage = coverage,
                         sortedBam = Bam,
                         method = "mrn",
                         threshold = 1,
                         cores = 1 )
clusters


###################################################
### code chunk number 14: wavCluster_vignette.Rnw:199-205 (eval = FALSE)
###################################################
## clusters <- getClusters( highConfSub = highConfSub,
##                          coverage = coverage,
##                          sortedBam = Bam,
##                          method = "mrn",
##                          cores = 1 )
## clusters


###################################################
### code chunk number 15: wavCluster_vignette.Rnw:219-228
###################################################
require(BSgenome.Hsapiens.UCSC.hg19)

wavclusters <- filterClusters( clusters = clusters, 
                   highConfSub = highConfSub,
                   coverage = coverage, 
                   model = model, 
                   genome = Hsapiens, 
                   refBase = "T", 
                   minWidth = 12)


###################################################
### code chunk number 16: wavCluster_vignette.Rnw:231-232
###################################################
wavclusters


###################################################
### code chunk number 17: wavCluster_vignette.Rnw:282-286 (eval = FALSE)
###################################################
## exportHighConfSub( highConfSub = highConfSub, 
##                    filename = "hcTC.bed", 
##                    trackname = "hcTC", 
##                    description = "hcTC" )


###################################################
### code chunk number 18: wavCluster_vignette.Rnw:291-295 (eval = FALSE)
###################################################
## exportClusters( clusters = wavclusters, 
##                 filename = "wavClusters.bed", 
##                 trackname = "wavClusters", 
##                 description = "wavClusters" )


###################################################
### code chunk number 19: wavCluster_vignette.Rnw:299-300 (eval = FALSE)
###################################################
## exportCoverage( coverage = coverage, filename = "coverage.bigWig" )


###################################################
### code chunk number 20: wavCluster_vignette.Rnw:308-309 (eval = FALSE)
###################################################
## txDB <- makeTranscriptDbFromUCSC(genome = "hg19", tablename = "ensGene")


###################################################
### code chunk number 21: wavCluster_vignette.Rnw:312-316 (eval = FALSE)
###################################################
## annotateClusters( clusters = wavclusters, 
##               txDB = txDB, 
##               plot = TRUE, 
##               verbose = TRUE)


###################################################
### code chunk number 22: wavCluster_vignette.Rnw:327-336 (eval = FALSE)
###################################################
## getMetaGene( clusters = wavclusters, 
##              txDB = txDB, 
##              upstream = 1e3, 
##              downstream = 1e3, 
##              nBins = 40, 
##              nBinsUD = 10, 
##              minLength = 1, 
##              plot = TRUE,
##              verbose = TRUE ) 


###################################################
### code chunk number 23: wavCluster_vignette.Rnw:343-351 (eval = FALSE)
###################################################
## getMetaTSS( sortedBam = Bam, 
##             txDB = txDB, 
##             upstream = 1e3, 
##             downstream = 1e3, 
##             nBins = 40, 
##             unique = FALSE, 
##             plot = TRUE, 
##             verbose = TRUE ) 


###################################################
### code chunk number 24: wavCluster_vignette.Rnw:363-364
###################################################
plotSizeDistribution( clusters = wavclusters, col = "skyblue2" )


###################################################
### code chunk number 25: wavCluster_vignette.Rnw:371-372
###################################################
plotSizeDistribution( clusters = wavclusters, showCov = TRUE, col = "skyblue2" )


###################################################
### code chunk number 26: session
###################################################
sessionInfo()