Merge branch '47-move.inference' into 'develop'

Resolve "Move inference to start of pipeline"

Closes #47

See merge request !28

.gitlab-ci.yml

@@ -32,9 +32,21 @@ parseMetadata:
   - singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/parseMeta.py -r Replicate_RID -m "./test_data/meta/metaTest.csv" -p repRID
   - singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/parseMeta.py -r Replicate_RID -m "./test_data/meta/metaTest.csv" -p ends
   - singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/parseMeta.py -r Replicate_RID -m "./test_data/meta/metaTest.csv" -p endsManual
-  - singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/parseMeta.py -r Replicate_RID -m "./test_data/meta/metaTest.csv" -p stranded -e se
+  - singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/parseMeta.py -r Replicate_RID -m "./test_data/meta/metaTest.csv" -p stranded
   - singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/parseMeta.py -r Replicate_RID -m "./test_data/meta/metaTest.csv" -p species
 
+inferMetadata:
+  stage: unit
+  script:
+  - >
+    align=$(echo $(grep "Overall alignment rate" ./test_data/meta/Q-Y5JA_1M.se.alignSummary.txt | cut -f2 -d ':' | cut -f2 -d ' ' | tr -d '%')) &&
+    if [[ ${align} == "" ]]; then exit 1; fi
+  - >
+    singularity run 'docker://bicf/rseqc3.0:2.0.0' infer_experiment.py -r "/project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/bed/genome.bed" -i "./test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.bam" 1>> Q-Y5JA_1M.se.inferMetadata.log &&
+    ended=`singularity run 'docker://bicf/python3:1.3' python3 ./workflow/scripts/inferMeta.sh endness Q-Y5JA_1M.se.inferMetadata.log` &&
+    if [[ ${ended} == "" ]]; then exit 1; fi
+  - pytest -m inferMetadata
+
 getRef:
   stage: unit
   script:
@@ -50,14 +62,20 @@ trimData:
   - singularity run 'docker://bicf/trimgalore:1.1' trim_galore --gzip -q 25 --illumina --length 35 --paired --basename Q-Y5JA_1M.pe -j 20 ./test_data/fastq/small/Q-Y5JA_1M.R1.fastq.gz ./test_data/fastq/small/Q-Y5JA_1M.R2.fastq.gz
   - pytest -m trimData
 
+downsampleData:
+  stage: unit
+  script:
+  - singularity exec 'docker://bicf/seqtk:2.0.0' seqtk sample -s100 ./test_data/fastq/small/Q-Y5JA_1M_trimmed.fq.gz 1000 1> sampled.1.fq
+  - pytest -m downsampleData
+
 alignData:
   stage: unit
   script:
-  - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' hisat2 -p 20 --add-chrname --un-gz Q-Y5JA_1M.se.unal.gz -S Q-Y5JA_1M.se.sam -x /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/hisat2/genome --rna-strandness F -U ./test_data/fastq/small/Q-Y5JA_1M_trimmed.fq.gz --summary-file ${repRID}.alignSummary.txt --new-summary
+  - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' hisat2 -p 20 --add-chrname --un-gz Q-Y5JA_1M.se.unal.gz -S Q-Y5JA_1M.se.sam -x /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/hisat2/genome --rna-strandness F -U ./test_data/fastq/small/Q-Y5JA_1M_trimmed.fq.gz --summary-file Q-Y5JA_1M.se.alignSummary.txt --new-summary
   - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' samtools view -1 -@ 20 -F 4 -F 8 -F 256 -o Q-Y5JA_1M.se.bam Q-Y5JA_1M.se.sam
   - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' samtools sort -@ 20 -O BAM -o Q-Y5JA_1M.se.sorted.bam Q-Y5JA_1M.se.bam
   - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' samtools index -@ 20 -b Q-Y5JA_1M.se.sorted.bam Q-Y5JA_1M.se.sorted.bam.bai
-  - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' hisat2 -p 20 --add-chrname --un-gz Q-Y5JA_1M.pe.unal.gz -S Q-Y5JA_1M.pe.sam -x /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/hisat2/genome --rna-strandness FR --no-mixed --no-discordant -1 ./test_data/fastq/small/Q-Y5JA_1M_R1_val_1.fq.gz -2 ./test_data/fastq/small/Q-Y5JA_1M_R2_val_2.fq.gz --summary-file ${repRID}.alignSummary.txt --new-summary
+  - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' hisat2 -p 20 --add-chrname --un-gz Q-Y5JA_1M.pe.unal.gz -S Q-Y5JA_1M.pe.sam -x /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/hisat2/genome --rna-strandness FR --no-mixed --no-discordant -1 ./test_data/fastq/small/Q-Y5JA_1M_R1_val_1.fq.gz -2 ./test_data/fastq/small/Q-Y5JA_1M_R2_val_2.fq.gz --summary-file Q-Y5JA_1M.pe.alignSummary.txt --new-summary
   - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' samtools view -1 -@ 20 -F 4 -F 8 -F 256 -o Q-Y5JA_1M.pe.bam Q-Y5JA_1M.pe.sam
   - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' samtools sort -@ 20 -O BAM -o Q-Y5JA_1M.pe.sorted.bam Q-Y5JA_1M.pe.bam
   - singularity run 'docker://bicf/gudmaprbkaligner:2.0.0' samtools index -@ 20 -b Q-Y5JA_1M.pe.sorted.bam Q-Y5JA_1M.pe.sorted.bam.bai
@@ -69,36 +87,39 @@ dedupData:
   - singularity run 'docker://bicf/gudmaprbkdedup:2.0.0' java -jar /picard/build/libs/picard.jar MarkDuplicates I=./test_data/bam/small/Q-Y5JA_1M.se.sorted.bam O=Q-Y5JA_1M.se.deduped.bam M=Q-Y5JA_1M.se.deduped.Metrics.txt REMOVE_DUPLICATES=true
   - singularity run 'docker://bicf/gudmaprbkdedup:2.0.0' samtools sort -@ 20 -O BAM -o Q-Y5JA_1M.se.sorted.deduped.bam ./test_data/bam/small/Q-Y5JA_1M.se.deduped.bam
   - singularity run 'docker://bicf/gudmaprbkdedup:2.0.0' samtools index -@ 20 -b ./test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.bam Q-Y5JA_1M.se.sorted.deduped.bam.bai
-  - for i in {"chr8","chr4","chrY"}; do 
+  - >
+    for i in {"chr8","chr4","chrY"}; do 
       echo "samtools view -b Q-Y5JA_1M.se.sorted.deduped.bam ${i} > Q-Y5JA_1M.se.sorted.deduped.${i}.bam; samtools index -@ 20 -b Q-Y5JA_1M.se.sorted.deduped.${i}.bam Q-Y5JA_1M.se.sorted.deduped.${i}.bam.bai;";
       done | singularity run 'docker://bicf/gudmaprbkdedup:2.0.0' parallel -j 20 -k
   - pytest -m dedupData
-  
-makeBigWig:
-  stage: unit
-  script:
-  - singularity run 'docker://bicf/deeptools3.3:2.0.0' bamCoverage -p 20 -b ./test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.bam -o Q-Y5JA_1M.se.bw
-  - pytest -m makeBigWig
 
-makeFeatureCounts:
+countData:
   stage: unit
   script:
   - singularity run 'docker://bicf/subread2:2.0.0' featureCounts -R SAM -p -G /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/genome.fna -T 20 -s 1 -a /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/genome.gtf -o Q-Y5JA_1M.se.featureCounts -g 'gene_name' --primary --ignoreDup -B ./test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.bam 
   - singularity run 'docker://bicf/subread2:2.0.0' Rscript ./workflow/scripts/calculateTPM.R --count Q-Y5JA_1M.se.featureCounts
   - pytest -m makeFeatureCounts
 
+makeBigWig:
+  stage: unit
+  script:
+  - singularity run 'docker://bicf/deeptools3.3:2.0.0' bamCoverage -p 20 -b ./test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.bam -o Q-Y5JA_1M.se.bw
+  - pytest -m makeBigWig
+
 fastqc:
   stage: unit
   script:
   - singularity run 'docker://bicf/fastqc:2.0.0' ./test_data/fastq/small/Q-Y5JA_1M.R1.fastq.gz -o .
   - pytest -m fastqc
 
-inferMetadata:
+dataQC:
   stage: unit
   script:
+  - echo -e  "geneID\tchrom\ttx_start\ttx_end\tTIN" > Q-Y5JA_1M.se.sorted.deduped.tin.xls
   - for i in {"chr8","chr4","chrY"}; do
-    echo " tin.py -i /project/BICF/BICF_Core/shared/gudmap/test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.${i}.bam -r /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/bed/genome.bed; cat Q-Y5JA_1M.se.sorted.deduped.${i}.tin.xls | tr -s \"\\w\" \"\\t\" | grep -P \"\\t${i}\\t\";"; done | singularity run 'docker://bicf/rseqc3.0:2.0.0' parallel -j 20 -k > Q-Y5JA_1M.se.sorted.deduped.tin.xls
-  - pytest -m inferMetadata
+    echo "tin.py -i /project/BICF/BICF_Core/shared/gudmap/test_data/bam/small/Q-Y5JA_1M.se.sorted.deduped.${i}.bam -r /project/BICF/BICF_Core/shared/gudmap/references/GRCh38.p12.v31/bed/genome.bed; cat Q-Y5JA_1M.se.sorted.deduped.${i}.tin.xls | tr -s \"\\w\" \"\\t\" | grep -P \"\\t${i}\\t\";"; done | singularity run 'docker://bicf/rseqc3.0:2.0.0' parallel -j 20 -k >> Q-Y5JA_1M.se.sorted.deduped.tin.xls
+  - pytest -m dataQC
+
 
 integration_se:
   stage: integration

README.md

@@ -9,7 +9,7 @@ RNA-Seq Analytic Pipeline for GUDMAP/RBK
 
 Introduction
 ------------
-This pipeline was created to be a standard mRNA-sequencing analysis pipeline which integrates with the GUDMAP and RBK consortium data-hub.
+This pipeline was created to be a standard mRNA-sequencing analysis pipeline which integrates with the GUDMAP and RBK consortium data-hub. It is designed to run on the HPC cluster ([BioHPC](https://portal.biohpc.swmed.edu)) at UT Southwestern Medical Center (in conjunction with the standard nextflow profile: config `biohpc.config`)
 
 ![flowchart](docs/RNA-Seq%20Pipeline%20Design%20Flowchart.jpg "Flowchart")
 
@@ -21,18 +21,19 @@ This pipeline is also capable of being run on AWS. To do so:
   * Replace workDir with the S3 bucket generated
   * Change region if different
   * Change queue to the aws batch queue generated
-* The user must have awscli configured with an appropriate authentication (with ```aws configure``` and access keys) in the environment which nextflow will be run
-* Add ```-profile ``` with the name aws config which was customized
+* The user must have awscli configured with an appropriate authentication (with `aws configure` and access keys) in the environment which nextflow will be run
+* Add `-profile` with the name aws config which was customized
 
 To Run:
 -------
 * Available parameters:
-  * ```--deriva``` active **credential.json** file from [deriva-auth](https://github.com/informatics-isi-edu/gudmap-rbk/wiki/Uploading-files-via-Deriva-client-tools#from-a-remote-server)
-  * ```--bdbag``` active **cookies.txt** file from [deriva-auth](https://github.com/informatics-isi-edu/gudmap-rbk/wiki/Uploading-files-via-Deriva-client-tools#from-a-remote-server)
-  * ```--repRID``` mRNA-seq replicate RID
-  * ```--refMoVersion``` mouse reference version ***(optional)***
-  * ```--refHuVersion``` human reference version ***(optional)***
-  * ```-profile``` config profile to use: standard = local processes on BioHPC (default), biohpc = BioHPC cluster, aws_ondemand = AWS Batch on-demand instant requests, aws_spot = AWS Batch spot instance requests ***(optional)***
+  * `--deriva` active **credential.json** file from [deriva-auth](https://github.com/informatics-isi-edu/gudmap-rbk/wiki/Uploading-files-via-Deriva-client-tools#from-a-remote-server)
+  * `--bdbag` active **cookies.txt** file from [deriva-auth](https://github.com/informatics-isi-edu/gudmap-rbk/wiki/Uploading-files-via-Deriva-client-tools#from-a-remote-server)
+  * `--repRID` mRNA-seq replicate RID
+  * `--refMoVersion` mouse reference version ***(optional)***
+  * `--refHuVersion` human reference version ***(optional)***
+  * `--refERCCVersion` human reference version ***(optional)***
+  * `-profile` config profile to use: standard = processes on BioHPC cluster, aws_ondemand = AWS Batch on-demand instant requests, aws_spot = AWS Batch spot instance requests ***(optional)***
 * NOTES:
   * once deriva-auth is run and authenticated, the two files above are saved in ```~/.deriva/``` (see official documents from [deriva](https://github.com/informatics-isi-edu/deriva-client#installer-packages-for-windows-and-macosx) on the lifetime of the credentials)
   * reference version consists of Genome Reference Consortium version, patch release and GENCODE annotation release # (leaving the params blank will use the default version tied to the pipeline version)
@@ -96,4 +97,4 @@ Please cite in publications: Pipeline was developed by BICF from funding provide
 Pipeline Directed Acyclic Graph
 -------------------------------
 
-![dag](docs/dag.png "DAG")
+![dag](docs/dag.png "DAG")
\ No newline at end of file

workflow/conf/aws_ondemand.config

-workDir = 's3://gudmap.rbk/work'
+workDir = 's3://'
 aws.client.storageEncryption = 'AES256'
 aws {
   region = 'us-east-2'
@@ -9,29 +9,7 @@ aws {
 
 process {
   executor = 'awsbatch'
-  queue = 'highpriority-3278a8b0-1fc8-11ea-b1ac-021e2396e2cc'
+  queue = 'highpriority-'
   cpus = 1
   memory = '2 GB'
-
-  withName:parseMetadata {
-    cpus = 5
-  }
-  withName:getRef {
-    cpus = 8
-  }
-  withName:trimData {
-    cpus = 8
-    memory = '2 GB'
-  }
-  withName:alignData {
-    cpus = 50
-    memory = '5 GB'
-  }
-  withName:dedupData {
-    cpus = 2
-    memory = '20 GB'
-  }
-  withName:fastqc {
-    memory = '5 GB'
-  }
-}
+}
\ No newline at end of file

workflow/conf/aws_spot.config

-workDir = 's3://gudmap.rbk/work'
+workDir = 's3://'
 aws.client.storageEncryption = 'AES256'
 aws {
   region = 'us-east-2'
@@ -9,29 +9,7 @@ aws {
 
 process {
   executor = 'awsbatch'
-  queue = 'default-3278a8b0-1fc8-11ea-b1ac-021e2396e2cc'
+  queue = 'default-'
   cpus = 1
   memory = '2 GB'
-
-  withName:parseMetadata {
-    cpus = 5
-  }
-  withName:getRef {
-    cpus = 8
-  }
-  withName:trimData {
-    cpus = 8
-    memory = '2 GB'
-  }
-  withName:alignData {
-    cpus = 50
-    memory = '5 GB'
-  }
-  withName:dedupData {
-    cpus = 2
-    memory = '20 GB'
-  }
-  withName:fastq  {
-    memory = '5 GB'
-  }
 }

workflow/conf/biohpc.config

@@ -3,34 +3,49 @@ process {
   queue = 'super'
   clusterOptions = '--hold'
 
-  withName:getBag {
+  withName: getBag {
     executor = 'local'
   }
-  withName:getData {
+  withName: getData {
     executor = 'local'
   }
-  withName:parseMetadata {
+  withName: parseMetadata {
     executor = 'local'
   }
-  withName:getRef {
+  withName: trimData {
+    queue = 'super'
+  }
+  withName: getRefInfer {
+    executor = 'local'
+  }
+  withName: downsampleData {
     executor = 'local'
   }
-  withName:trimData {
+  withName: alignSampleData {
+    queue = 'super'
+  }
+  withName: inferMetadata {
     queue = 'super'
   }
-  withName:alignData {
+  withName: getRef {
+    executor = 'local'
+  }
+  withName: alignData {
     queue = '256GB,256GBv1'
   }
   withName: dedupData {
     queue = 'super'
   }
-  withName:fastqc {
+  withName: countData {
+    queue = 'super'
+  }
+  withName: makeBigWig {
     queue = 'super'
   }
-  withName:inferMetadata {
+  withName: fastqc {
     queue = 'super'
   }
-  withName: makeBigWig {
+  withName: dataQC {
     queue = 'super'
   }
 }

workflow/nextflow.config

@@ -20,30 +20,42 @@ process {
   withName: parseMetadata {
     container = 'bicf/python3:1.3'
   }
-  withName: getRef {
-    container = 'bicf/awscli:1.1'
-  }
   withName: trimData {
     container = 'bicf/trimgalore:1.1'
   }
+  withName: getRefInfer {
+    container = 'bicf/awscli:1.1'
+  }
+  withName: downsampleData {
+    container = 'bicf/seqtk:2.0.0'
+  }
+  withName: alignSampleData {
+    container = 'bicf/gudmaprbkaligner:2.0.0'
+  }
+  withName: inferMetadata {
+    container = 'bicf/rseqc3.0:2.0.0'
+  }
+  withName: getRef {
+    container = 'bicf/awscli:1.1'
+  }
   withName: alignData {
     container = 'bicf/gudmaprbkaligner:2.0.0'
   }
   withName: dedupData {
     container = 'bicf/gudmaprbkdedup:2.0.0'
   }
+  withName: countData {
+    container = 'bicf/subread2:2.0.0'
+  }
   withName: makeBigWig {
     container = 'bicf/deeptools3.3:2.0.0'
   }
   withName: fastqc {
     container = 'bicf/fastqc:2.0.0'
   }
-  withName:inferMetadata{
+  withName: dataQC {
     container = 'bicf/rseqc3.0:2.0.0'
   }
-  withName: makeFeatureCounts {
-    container = 'bicf/subread2:2.0.0'
-  }
 }
 
 trace {

workflow/scripts/parseMeta.py

@@ -10,7 +10,6 @@ def get_args():
     parser.add_argument('-r', '--repRID',help="The replicate RID.",required=True)
     parser.add_argument('-m', '--metaFile',help="The metadata file to extract.",required=True)
     parser.add_argument('-p', '--parameter',help="The parameter to extract.",required=True)
-    parser.add_argument('-e', '--endsManual',help="The endness.",required=False)
     args = parser.parse_args()
     return args
 
@@ -56,12 +55,9 @@ def main():
     # Get strandedness metadata from 'Experiment Settings.csv'
     if (args.parameter == "stranded"):
         if (metaFile.Has_Strand_Specific_Information.unique() == "yes"):
-            if (args.endsManual=="se"):
-                stranded = "--rna-strandness F"
-            elif (args.endsManual=="pe"):
-                stranded = "--rna-strandness FR"
+            stranded = "stranded"
         elif (metaFile.Has_Strand_Specific_Information.unique() == "no"):
-            stranded = ""
+            stranded = "unstranded"
         else:
             print("Stranded metadata not match expected options: " + metaFile.Has_Strand_Specific_Information.unique())
             exit(1)

workflow/tests/test_dataQC.py

+#!/usr/bin/env python3
+
+import pytest
+import pandas as pd
+from io import StringIO
+import os
+
+test_output_path = os.path.dirname(os.path.abspath(__file__)) + \
+                '/../../'
+
+@pytest.mark.dataQC
+def test_dataQC():
+    assert os.path.exists(os.path.join(test_output_path, 'Q-Y5JA_1M.se.sorted.deduped.tin.xls'))
+

workflow/tests/test_downsampleData.py

+#!/usr/bin/env python3
+
+import pytest
+import pandas as pd
+from io import StringIO
+import os
+
+test_output_path = os.path.dirname(os.path.abspath(__file__)) + \
+                '/../../'
+
+@pytest.mark.downsampleData
+def test_downsampleData():
+    assert os.path.exists(os.path.join(test_output_path, 'sampled.1.fq'))
\ No newline at end of file

workflow/tests/test_inferMetadata.py

@@ -10,5 +10,5 @@ test_output_path = os.path.dirname(os.path.abspath(__file__)) + \
 
 @pytest.mark.inferMetadata
 def test_inferMetadata():
-    assert os.path.exists(os.path.join(test_output_path, 'Q-Y5JA_1M.se.sorted.deduped.tin.xls'))
+    assert os.path.exists(os.path.join(test_output_path, 'Q-Y5JA_1M.se.inferMetadata.log'))