##################################################################################
#Andy Rampersaud, 07.19.17
#This 01_Pipeline_Setup.sh file is used to initialize all variables used in the pipeline
#The purpose of this file is to streamline the pipeline so that 
#	1. Variables are initialized from a single file
#	2. Minimize redundancy between pipeline steps
#This file is organized by 
#	1. "Global Variables"  		= variables used by multiple steps
#	2. "Step-specific Variables"  	= variables used by a particular step
##################################################################################
#---------------------------------------------------------------------------------
#"Global Variables"  = variables used by multiple steps
#---------------------------------------------------------------------------------
Dataset_DIR=/restricted/projectnb/waxmanlab/kkarri/tabula_muris/smartseq_fastq/test_samples
BU_User="kkarri"
Dataset_Label="SRR"
#Need this dir that contains the various GTF files
GTF_Files_DIR=/restricted/projectnb/waxmanlab/kkarri/tabula_muris/GTF_Files
#---------------------------------------------------------------------------------
##################################################################################
#---------------------------------------------------------------------------------
#"Step-specific Variables"  = variables used by a particular step
#---------------------------------------------------------------------------------
#TopHat_Paired_End
#Need this dir for the Bowtie2 Index
Bowtie2Index_DIR=/restricted/projectnb/waxmanlab/kkarri/tabula_muris/Bowtie2
#---------------------------------------------------------------------------------
#Set this variable to the appropriate value <0, 1, or 2> 
#Refer to the Parameter_Descriptions.txt in the TopHat_Paired_End folder for more detailed information
#0="fr-unstranded"
#1="fr-firststrand"
#2="fr-secondstrand"
STRANDEDNESS="0"
#---------------------------------------------------------------------------------
#This script will do the arithmetic to calculate your (--mate-inner-dist)
#Input the Bioanalyzer length (bp) (directly from Bioanalyzer tracings)
BIOANALYZER_LEN=277
#Note: Lengths directly from Bioanalyzer tracings still include adaptors (ligated to both ends)
#Usually the adaptor length is 60bp (60bp on each end: 120bp total)
#Input the adaptor length (bp) for one end:
ADAPTOR_LEN=60
#Input the read length (bp)
#Should have this from the Read_Length job
READ_LEN=99	
#---------------------------------------------------------------------------------
#Set this variable to the appropriate value <0, 1, or 2> 
#Refer to the Parameter_Descriptions.txt in the CollectRnaSeqMetrics folder for more detailed information
#0="NONE"
#1="FIRST_READ_TRANSCRIPTION_STRAND"
#2="SECOND_READ_TRANSCRIPTION_STRAND"
STRAND_SPECIFICITY="2"
#---------------------------------------------------------------------------------
#Extract_Counts
#---------------------------------------------------------------------------------
#Set this variable to the appropriate value <0, 1, or 2> 
#Refer to the Parameter_Descriptions.txt in the Extract_Counts folder for more detailed information
#0="union"
#1="intersection-strict"
#2="intersection-nonempty"
MODE="2"
#---------------------------------------------------------------------------------
#Set this variable to the appropriate value <0, 1, or 2> 
#Refer to the Parameter_Descriptions.txt in the Extract_Counts folder for more detailed information
#0="no"
#1="yes"
#2="reverse"
STRANDEDNESS_HTSeq="0"
#---------------------------------------------------------------------------------
#Extract_Counts
#---------------------------------------------------------------------------------
#Set this variable to the appropriate value <0, 1, or 2> 
#Refer to the Parameter_Descriptions.txt in the Extract_Counts folder for more detailed information
#0="unstranded"
#1="stranded"
#2="reversely stranded"
STRANDEDNESS_featureCount="0"
#---------------------------------------------------------------------------------
#09b_DiffExp_1
#Note: 
#	Modify setup_DiffExp.sh (located in job-specific folder) as needed
#	Modify Condition_*.txt files as needed
#---------------------------------------------------------------------------------
#	For users with multiple differential gene expression (DE) comparisons:
#	Each comparison is run in a separate folder; for instance (regardless of the number of replicates per condition):
#	Condition_A vs. Condition_B => calculated in 09*_DiffExp_1
#	Condition_C vs. Condition_D => calculated in 09*_DiffExp_2
#	Condition_E vs. Condition_F => calculated in 09*_DiffExp_3
#	etc ...
#	The user needs to copy a 09*_DiffExp_* folder and use it as a template for any additional DE comparisons 
#	Then the pipeline or individual step can be run as usual
#---------------------------------------------------------------------------------
#Note:
#Three types of DiffExp folders are provided: 
#	(1) HTSeq, (2) featureCounts and (3) lncRNA_featureCounts. 
#	The recommended option is to use featureCounts, and not HTSeq, for counting of RefSeq genes, and to always use lncRNA_featureCounts when counting mouse liver datasets. 
#featureCounts vs. HTSeq: 
#	For two RefSeq genes that show extensive overlap, leaving less than 130 bp of unique sequence, featureCounts counts the overlapping region reads for BOTH RefSeq genes. HTSeq excludes from the count all such overlapping reads, and may thus seriously undercount such a gene. Omit the HTSeq folder unless you have a specific need for that output. 
#	featureCounts, but not HTSeq, counts non-unique reads for any genomic region that is covered by two overlapping RefSeq genes that show extensive overlap . The recommended option is to  use featureCounts folders only (omit HTSeq folders).
#lncRNAs: 
#	This counting method is based on lncRNA gene structures defined for mouse liver. There is, in general, no need to perform DE analysis for lncRNAs for non-liver samples.
#---------------------------------------------------------------------------------
#bamCorrelate
zMin=0.6
zMax=1.0
#---------------------------------------------------------------------------------
#UCSC_BigWig
#---------------------------------------------------------------------------------
#Set this variable to the appropriate value <0, 1, or 2> 
#Refer to the Parameter_Descriptions.txt in the UCSC_BigWig folder for more detailed information
#Paired-end options:
#0="none"
#1="1++,1--,2+-,2-+"
#2="1+-,1-+,2++,2--"
STRAND_RULE="0"
#---------------------------------------------------------------------------------
#11_Generate_Tracks
#Notes: 
#	Please refer to README_Generate_Tracks.sh (located in job-specific folder)
#	Modify Sample_Labels.txt (located in job-specific folder) as needed
#---------------------------------------------------------------------------------
##################################################################################
#DO NOT EDIT CODE BELOW THIS LINE
##################################################################################
#---------------------------------------------------------------------------------
#Need if statements to process variables with <0, 1, or 2>
#For multiple if statements - use elif
#Need to re-assign these variables using the user defined value and convert them to valid arguments for its respective program
#---------------------------------------------------------------------------------
#TopHat_Paired_End: STRANDEDNESS
if [ ${STRANDEDNESS} -eq 0 ]
then
STRANDEDNESS="fr-unstranded"
elif [ ${STRANDEDNESS} -eq 1 ]
then
STRANDEDNESS="fr-firststrand"
elif [ ${STRANDEDNESS} -eq 2 ]
then
STRANDEDNESS="fr-secondstrand"
fi
#echo ${STRANDEDNESS}
#---------------------------------------------------------------------------------
#CollectRnaSeqMetrics: STRAND_SPECIFICITY
if [ ${STRAND_SPECIFICITY} -eq 0 ]
then
STRAND_SPECIFICITY="NONE"
elif [ ${STRAND_SPECIFICITY} -eq 1 ]
then
STRAND_SPECIFICITY="FIRST_READ_TRANSCRIPTION_STRAND"
elif [ ${STRAND_SPECIFICITY} -eq 2 ]
then
STRAND_SPECIFICITY="SECOND_READ_TRANSCRIPTION_STRAND"
fi
#echo ${STRAND_SPECIFICITY}
#---------------------------------------------------------------------------------
#Extract_Counts: MODE
if [ ${MODE} -eq 0 ]
then
MODE="union"
elif [ ${MODE} -eq 1 ]
then
MODE="intersection-strict"
elif [ ${MODE} -eq 2 ]
then
MODE="intersection-nonempty"
fi
#echo ${MODE}
#---------------------------------------------------------------------------------
#Extract_Counts: STRANDEDNESS_HTSeq
if [ ${STRANDEDNESS_HTSeq} -eq 0 ]
then
STRANDEDNESS_HTSeq="no"
elif [ ${STRANDEDNESS_HTSeq} -eq 1 ]
then
STRANDEDNESS_HTSeq="yes"
elif [ ${STRANDEDNESS_HTSeq} -eq 2 ]
then
STRANDEDNESS_HTSeq="reverse"
fi
#echo ${STRANDEDNESS_HTSeq}
#---------------------------------------------------------------------------------
#Extract_Counts: STRANDEDNESS_featureCount
if [ ${STRANDEDNESS_featureCount} -eq 0 ]
then
STRANDEDNESS_featureCount="0"
elif [ ${STRANDEDNESS_featureCount} -eq 1 ]
then
STRANDEDNESS_featureCount="1"
elif [ ${STRANDEDNESS_featureCount} -eq 2 ]
then
STRANDEDNESS_featureCount="2"
fi
#echo ${STRANDEDNESS_featureCount}
#---------------------------------------------------------------------------------
#UCSC_BigWig: STRAND_RULE
if [ ${STRAND_RULE} -eq 0 ]
then
STRAND_RULE="none"
elif [ ${STRAND_RULE} -eq 1 ]
then
STRAND_RULE="1++,1--,2+-,2-+"
elif [ ${STRAND_RULE} -eq 2 ]
then
STRAND_RULE="1+-,1-+,2++,2--"
fi
#echo ${STRAND_RULE}
#---------------------------------------------------------------------------------
##################################################################################
#---------------------------------------------------------------------------------
#"Additional Variables" = variables not user-specified
#---------------------------------------------------------------------------------
#Extract_Counts
#We always count by gene_id when using HTSeq
FEATURE_ID="gene_id"
#---------------------------------------------------------------------------------
#Location of Sample_Labels.txt (text file indicating the samples to process):
Sample_Labels_DIR=${Dataset_DIR}/Scripts/00_Setup_Pipeline
#---------------------------------------------------------------------------------
#Need to get the current dir
#This variable is used when creating the job name for qsub commands
SCRIPT_DIR=$(pwd)
#---------------------------------------------------------------------------------
#Do the calculation(s):
FRAGMENT_LEN=$(echo "scale=4;$BIOANALYZER_LEN-(2*$ADAPTOR_LEN)" | bc)
DISTANCE_BT_READ_PAIR=$(echo "scale=4;$FRAGMENT_LEN-(2*$READ_LEN)" | bc)
#---------------------------------------------------------------------------------
#For read mapping purposes we want the full gene body information (use the RefSeq_GeneBody.gtf)
#Presence of splice junctions should be based on the full gene structure for any isoform of a gene symbol
ANNOTATION_FILE=RefSeq_GeneBody.gtf
#---------------------------------------------------------------------------------
#Webserver location(s) for hosting files
#Note: We are using the waxman-server mount point on the SCC
VM_DIR_FASTQC=/net/waxman-server/mnt/data/waxmanlabvm_home/waxmanlab/FASTQC/${Dataset_Label}
VM_DIR_UCSC=/net/waxman-server/mnt/data/waxmanlabvm_home/${BU_User}/${Dataset_Label}
#---------------------------------------------------------------------------------
#Time hour limit
#On SCC a 12-hour runtime limit is enforced on all jobs, unless specified explicitly. 
#A runtime limit can be specified in the format "hh:mm:ss"
TIME_LIMIT="96:00:00"
#---------------------------------------------------------------------------------
##################################################################################