#v2.1
/data/results/tools/r/R-4.0.0/bin/R
library(metaseqR2)

#@hg19
- CNTR untreated vs CNTR treated
- CNTR untreated vs Cre treated
- CNTR treated vs Cre treated

the.path <- "/home/reczko/bak/doc/fleming/kolliarakis/VKlab/run514"
the.contrasts.1 <- c(
	"Ctrl_tr_vs_Ctrl_untr",
	"Cre_tr_vs_Ctrl_untr",
	"Cre_tr_vs_Ctrl_tr"	
)

buildDir <- file.path(tempdir(),"test_anndb")
dir.create(buildDir)

# The location of the custom database
myDb <- file.path(buildDir,"testann.sqlite")

# Since we are using Ensembl, we can also ask for a version
organisms <- list(hg19=105)
sources <- ifelse(.Platform$OS.type=="unix",c("ensembl","refseq"),"ensembl")

# If the example is not running in a multicore system, rc is ignored
buildAnnotationDatabase(organisms,sources,forceDownload=FALSE,db=myDb,rc=0.5)


# Opening metaseqR2 SQLite database /tmp/Rtmp8hyZU4/test_anndb/testann.sqlite

result <- metaseqr2(
    sampleList=file.path(the.path, "targets.txt"),
    fileType = "bam",
    contrast=the.contrasts.1,
#    org="mm10",
#    annotation="download",
    org="hg19",
#annotation="embedded",
#localDb ="/data/results/tools/rnaseq/metaseqr/mm10/annotation.sqlite",
localDb ="/tmp/Rtmp8hyZU4/test_anndb/testann.sqlite",
    refdb = "ensembl",
    transLevel = "gene",
    countType="utr",
    normalization="deseq",
    statistics=c("deseq","deseq2","edger","noiseq","limma","nbpseq","absseq","dss"),
    adjustMethod = "fdr",
    metaP="pandora",
    figFormat=c("png","pdf", "jpg"),
    exportWhere=file.path(the.path, "metaseqR2_run514b"),
    restrictCores=0.5,
    qcPlots=c("mds","biodetection","countsbio","saturation","readnoise",
              "rnacomp","gcbias",
              "pairwise","filtered","correl","boxplot","lengthbias","meandiff",
              "meanvar","boxplot", "filtered", "biodist","volcano","mastat"
              #, "deheatmap"
              ),
    exonFilters=NULL,
    geneFilters=list(
      #length=list(length=500),
      avgReads=list(averagePerBp=100, quantile=0.25),
      expression=list(median=TRUE,
                     mean=FALSE,
                     quantile=NA,
                     known=NA,
                     custom=NA),
      biotype=getDefaults("biotypeFilter", "mm10"),
      presence = list(
        frac=0.5,
        minCount=1,
        perCondition=TRUE
      )
    ),
    pcut=0.05,
    outList = TRUE,
    exportWhat=c("annotation","p_value","adj_p_value","meta_p_value",
                 "adj_meta_p_value","fold_change","stats","counts","flags"),
    exportScale=c("natural","log2","rpgm"),
    exportValues="normalized",
    exportStats = c("mean", "median", "sd", "mad", "cv"),
    exportCountsTable=TRUE,
    #reportTop=0.05,
    saveGeneModel = TRUE,
    createTracks=TRUE,
    overwrite=TRUE,
    trackInfo=list(
      stranded=TRUE,
      normTo=1e+8,
      #urlBase="http://epigenomics.fleming.gr/home/tzanos/public_html/DK/HuR/metaseqR2_DK_pandora_tracks/tracks",
      hubInfo=list(
        name="VKHub",
        shortLabel="VK Hub",
        longLabel="VK hub long",
        email="reczko@fleming.gr"
      )
    )
  )

#@ exclude VK272-2_CNTR_treated	/home/reczko/bak/doc/fleming/kolliarakis/VKlab/run514/VK272-2_CNTR_treated.bam	Ctrl_tr	single	yes

result <- metaseqr2(
    sampleList=file.path(the.path, "targets2.txt"),
    fileType = "bam",
    contrast=the.contrasts.1,
#    org="mm10",
#    annotation="download",
    org="hg19",
#annotation="embedded",
#localDb ="/data/results/tools/rnaseq/metaseqr/mm10/annotation.sqlite",
localDb ="/tmp/RtmpiZTU6I/test_anndb/testann.sqlite",
    refdb = "ensembl",
    transLevel = "gene",
    countType="utr",
    normalization="deseq",
    statistics=c("deseq","deseq2","edger","noiseq","limma","nbpseq","absseq","dss"),
    adjustMethod = "fdr",
    metaP="pandora",
    figFormat=c("png","pdf", "jpg"),
    exportWhere=file.path(the.path, "metaseqR2_run514c"),
    restrictCores=0.5,
    qcPlots=c("mds","biodetection","countsbio","saturation","readnoise",
              "rnacomp","gcbias",
              "pairwise","filtered","correl","boxplot","lengthbias","meandiff",
              "meanvar","boxplot", "filtered", "biodist","volcano","mastat"
              #, "deheatmap"
              ),
    exonFilters=NULL,
    geneFilters=list(
      #length=list(length=500),
      avgReads=list(averagePerBp=100, quantile=0.25),
      expression=list(median=TRUE,
                     mean=FALSE,
                     quantile=NA,
                     known=NA,
                     custom=NA),
      biotype=getDefaults("biotypeFilter", "mm10"),
      presence = list(
        frac=0.5,
        minCount=1,
        perCondition=TRUE
      )
    ),
    pcut=0.05,
    outList = TRUE,
    exportWhat=c("annotation","p_value","adj_p_value","meta_p_value",
                 "adj_meta_p_value","fold_change","stats","counts","flags"),
    exportScale=c("natural","log2","rpgm"),
    exportValues="normalized",
    exportStats = c("mean", "median", "sd", "mad", "cv"),
    exportCountsTable=TRUE,
    #reportTop=0.05,
    saveGeneModel = TRUE,
    createTracks=TRUE,
    overwrite=TRUE,
    trackInfo=list(
      stranded=TRUE,
      normTo=1e+8,
      #urlBase="http://epigenomics.fleming.gr/home/tzanos/public_html/DK/HuR/metaseqR2_DK_pandora_tracks/tracks",
      hubInfo=list(
        name="VKHub",
        shortLabel="VK Hub",
        longLabel="VK hub long",
        email="reczko@fleming.gr"
      )
    )
  )