### R code from vignette source 'CAGEr.Rnw'

###################################################
### code chunk number 1: setup
###################################################
options(width = 60)
olocale=Sys.setlocale(locale="C")


###################################################
### code chunk number 2: CAGEr.Rnw:137-138
###################################################
library(CAGEr)


###################################################
### code chunk number 3: CAGEr.Rnw:144-145
###################################################
library(BSgenome.Drerio.UCSC.danRer7)


###################################################
### code chunk number 4: CAGEr.Rnw:166-169
###################################################
inputDir <- system.file("extdata", package = "CAGEr")
pathsToInputFiles <- list.files(inputDir, full.names = TRUE)
basename(pathsToInputFiles)


###################################################
### code chunk number 5: CAGEr.Rnw:177-181
###################################################
myCAGEset <- new("CAGEset", genomeName = "BSgenome.Drerio.UCSC.danRer7",
		inputFiles = pathsToInputFiles, inputFilesType = "ctss",
		sampleLabels = c("zf_30p_dome", "zf_high", 
		"zf_prim6_rep1", "zf_prim6_rep2", "zf_unfertilized_egg"))


###################################################
### code chunk number 6: CAGEr.Rnw:185-186
###################################################
myCAGEset


###################################################
### code chunk number 7: CAGEr.Rnw:208-218
###################################################
TSS.df <- read.table(system.file("extdata/Zf.unfertilized.egg.chr17.ctss", 
                                package = "CAGEr"))
# make sure the column names are as required
colnames(TSS.df) <- c("chr", "pos", "strand", "zf_unfertilized_egg")
# make sure the column classes are as required
TSS.df$chr <- as.character(TSS.df$chr)
TSS.df$pos <- as.integer(TSS.df$pos)
TSS.df$strand <- as.character(TSS.df$strand)
TSS.df$zf_unfertilized_egg <- as.integer(TSS.df$zf_unfertilized_egg)
head(TSS.df)


###################################################
### code chunk number 8: CAGEr.Rnw:222-224
###################################################
myCAGEset.coerced <- as(TSS.df, "CAGEset")
myCAGEset.coerced


###################################################
### code chunk number 9: CAGEr.Rnw:229-230
###################################################
getCTSS(myCAGEset)


###################################################
### code chunk number 10: CAGEr.Rnw:235-237
###################################################
ctss <- CTSStagCount(myCAGEset)
head(ctss)


###################################################
### code chunk number 11: CAGEr.Rnw:241-242
###################################################
sampleLabels(myCAGEset)


###################################################
### code chunk number 12: CAGEr.Rnw:253-254 (eval = FALSE)
###################################################
## corr.m <- plotCorrelation(myCAGEset, samples = "all", method = "pearson")


###################################################
### code chunk number 13: CAGEr.Rnw:273-276
###################################################
mergeSamples(myCAGEset, mergeIndex = c(3,2,4,4,1), 
			mergedSampleLabels = c("zf_unfertilized_egg", 
			"zf_high", "zf_30p_dome", "zf_prim6"))


###################################################
### code chunk number 14: CAGEr.Rnw:288-289
###################################################
librarySizes(myCAGEset)


###################################################
### code chunk number 15: CAGEr.Rnw:301-302 (eval = FALSE)
###################################################
## plotReverseCumulatives(myCAGEset, fitInRange = c(5, 1000), onePlot = TRUE)


###################################################
### code chunk number 16: CAGEr.Rnw:320-322
###################################################
normalizeTagCount(myCAGEset, method = "powerLaw", 
		fitInRange = c(5, 1000), alpha = 1.2, T = 5*10^4)


###################################################
### code chunk number 17: CAGEr.Rnw:341-342 (eval = FALSE)
###################################################
## exportCTSStoBedGraph(myCAGEset, values = "normalized", oneFile = TRUE)


###################################################
### code chunk number 18: CAGEr.Rnw:373-376
###################################################
clusterCTSS(object = myCAGEset, threshold = 1, thresholdIsTpm = TRUE, 
		nrPassThreshold = 1, method = "distclu", maxDist = 20, 
		removeSingletons = TRUE, keepSingletonsAbove = 5)


###################################################
### code chunk number 19: CAGEr.Rnw:382-384
###################################################
tc <- tagClusters(myCAGEset, sample = "zf_unfertilized_egg")
head(tc)


###################################################
### code chunk number 20: CAGEr.Rnw:421-423
###################################################
cumulativeCTSSdistribution(myCAGEset, clusters = "tagClusters")
quantilePositions(myCAGEset, clusters = "tagClusters", qLow = 0.1, qUp = 0.9)


###################################################
### code chunk number 21: CAGEr.Rnw:427-430
###################################################
tc <- tagClusters(myCAGEset, sample = "zf_unfertilized_egg",
		returnInterquantileWidth = TRUE,  qLow = 0.1, qUp = 0.9)
head(tc)


###################################################
### code chunk number 22: CAGEr.Rnw:434-436 (eval = FALSE)
###################################################
## exportToBed(object = myCAGEset, what = "tagClusters", 
## 		qLow = 0.1, qUp = 0.9, oneFile = TRUE)


###################################################
### code chunk number 23: CAGEr.Rnw:450-452 (eval = FALSE)
###################################################
## plotInterquantileWidth(myCAGEset, clusters = "tagClusters", 
## 		tpmThreshold = 3, qLow = 0.1, qUp = 0.9)


###################################################
### code chunk number 24: CAGEr.Rnw:470-472
###################################################
aggregateTagClusters(myCAGEset, tpmThreshold = 5, 
		qLow = 0.1, qUp = 0.9, maxDist = 100)


###################################################
### code chunk number 25: CAGEr.Rnw:477-479
###################################################
consensusCl <- consensusClusters(myCAGEset)
head(consensusCl)


###################################################
### code chunk number 26: CAGEr.Rnw:494-496
###################################################
getExpressionProfiles(myCAGEset, what = "consensusClusters", tpmThreshold = 10, 
		nrPassThreshold = 1, method = "som", xDim = 4, yDim = 2)


###################################################
### code chunk number 27: CAGEr.Rnw:500-501 (eval = FALSE)
###################################################
## plotExpressionProfiles(myCAGEset, what = "consensusClusters")


###################################################
### code chunk number 28: CAGEr.Rnw:513-516
###################################################
class3_1 <- extractExpressionClass(myCAGEset, 
		what = "consensusClusters", which = "3_1")
head(class3_1)


###################################################
### code chunk number 29: CAGEr.Rnw:521-523 (eval = FALSE)
###################################################
## exportToBed(myCAGEset, what = "consensusClusters", 
## 		colorByExpressionProfile = TRUE)


###################################################
### code chunk number 30: CAGEr.Rnw:545-546
###################################################
cumulativeCTSSdistribution(myCAGEset, clusters = "consensusClusters")


###################################################
### code chunk number 31: CAGEr.Rnw:549-551
###################################################
scoreShift(myCAGEset, groupX = "zf_unfertilized_egg", groupY = "zf_prim6",
		testKS = TRUE, useTpmKS = FALSE)


###################################################
### code chunk number 32: CAGEr.Rnw:575-579
###################################################
shifting.promoters <- getShiftingPromoters(myCAGEset, 
		tpmThreshold = 5, scoreThreshold = 0.6,
		fdrThreshold = 0.01)
head(shifting.promoters)


###################################################
### code chunk number 33: CAGEr.Rnw:614-617
###################################################
data(FANTOM5humanSamples)
head(FANTOM5humanSamples)
nrow(FANTOM5humanSamples)


###################################################
### code chunk number 34: CAGEr.Rnw:635-638
###################################################
astrocyteSamples <- FANTOM5humanSamples[grep("Astrocyte", 
				FANTOM5humanSamples[,"description"]),]
astrocyteSamples


###################################################
### code chunk number 35: CAGEr.Rnw:641-644
###################################################
data(FANTOM5mouseSamples)
head(FANTOM5mouseSamples)
nrow(FANTOM5mouseSamples)


###################################################
### code chunk number 36: CAGEr.Rnw:648-649
###################################################
astrocyteSamples[,"sample"]


###################################################
### code chunk number 37: CAGEr.Rnw:652-655
###################################################
astrocyteCAGEset <- importPublicData(source = "FANTOM5", dataset = "human", 
					sample = astrocyteSamples[1:3,"sample"])
astrocyteCAGEset


###################################################
### code chunk number 38: CAGEr.Rnw:663-665 (eval = FALSE)
###################################################
## source("http://bioconductor.org/biocLite.R")
## biocLite("FANTOM3and4CAGE")


###################################################
### code chunk number 39: CAGEr.Rnw:668-674
###################################################
library(FANTOM3and4CAGE)

data(FANTOMhumanSamples)
head(FANTOMhumanSamples)
data(FANTOMmouseSamples)
head(FANTOMmouseSamples)


###################################################
### code chunk number 40: CAGEr.Rnw:678-682
###################################################
kidneyCAGEset <- importPublicData(source = "FANTOM3and4",
			dataset = "FANTOMtissueCAGEhuman", 
			group = "kidney", sample = c("kidney", "malignancy"))
kidneyCAGEset


###################################################
### code chunk number 41: CAGEr.Rnw:685-691
###################################################
mixedCAGEset <- importPublicData(source = "FANTOM3and4",
		dataset = c("FANTOMtissueCAGEmouse", "FANTOMtissueCAGEmouse", 
		"FANTOMtimecourseCAGEmouse"), group = c("liver", "liver", 
		"liver_under_constant_darkness"), 
		sample = c("cloned_mouse", "control_mouse", "4_hr"))
mixedCAGEset


###################################################
### code chunk number 42: CAGEr.Rnw:700-703 (eval = FALSE)
###################################################
## library(ENCODEprojectCAGE)
## 
## data(ENCODEhumanCellLinesSamples)


###################################################
### code chunk number 43: CAGEr.Rnw:707-711 (eval = FALSE)
###################################################
## ENCODEset <- importPublicData(source = "ENCODE", 
## 	dataset = c("A549", "H1-hESC", "IMR90"), 
## 	group = c("cell", "cell", "cell"), sample = c("A549_cell_rep1", 
## 	"H1-hESC_cell_rep1", "IMR90_cell_rep1"))


###################################################
### code chunk number 44: CAGEr.Rnw:720-723 (eval = FALSE)
###################################################
## library(ZebrafishDevelopmentalCAGE)
## 
## data(ZebrafishSamples)


###################################################
### code chunk number 45: CAGEr.Rnw:727-730 (eval = FALSE)
###################################################
## zebrafishCAGEset <- importPublicData(source = "ZebrafishDevelopment", 
## 			dataset = "ZebrafishCAGE", group = "development", 
## 			sample = c("zf_64cells", "zf_prim6"))


###################################################
### code chunk number 46: CAGEr.Rnw:742-743
###################################################
sessionInfo()