RNA-Seq germline variant calling pipeline

.gitattributes

@@ -21,3 +21,5 @@ bin/check-job-alive !text !filter !merge !diff
 *.py !text !filter !merge !diff
 *.pl !text !filter !merge !diff
 *.conf !text !filter !merge !diff
+.fa filter=lfs diff=lfs merge=lfs -text
+*.fa filter=lfs diff=lfs merge=lfs -text

tests/tools/input/GRCh38.chr1_chr19.dict

@@ -0,0 +1,3 @@
+@HD	VN:1.0	SO:unsorted
+@SQ	SN:chr1	LN:248956422	M5:6aef897c3d6ff0c78aff06ac189178dd
+@SQ	SN:chr19	LN:58617616	M5:85f9f4fc152c58cb7913c06d6b98573a

tests/tools/input/GRCh38.chr1_chr19.fa.fai

@@ -0,0 +1,2 @@
+chr1	248956422	112	70	71
+chr19	58617616	252513180	70	71

tests/tools/input/chr1.interval_list

@@ -0,0 +1,19 @@
+@HD	VN:1.5	SO:coordinate
+@SQ	SN:chr1	LN:248956422	M5:6aef897c3d6ff0c78aff06ac189178dd	AS:38	UR:/seq/references/Homo_sapiens_assembly38/v0/Homo_sapiens_assembly38.fasta	SP:Homo sapiens
+@SQ	SN:chr19	LN:58617616	M5:85f9f4fc152c58cb7913c06d6b98573a	AS:38	UR:/seq/references/Homo_sapiens_assembly38/v0/Homo_sapiens_assembly38.fasta	SP:Homo sapiens
+chr1	10001	207666	+	. intersection ACGTmer
+chr1	257667	297968	+	. intersection ACGTmer
+chr1	347969	535988	+	. intersection ACGTmer
+chr1	585989	2702781	+	. intersection ACGTmer
+chr1	2746291	12954384	+	. intersection ACGTmer
+chr1	13004385	16799163	+	. intersection ACGTmer
+chr1	16849164	29552233	+	. intersection ACGTmer
+chr1	29553836	121976459	+	. intersection ACGTmer
+chr1	122026460	124977944	+	. intersection ACGTmer
+chr1	124978327	125130246	+	. intersection ACGTmer
+chr1	125131848	125171347	+	. intersection ACGTmer
+chr1	125173584	125184587	+	. intersection ACGTmer
+chr1	143184588	223558935	+	. intersection ACGTmer
+chr1	223608936	228558364	+	. intersection ACGTmer
+chr1	228608365	248946422	+	. intersection ACGTmer
+chr19	20002	208000	+	. intersection ACGTmer

tests/tools/input_json/gatk4_apply_bqsr.json

@@ -0,0 +1,5 @@
+{
+    "bam": "tests/tools/input/test_rnaseq_variant.bam",
+    "bam_index": "tests/tools/input/test_rnaseq_variant.bam.bai",
+    "recalibration_report": "tests/tools/input/test_rnaseq_variant.recal.txt"
+}

tests/tools/input_json/gatk4_base_recalibrator.json

@@ -0,0 +1,11 @@
+{
+    "bam": "tests/tools/input/test_rnaseq_variant.bam",
+    "bam_index": "tests/tools/input/test_rnaseq_variant.bam.bai",
+    "fasta": "tests/tools/input/GRCh38.chr1_chr19.fa",
+    "fasta_index": "tests/tools/input/GRCh38.chr1_chr19.fa.fai",
+    "dict": "tests/tools/input/GRCh38.chr1_chr19.dict",
+    "dbSNP_vcf":"tests/tools/input/Homo_sapiens_assembly38.dbsnp138.top5000.vcf",
+    "dbSNP_vcf_index": "tests/tools/input/Homo_sapiens_assembly38.dbsnp138.top5000.vcf.idx",
+    "known_indels_sites_VCFs": ["tests/tools/input/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz"],
+    "known_indels_sites_indices": ["tests/tools/input/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz.tbi"]
+}

tests/tools/input_json/gatk4_haplotype_caller.json

@@ -0,0 +1,10 @@
+{
+    "bam": "tests/tools/input/test_rnaseq_variant.bam",
+    "bam_index": "tests/tools/input/test_rnaseq_variant.bam.bai",
+    "fasta": "tests/tools/input/GRCh38.chr1_chr19.fa", 
+    "fasta_index": "tests/tools/input/GRCh38.chr1_chr19.fa.fai",
+    "dict": "tests/tools/input/GRCh38.chr1_chr19.dict",
+    "dbSNP_vcf": "tests/tools/input/Homo_sapiens_assembly38.dbsnp138.top5000.vcf",
+    "dbSNP_vcf_index": "tests/tools/input/Homo_sapiens_assembly38.dbsnp138.top5000.vcf.idx",
+    "interval_list": "tests/tools/input/chr1.interval_list"
+}

tests/tools/input_json/gatk4_split_n_cigar_reads.json

@@ -0,0 +1,9 @@
+{
+    "bam": "tests/tools/input/test.bam",
+    "bam_index": "tests/tools/input/test.bam.bai",
+    "fasta": "tests/tools/input/GRCh38.chr1_chr19.fa",
+    "fasta_index": "tests/tools/input/GRCh38.chr1_chr19.fa.fai",
+    "dict": "tests/tools/input/GRCh38.chr1_chr19.dict",
+    "prefix": "split",
+    "interval_list": "tests/tools/input/wgs_calling_regions.hg38.interval_list"
+}

tests/tools/input_json/gatk4_variant_filtration.json

@@ -0,0 +1,7 @@
+{
+    "vcf": "tests/tools/input/test1.vcf.gz",
+    "vcf_index": "tests/tools/input/test1.vcf.gz.tbi",
+    "fasta": "tests/tools/input/GRCh38.chr1_chr19.fa",
+    "fasta_index": "tests/tools/input/GRCh38.chr1_chr19.fa.fai",
+    "dict": "tests/tools/input/GRCh38.chr1_chr19.dict"
+}

tests/tools/input_json/picard_merge_vcfs.json

@@ -0,0 +1,5 @@
+{
+    "vcfs": ["tests/tools/input/test1.vcf.gz", "tests/tools/input/test2.vcf.gz"],
+    "vcfs_indexes": ["tests/tools/input/test1.vcf.gz.tbi", "tests/tools/input/test2.vcf.gz.tbi"],
+    "output_vcf_name": "test.vcf.gz"
+}

tests/tools/input_json/sambamba_merge.json

@@ -0,0 +1,4 @@
+{
+    "bams": ["tests/tools/input/possorted_genome_bam.bam", "tests/tools/input/test.bwa_aln_pe.bam"],
+    "prefix": "test"
+}

tests/tools/test_gatk4.yaml

@@ -0,0 +1,49 @@
+- name: gatk4_split_n_cigar_reads
+  tags:
+    - miniwdl
+    - gatk4
+  command: >-
+    miniwdl run --verbose -d test-output/. --task split_n_cigar_reads -i tests/tools/input_json/gatk4_split_n_cigar_reads.json tools/gatk4.wdl
+  files:
+    - path: test-output/out/split_n_reads_bam/split.bam
+    - path: test-output/out/split_n_reads_bam_index/split.bam.bai
+    - path: test-output/out/split_n_reads_bam_md5/split.bam.md5
+
+- name: gatk4_base_recalibrator
+  tags:
+    - miniwdl
+    - gatk4
+  command: >-
+    miniwdl run --verbose -d test-output/. --task base_recalibrator -i tests/tools/input_json/gatk4_base_recalibrator.json tools/gatk4.wdl
+  files:
+    - path: test-output/out/recalibration_report/test_rnaseq_variant.recal.txt
+
+- name: gatk4_apply_bqsr
+  tags:
+    - miniwdl
+    - gatk4
+  command: >-
+    miniwdl run --verbose -d test-output/. --task apply_bqsr -i tests/tools/input_json/gatk4_apply_bqsr.json tools/gatk4.wdl
+  files:
+    - path: test-output/out/recalibrated_bam/test_rnaseq_variant.bqsr.bam
+    - path: test-output/out/recalibrated_bam_index/test_rnaseq_variant.bqsr.bam.bai
+
+
+- name: gatk4_haplotype_caller
+  tags:
+    - miniwdl
+    - gatk4
+  command: >-
+    miniwdl run --verbose -d test-output/. --task haplotype_caller -i tests/tools/input_json/gatk4_haplotype_caller.json tools/gatk4.wdl
+  files:
+    - path: test-output/out/vcf/test_rnaseq_variant.vcf.gz
+
+- name: gatk4_variant_filtration
+  tags:
+    - miniwdl
+    - gatk4
+  command: >-
+    miniwdl run --verbose -d test-output/. --task variant_filtration -i tests/tools/input_json/gatk4_variant_filtration.json tools/gatk4.wdl
+  files:
+    - path: test-output/out/vcf_filtered/test1.filtered.vcf.gz
+    - path: test-output/out/vcf_filtered_index/test1.filtered.vcf.gz.tbi

tests/tools/test_picard.yaml

@@ -105,3 +105,32 @@
   files:
     - path: test-output/out/quality_score_distribution_txt/test.bwa_aln_pe.QualityScoreDistribution.txt
     - path: test-output/out/quality_score_distribution_pdf/test.bwa_aln_pe.QualityScoreDistribution.pdf
+
+- name: picard_merge_vcfs
+  tags:
+    - miniwdl
+    - picard
+  command: >-
+    miniwdl run -d test-output/. --task merge_vcfs -i tests/tools/input_json/picard_merge_vcfs.json tools/picard.wdl
+  files:
+    - path: test-output/out/output_vcf/test.vcf.gz
+
+- name: picard_scatter_interval_list
+  tags:
+    - miniwdl
+    - picard
+  command: >-
+    miniwdl run -d test-output/. --task scatter_interval_list tools/picard.wdl interval_list="tests/tools/input/wgs_calling_regions.hg38.interval_list" scatter_count=3
+  files:
+    - path: test-output/out/interval_lists_scatter/0/1scattered.interval_list
+    - path: test-output/out/interval_lists_scatter/1/2scattered.interval_list
+    - path: test-output/out/interval_lists_scatter/2/3scattered.interval_list
+
+- name: picard_create_sequence_dictionary
+  tags:
+    - miniwdl
+    - picard
+  command: >-
+    miniwdl run -d test-output/. --task create_sequence_dictionary tools/picard.wdl fasta="tests/tools/input/GRCh38.chrY_chrM.fa.gz" outfile_name="GRCh38.chrY_chrM.dict"
+  files:
+    - path: test-output/out/dictionary/GRCh38.chrY_chrM.dict

tests/tools/test_sambamba.yaml

@@ -0,0 +1,46 @@
+- name: sambamba_index
+  tags:
+    - miniwdl
+    - sambamba
+  command: >-
+    miniwdl run -d test-output/. --task index tools/sambamba.wdl bam="tests/tools/input/test.bwa_aln_pe.bam"
+  files:
+    - path: test-output/out/bam_index/test.bwa_aln_pe.bam.bai
+
+- name: sambamba_merge
+  tags:
+    - miniwdl
+    - sambamba
+  command: >-
+    miniwdl run -d test-output/. --task merge -i tests/tools/input_json/sambamba_merge.json tools/sambamba.wdl
+  files:
+    - path: test-output/out/merged_bam/test.bam
+
+- name: sambamba_sort
+  tags:
+    - miniwdl
+    - sambamba
+  command: >-
+    miniwdl run -d test-output/. --task sort tools/sambamba.wdl bam="tests/tools/input/test.bwa_aln_pe.bam"
+  files:
+    - path: test-output/out/sorted_bam/test.bwa_aln_pe.sorted.bam
+
+- name: sambamba_flagstat
+  tags:
+    - miniwdl
+    - sambamba
+  command: >-
+    miniwdl run -d test-output/. --task flagstat tools/sambamba.wdl bam="tests/tools/input/test.bwa_aln_pe.bam"
+  files:
+    - path: test-output/out/flagstat_report/test.bwa_aln_pe.flagstat.txt
+
+- name: sambamba_markdup
+  tags:
+    - miniwdl
+    - sambamba
+  command: >-
+    miniwdl run -d test-output/. --task markdup tools/sambamba.wdl bam="tests/tools/input/test.bwa_aln_pe.bam"
+  files:
+    - path: test-output/out/duplicate_marked_bam/test.bwa_aln_pe.markdup.bam
+    - path: test-output/out/duplicate_marked_bam_index/test.bwa_aln_pe.markdup.bam.bai
+    - path: test-output/out/markdup_log/test.bwa_aln_pe.markdup_log.txt

tests/tools/test_samtools.yaml

@@ -82,3 +82,12 @@
     - path: test-output/out/collated_bam/test.bwa_aln_pe.collated.bam
     - path: test-output/out/read_one_fastq_gz/test.bwa_aln_pe.R1.fastq.gz
     - path: test-output/out/read_two_fastq_gz/test.bwa_aln_pe.R2.fastq.gz
+
+- name: samtools_faidx
+  tags:
+    - miniwdl
+    - samtools
+  command: >-
+    miniwdl run -d test-output/. --task faidx tools/samtools.wdl fasta="tests/tools/input/test.fa"
+  files:
+    - path: test-output/out/fasta_index/test.fa.fai