## ----knitr, echo=FALSE, results="hide"-----------------------------------
library("knitr")
opts_chunk$set(tidy=FALSE,dev="png",fig.show="hide",
               fig.width=6.5,fig.height=5.5,fig.keep="high",
               message=FALSE)

## ----getVers, echo=FALSE, results="hide"---------------------------------
#Get the versions, if available. otherwise fall back to the static version.

suppressWarnings({
  QoRTsStaticVer = "1.0.1" #REPLACE_THIS_LINE_WITH_QORTS_VERSION_NUMBER
  QoRTsVer <- if( requireNamespace("QoRTs", quietly=TRUE) ){ packageVersion("QoRTs") } else { QoRTsStaticVer }
  JSStaticVer = "1.0.0" #REPLACE_THIS_LINE_WITH_JUNCTIONSEQ_VERSION_NUMBER
  JSVer <- if( requireNamespace("JunctionSeq", quietly=TRUE) ){ packageVersion("JunctionSeq") } else { JSStaticVer }
})


## ----setup,echo=FALSE,results="hide"-------------------------------------
options(width=80, signif=3, digits=3, prompt=" ", continue=" ")
set.seed(0xdada)
#suppressWarnings({
#      library("JunctionSeq")
#})

## ----decoder-------------------------------------------------------------
#Read the decoder:
decoder.file <- system.file("extdata/annoFiles/decoder.bySample.txt",
                             package="JctSeqData");
decoder <- read.table(decoder.file,
                      header=TRUE,
                      stringsAsFactors=FALSE);
print(decoder);

## ----gffdata-------------------------------------------------------------
#The original gtf file, from Ensembl 
#  (all but a few genes were removed, to save space):
anno.original.gtf.file <- system.file("extdata/annoFiles/anno-original.gtf.gz",
                        package="JctSeqData");
head(read.table(anno.original.gtf.file,sep='\t'));

#Modified gtf file, with a few genes 
#  changed to make them into a better test dataset:
anno.gtf.file <- system.file("extdata/annoFiles/anno.gtf.gz",
                        package="JctSeqData");
head(read.table(anno.gtf.file,sep='\t'));

#"Flattened" gff file, for DEXSeq:
DEX.gff.file <- system.file("extdata/annoFiles/DEXSeq.flat.gff.gz",
                        package="JctSeqData");
head(read.table(DEX.gff.file,sep='\t'));

#"Flattened" gff file, for JunctionSeq
#  (w/o novel splice junctions):
JS.gff.file <- system.file("extdata/annoFiles/JunctionSeq.flat.gff.gz",
                        package="JctSeqData");
head(read.table(JS.gff.file,sep='\t'));

#"Flattened" gff file, for JunctionSeq
#  (w/ novel splice junctions):
JS.novel.gff.file <- 
    system.file("extdata/annoFiles/withNovel.forJunctionSeq.gff.gz",
                        package="JctSeqData");
head(read.table(JS.novel.gff.file,sep='\t'));

## ----otherAnno-----------------------------------------------------------
#rn6 chrom.sizes file, from UCSC:
rn6.chrom.sizes <- system.file("extdata/annoFiles/rn6.chrom.sizes",
                        package="JctSeqData");
head(read.table(rn6.chrom.sizes,sep='\t'));


#stripped-down chrom.sizes file
# only includes chr14, chrX, chrM, and loose chr14 contigs:
chrom.sizes <- system.file("extdata/annoFiles/chrom.sizes",
                        package="JctSeqData");
head(read.table(chrom.sizes,sep='\t'));

#mapping of ensembl id to gene symbol:
ensid.2.symbol.file <- system.file("extdata/annoFiles/ensid.2.symbol.txt",
                        package="JctSeqData");
head(read.table(ensid.2.symbol.file,sep='\t',header=TRUE));

## ----DESeqCounts---------------------------------------------------------
#Recall the decoder:
print(decoder);

#Gene level counts:
gene.count.files <- system.file(paste0("extdata/cts/",
                                       decoder$sample.ID,
                                       "/QC.geneCounts.formatted.for.DESeq.txt.gz"
                                       ),
                        package="JctSeqData");

#One of the count files:
read.table(gene.count.files[1])[1:10,]

## ----DEXSeqCounts--------------------------------------------------------
#Exon level counts:
exon.count.files <- system.file(paste0("extdata/cts/",
                                       decoder$sample.ID,
                                       "/QC.exonCounts.formatted.for.DEXSeq.txt.gz"
                                       ),
                        package="JctSeqData");

#Part of One of the count files:
read.table(exon.count.files[1])[110:130,]

## ----JSCounts------------------------------------------------------------
#JunctionSeq counts:
JS.count.files <- system.file(
    paste0("extdata/cts/",
           decoder$sample.ID,
           "/QC.spliceJunctionAndExonCounts.forJunctionSeq.txt.gz"
           ),
     package="JctSeqData");

#Part of one of the count files:
read.table(JS.count.files[1])[526:552,]

## ----JSCounts2-----------------------------------------------------------
#JunctionSeq counts:
JS.novel.count.files <- system.file(
    paste0("extdata/cts/",
           decoder$sample.ID,
           "/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz"
           ),
     package="JctSeqData");

#Part of one of the count files:
read.table(JS.novel.count.files[1])[526:553,]

## ----loadFullData--------------------------------------------------------
data(fullExampleDataSet,package="JctSeqData");

## ----loadTinyData--------------------------------------------------------
data(exampleDataSet,package="JctSeqData");

## ----loadExampleDataSetForSerialization----------------------------------
#Read the decoder:
decoder.file <- system.file("extdata/annoFiles/decoder.bySample.txt",
                             package="JctSeqData");
decoder <- read.table(decoder.file,
                      header=TRUE,
                      stringsAsFactors=FALSE);
#Here are the full-size count and gff files:
gff.file.FULL <- system.file("extdata/cts/withNovel.forJunctionSeq.gff.gz",
                        package="JctSeqData");
countFiles.FULL <- system.file(paste0("extdata/cts/",
      decoder$sample.ID,
      "/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz"),
      package="JctSeqData");
#Here are the "tiny" subset count and gff files:
gff.file.TINY <- system.file("extdata/tiny/withNovel.forJunctionSeq.gff.gz",
                        package="JctSeqData");
countFiles.TINY <- system.file(paste0("extdata/tiny/",
      decoder$sample.ID,
      "/QC.spliceJunctionAndExonCounts.withNovel.forJunctionSeq.txt.gz"),
      package="JctSeqData");

## ----generateFullSerializedExData, eval = FALSE--------------------------
## jscs2 <- runJunctionSeqAnalyses(sample.files = countFiles,
##            sample.names = decoder$sample.ID,
##            condition=factor(decoder$group.ID),
##            flat.gff.file = gff.file);

## ----saveSerializedExData1, eval = FALSE---------------------------------
## save(jscs2, file = "outputData/JctSeqData/data/fullExampleDataSet.RData");

## ----generateTinySerializedExData, eval = FALSE--------------------------
## jscs <- runJunctionSeqAnalyses(sample.files = countFiles.TINY,
##            sample.names = decoder$sample.ID,
##            condition=factor(decoder$group.ID),
##            flat.gff.file = gff.file.TINY);

## ----saveSerializedExData2, eval = FALSE---------------------------------
## save(jscs, file = "outputData/JctSeqData/data/tinyExampleDataSet.RData");