Update by dgwoo@unist.ac.kr

Signed-off-by: Jaewoong Lee <230@fumire.moe>

Author: Fumire

02_Data_pre-processing_for_variant_discovery/02_1_BWA.bash

@@ -3,4 +3,3 @@ set -euo pipefail
 IFS=$'\n\t'
 
 python3 02_1_BWA.py /BiO/Store/Standard-Pipeline/cn95N_S0_L009_R1_001.fastq.gz /BiO/Store/Standard-Pipeline/cn95N_S0_L009_R2_001.fastq.gz .
-python3 02_1_BWA.py /BiO/Store/Standard-Pipeline/cn95P_S0_L009_R1_001.fastq.gz /BiO/Store/Standard-Pipeline/cn95P_S0_L009_R2_001.fastq.gz .

02_Data_pre-processing_for_variant_discovery/02_1_BWA.py

@@ -1,66 +1,70 @@
-#!/usr/bin/env python3
-"""
-02_1_BWA.py: Mapping with BWA
-"""
-import argparse
-import configparser
-import os
-import subprocess
-
-parser = argparse.ArgumentParser()
-
-parser.add_argument("input", help="Input FASTQ file", nargs=2)
-parser.add_argument("output", help="Output directory", default=os.getcwd())
-parser.add_argument("-c", "--config", help="config INI file", default="../config.ini")
-
-parser.add_argument("-n", "--dryrun", help="Don't actually run any recipe; just make .SH only", default=False, action="store_true")
-
-args = parser.parse_args()
-
-config = configparser.ConfigParser(interpolation=configparser.ExtendedInterpolation())
-config.read(args.config)
-
-args.input.sort()
-name = args.input[0].split("/")[-1].split("_")[0]
-args.output = os.path.realpath(args.output)
-
-# Mapping
-with open(f"BWA_{name}.sh", "w") as sh:
-    sh.write("#!/bin/bash\n")
-    sh.write(f"{config['TOOLS']['bwa']} mem -M -t {config['DEFAULT']['threads']} -R '@RG\\tID:{name}\\tPL:ILLUMINA\\tLB:{name}\\tSM:{name}\\tCN:UNIST' -v 3 {config['REFERENCES']['fasta']} {args.input[0]} {args.input[0]} | {config['TOOLS']['samtools']} view --bam --with-header --threads {config['DEFAULT']['threads']} --reference {config['REFERENCES']['fasta']} --output {args.output}/{name}.bam")
-
-if not args.dryrun:
-    mapping_job_id = subprocess.check_output(f"sbatch --chdir=$(realpath .) --cpus-per-task={config['DEFAULT']['threads']} --error='%x-%A.txt' --job-name='BWA_{name}' --mem={config['DEFAULT']['memory']}G --output='%x-%A.txt' --export=ALL BWA_{name}.sh", encoding="utf-8", shell=True).split()[-1]
-
-# Sort
-with open(f"Sort_{name}.sh", "w") as sh:
-    sh.write("#!/bin/bash\n")
-    sh.write(f"{config['TOOLS']['samtools']} sort -l 9 --threads {config['DEFAULT']['threads']} -m {int(config['DEFAULT']['memory']) // int(config['DEFAULT']['threads'])}G --reference {config['REFERENCES']['fasta']} --write-index -o {args.output}/{name}.Sort.bam {args.output}/{name}.bam")
-
-if not args.dryrun:
-    sorting_job_id = subprocess.check_output(f"sbatch --dependency=afterok:{mapping_job_id} --chdir=$(realpath .) --cpus-per-task={config['DEFAULT']['threads']} --error='%x-%A.txt' --job-name='Sort_{name}' --mem={config['DEFAULT']['memory']}G --output='%x-%A.txt' --export=ALL Sort_{name}.sh", encoding="utf-8", shell=True).split()[-1]
-
-# Mark Duplicates
-with open(f"MarkDup_{name}.sh", "w") as sh:
-    sh.write("#!/bin/bash\n")
-    sh.write(f"{config['TOOLS']['gatk']} MarkDuplicatesSpark --input {args.output}/{name}.Sort.bam --output {args.output}/{name}.Sort.MarkDuplicates.bam --reference {config['REFERENCES']['fasta']} --metrics-file {args.output}/{name}.Sort.MarkDuplicates.metrics --duplicate-tagging-policy 'OpticalOnly' -- --spark-master 'local[{config['DEFAULT']['threads']}]' --spark-verbosity 'INFO'")
-
-if not args.dryrun:
-    markduplicates_job_id = subprocess.check_output(f"sbatch --dependency=afterok:{sorting_job_id} --chdir=$(realpath .) --cpus-per-task={config['DEFAULT']['threads']} --error='%x-%A.txt' --job-name='MarkDup_{name}' --mem={config['DEFAULT']['memory']}G --output='%x-%A.txt' --export=ALL MarkDup_{name}.sh", encoding="utf-8", shell=True).split()[-1]
-
-# Base Quality Score Recalibration (BQSR)
-with open(f"BQSR_{name}.sh", "w") as sh:
-    sh.write("#!/bin/bash\n")
-    sh.write(f"{config['TOOLS']['gatk']} BaseRecalibrator --input {args.output}/{name}.Sort.MarkDuplicates.bam --reference {config['REFERENCES']['fasta']} --output {args.output}/{name}.Sort.MarkDuplicates.BQSR.table --create-output-bam-index true")
-    for site in config['REFERENCES']['sites'].split(" "):
-        sh.write(f" --known-sites {site}")
-
-if not args.dryrun:
-    BQSR_job_id = subprocess.check_output(f"sbatch --dependency=afterok:{markduplicates_job_id} --chdir=$(realpath .) --cpus-per-task={config['DEFAULT']['threads']} --error='%x-%A.txt' --job-name='BQSR_{name}' --mem={config['DEFAULT']['memory']}G --output='%x-%A.txt' --export=ALL BQSR_{name}.sh", encoding="utf-8", shell=True).split()[-1]
-
-with open(f"ApplyBQSR_{name}.sh", "w") as sh:
-    sh.write("#!/bin/bash\n")
-    sh.write(f"{config['TOOLS']['gatk']} ApplyBQSR --bqsr-recal-file {args.output}/{name}.Sort.MarkDuplicates.BQSR.table --input {args.output}/{name}.Sort.MarkDuplicates.bam --output {args.output}/{name}.Sort.MarkDuplicates.BQSR.bam --reference {config['REFERENCES']['fasta']} --create-output-bam-index true")
-
-if not args.dryrun:
-    ApplyBQSR_job_id = subprocess.check_output(f"sbatch --dependency=afterok:{BQSR_job_id} --chdir=$(realpath .) --cpus-per-task={config['DEFAULT']['threads']} --error='%x-%A.txt' --job-name='ApplyBQSR_{name}' --mem={config['DEFAULT']['memory']}G --output='%x-%A.txt' --export=ALL ApplyBQSR_{name}.sh", encoding="utf-8", shell=True).split()[-1]
+#!/usr/bin/env python3
+"""
+02_1_BWA.py: Mapping with BWA
+"""
+import argparse
+import os
+import sys
+from pipeline_utils import PipelineManagerBase
+
+sys.path.append(os.path.dirname(os.path.abspath(os.path.dirname(__file__))))
+
+
+class PipelineManager(PipelineManagerBase):
+    def __init__(self, input_files, output, config_file, dryrun):
+        super().__init__(config_file, dryrun, output_dir=output)
+        self.input_files = sorted(input_files)
+        self.name = os.path.basename(self.input_files[0]).split("_")[0]
+
+    def run_bwa(self, dependency_id=None):
+        command = f"{self.config['TOOLS']['bwa']} mem -M -t {self.config['DEFAULT']['threads']} -R '@RG\\tID:{self.name}\\tPL:ILLUMINA\\tLB:{self.name}\\tSM:{self.name}\\tCN:UNIST' -v 3 {self.config['REFERENCES']['fasta']} {self.input_files[0]} {self.input_files[0]} | {self.config['TOOLS']['samtools']} view --bam --with-header --threads {self.config['DEFAULT']['threads']} --reference {self.config['REFERENCES']['fasta']} --output {self.output_dir}/{self.name}.bam"
+        self.create_sh("BWA", command)
+        return self.submit_job("BWA")
+
+    def run_sort(self, dependency_id=None):
+        command = f"{self.config['TOOLS']['samtools']} sort -l 9 --threads {self.config['DEFAULT']['threads']} -m {int(self.config['DEFAULT']['memory']) // int(self.config['DEFAULT']['threads'])}G --reference {self.config['REFERENCES']['fasta']} --write-index -o {self.output_dir}/{self.name}.Sort.bam {self.output_dir}/{self.name}.bam"
+        self.create_sh("Sort", command)
+        return self.submit_job("Sort", dependency_id=dependency_id)
+
+    def run_mark_duplicates(self, dependency_id=None):
+        command = f"{self.config['TOOLS']['gatk']} MarkDuplicatesSpark --input {self.output_dir}/{self.name}.Sort.bam --output {self.output_dir}/{self.name}.Sort.MarkDuplicates.bam --reference {self.config['REFERENCES']['fasta']} --metrics-file {self.output_dir}/{self.name}.Sort.MarkDuplicates.metrics --duplicate-tagging-policy 'OpticalOnly' -- --spark-master 'local[{self.config['DEFAULT']['threads']}]' --spark-verbosity 'INFO'"
+        self.create_sh("MarkDup", command)
+        return self.submit_job("MarkDup", dependency_id=dependency_id)
+
+    def run_bqsr(self, dependency_id=None):
+        known_sites = " ".join([f"--known-sites {site}" for site in self.config["REFERENCES"]["sites"].split(" ")])
+        command = f"{self.config['TOOLS']['gatk']} BaseRecalibrator --input {self.output_dir}/{self.name}.Sort.MarkDuplicates.bam --reference {self.config['REFERENCES']['fasta']} --output {self.output_dir}/{self.name}.Sort.MarkDuplicates.BQSR.table --create-output-bam-index true {known_sites}"
+        self.create_sh("BQSR", command)
+        return self.submit_job("BQSR", dependency_id=dependency_id)
+
+    def run_apply_bqsr(self, dependency_id=None):
+        command = f"{self.config['TOOLS']['gatk']} ApplyBQSR --bqsr-recal-file {self.output_dir}/{self.name}.Sort.MarkDuplicates.BQSR.table --input {self.output_dir}/{self.name}.Sort.MarkDuplicates.bam --output {self.output_dir}/{self.name}.Sort.MarkDuplicates.BQSR.bam --reference {self.config['REFERENCES']['fasta']} --create-output-bam-index true"
+        self.create_sh("ApplyBQSR", command)
+        return self.submit_job("ApplyBQSR", dependency_id=dependency_id)
+
+
+def parse_arguments():
+    parser = argparse.ArgumentParser()
+
+    parser.add_argument("input", help="Input FASTQ file", nargs=2)
+    parser.add_argument("output", help="Output directory", default=os.getcwd())
+    parser.add_argument("-c", "--config", help="config INI file", default="../config.ini")
+    parser.add_argument("-n", "--dryrun", help="Don't actually run any recipe; just make .SH only", default=False, action="store_true")
+
+    return parser.parse_args()
+
+
+def main():
+    args = parse_arguments()
+
+    pipeline = PipelineManager(input_files=args.input, output=args.output, config_file=args.config, dryrun=args.dryrun)
+
+    mapping_job_id = pipeline.run_bwa()
+    sort_job_id = pipeline.run_sort(dependency_id=mapping_job_id)
+    mark_duplicates_job_id = pipeline.run_mark_duplicates(dependency_id=sort_job_id)
+    bqsr_job_id = pipeline.run_bqsr(dependency_id=mark_duplicates_job_id)
+    apply_bqsr_job_id = pipeline.run_apply_bqsr(dependency_id=bqsr_job_id)
+
+
+if __name__ == "__main__":
+    main()

03_Somatic_short_variant_discovery/03_1_Mutect2.bash

@@ -2,4 +2,4 @@
 set -euo pipefail
 #IFS=$'\n\t'
 
-python3 03_1_Mutect2.py /BiO/Research/Standard-Pipeline/02_Data_pre-processing_for_variant_discovery/cn95N.Sort.MarkDuplicates.BQSR.bam /BiO/Research/Standard-Pipeline/02_Data_pre-processing_for_variant_discovery/cn95P.Sort.MarkDuplicates.BQSR.bam cn95.maf
+python3 03_1_Mutect2.py /BiO/Research/Standard-Pipeline/02_Data_pre-processing_for_variant_discovery/cn95N.Sort.MarkDuplicates.BQSR.bam /BiO/Research/Standard-Pipeline/02_Data_pre-processing_for_variant_discovery/cn95P.Sort.MarkDuplicates.BQSR.bam cn95.maf -p