Answers.Body: <p>There are a couple ways to approach this, one is to create output files based on the <code>BAM</code> input file name - the other approach is to make a new folder (based on the input <code>BAM</code> filename containing result files with the same name for each input. I generally use the first method, but either work.</p>
  <p>Here is a script that may run the 5 <code>GATK</code> functions. It calculates a prefix to be used for derivative output files <code>BASE</code>, and should run your pipeline on all the <code>BAM</code> files within the <code>INPUT_FOLDER</code>.</p>
  <pre><code>#!/bin/bash
  
  set -e
  
  INPUT_FOLDER=&quot;/data/Continuum/WES/results/&quot;
  
  for BAM_FILE in &quot;<span class="math-container">$INPUT_FOLDER/"*.bam
  do
          BASE=$</span>(basename &quot;$BAM_FILE&quot; .bam) # NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup
  
          # function 1
          ./gatk CollectSequencingArtifactMetrics \
                  -I &quot;<span class="math-container">$BAM_FILE" \
                  -O "$</span>BASE&quot; \
                  --FILE_EXTENSION .txt \
                  -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
          # function 2
          ./gatk GetPileupSummaries \
                  -I &quot;<span class="math-container">$BAM_FILE" \
                  -O "$</span>BASE.targeted_sequencing.table&quot; \
                  -V af-only-gnomad.hg38.vcf.gz \
                  -L resources_broad_hg38_v0_wgs_calling_regions.hg38.interval_list \
                  -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
          # function 3
          ./gatk CalculateContamination \
                  -I &quot;<span class="math-container">$BASE.targeted_sequencing.table" \
                  -O "$</span>BASE.targeted_sequencing.contamination.table&quot; \
  
          # function 4
          ./gatk GetSampleName \
                  -I &quot;<span class="math-container">$BAM_FILE" \
                  -O "$</span>BASE.targeted_sequencing.sample_name&quot;
  
          # function 5
          ./gatk Mutect2 \
                  -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta \
                  -I &quot;<span class="math-container">$BAM_FILE" \
                  -O "$</span>BASE&quot;.vcf \
                  -tumor &quot;$BASE.targeted_sequencing.sample_name&quot; \
                  --af-of-alleles-not-in-resource 2.5e-06 \
                  --germline-resource af-only-gnomad.hg38.vcf.gz \
                  -pon 1000g_pon.hg38.vcf.gz
  done
  </code></pre>
  <p>It may need some adjustments - I did not test this, and used the <code>GATK</code> incantations you provided (thank you for including all of your steps!)</p>
  <p>It includes the addition of the command <code>set -e</code>, which depending on GATKs return status may stop the script if any command fails - otherwise the script will happily run through a failed command.</p>
  <p>You also may need to adjust the location of your output, and the directory the script is in upon execution, as you are running <code>GATK</code> from the current directory (<code>./gatk</code>)</p>
  <p>In order to make these long lines readable, I use the method of placing a <code>\</code> at the end of a line - this allows the command to wrap to the next line. Be sure that the <code>\</code> character is at the end of your line for this to function properly.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-21T03:09:41.763Z
Answers.Text: Wow, cool, thanks.
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-21T10:21:14.013Z
Answers.Text: OMG Thanks a million the script is running
Answers.UserDisplayName: Zizogolu
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-21T10:27:53.940Z
Answers.Text: Sorry @Scot, as the Mutect2 part of the script takes 15 hours in average for each sample, does this mean that for function 5, script will wait for the first sample to being completed by Mutect2 then goes for the second sample or Mutect2 will run of all of 18 samples in parallel? In the first case, the script will take 18 days?
Answers.UserDisplayName: Zizogolu
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-21T15:24:54.927Z
Answers.Text: Shell scripts are run line by line. If a line takes 15 hours, the interpreter will not move to the next line before that line is done. Loops are merely instructions to run a set of lines multiple times in series, so there's no parallel action happening here. You should ideally write a snakemake script with a rule for each GATK step then run a bunch of snakemake processes in parallel using whatever parallelization means are available to you.
Answers.UserDisplayName: Ram RS
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-21T17:54:03.227Z
Answers.Text: Hi @Angel - yes, the current version of this script runs sequentially, and is not parallelized. If function5 `Mutect2` takes a long time, then it may be worthwhile to pull it out of the loop, and run it upon completion of all other steps, and addition, explore if it can be run in parallel.
Answers.UserDisplayName: Scot
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-22T11:03:14.103Z
Answers.Text: You really don't want to use CAPS for shell variable names. Since, [by convention](https://stackoverflow.com/q/673055/1081936), global environment variables are capitalized, it is bad practice to also use caps for your own, local variables because this can lead to naming collisions and hard to find bugs.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2023-08-21T03:00:46.823Z
Answers.LastActivityDate: 2023-08-21T03:18:54.700Z
Answers.LastEditDate: 2023-08-21T03:18:54.700Z
Answers.UserDisplayName: Scot
Body: <p>I have some bam files in this directory</p>
  <pre><code>/data/Continuum/WES/results/
  </code></pre>
  <p>I want to run <code>GATK</code> mutation calling over bam files</p>
  <p>I googled and I realised for one function, I can do this</p>
  <pre><code>cd /data/Continuum/WES/vcf/
  
  for file in *.bam ; do ./gatk CollectSequencingArtifactMetrics -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup --FILE_EXTENSION .txt -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta done;
  </code></pre>
  <p>I have five functions to run but I do not know how to include all in the script</p>
  <p>These are my functions</p>
  <pre><code>./gatk CollectSequencingArtifactMetrics -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup --FILE_EXTENSION .txt -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
  
  ./gatk GetPileupSummaries -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.table -V af-only-gnomad.hg38.vcf.gz -L resources_broad_hg38_v0_wgs_calling_regions.hg38.interval_list -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
  
  ./gatk CalculateContamination -I NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.table -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.contamination.table
  
  ./gatk GetSampleName -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.sample_name
  
  ./gatk Mutect2 -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.vcf -tumor NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.sample_name --af-of-alleles-not-in-resource 2.5e-06 --germline-resource af-only-gnomad.hg38.vcf.gz -pon 1000g_pon.hg38.vcf.gz
  </code></pre>
  <p>I thought to use <code>|</code> to include five functions but I do not know how to change input and output for each iteration</p>
  <p>I tried this for the first function,</p>
  <pre><code>    FILES=&quot;/data/Continuum/WES/testAlignmentBROADGenome2/results/.MarkedDup.bam*&quot;
  for f in <span class="math-container">$FILES
  do
    ./gatk CollectSequencingArtifactMetrics -I $</span>f -O $f --FILE_EXTENSION .txt -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  done
  </code></pre>
  <p>but returned this error</p>
  <pre><code>SAMException: Cannot read non-existent file: file:///data/Continuum/WES/testAlignmentBROADGenome2/results/.MarkedDup.bam*
  </code></pre>
  <p>Anyone, please help me in develop this?</p>
  
Tags: <bash><gatk><wes><mutect2>
Title: A bash script for running on a bunch of bam files
ThreadText: # A bash script for running on a bunch of bam files
  
  (Asked by: username_0 on 2023-08-20)
  
  I have some bam files in this directory
  
  ```
  /data/Continuum/WES/results/
  
  ```
  
  I want to run `GATK` mutation calling over bam files
  
  I googled and I realised for one function, I can do this
  
  ```
  cd /data/Continuum/WES/vcf/
  
  for file in *.bam ; do ./gatk CollectSequencingArtifactMetrics -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup --FILE_EXTENSION .txt -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta done;
  
  ```
  
  I have five functions to run but I do not know how to include all in the script
  
  These are my functions
  
  ```
  ./gatk CollectSequencingArtifactMetrics -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup --FILE_EXTENSION .txt -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
  ./gatk GetPileupSummaries -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.table -V af-only-gnomad.hg38.vcf.gz -L resources_broad_hg38_v0_wgs_calling_regions.hg38.interval_list -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
  ./gatk CalculateContamination -I NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.table -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.contamination.table
  
  ./gatk GetSampleName -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.sample_name
  
  ./gatk Mutect2 -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta -I /data/Continuum/WES/results/NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.bam -O NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.vcf -tumor NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup.targeted_sequencing.sample_name --af-of-alleles-not-in-resource 2.5e-06 --germline-resource af-only-gnomad.hg38.vcf.gz -pon 1000g_pon.hg38.vcf.gz
  
  ```
  
  I thought to use `|` to include five functions but I do not know how to change input and output for each iteration
  
  I tried this for the first function,
  
  ```
      FILES="/data/Continuum/WES/testAlignmentBROADGenome2/results/.MarkedDup.bam*"
  for f in $FILES
  do
    ./gatk CollectSequencingArtifactMetrics -I $f -O $f --FILE_EXTENSION .txt -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  done
  
  ```
  
  but returned this error
  
  ```
  SAMException: Cannot read non-existent file: file:///data/Continuum/WES/testAlignmentBROADGenome2/results/.MarkedDup.bam*
  
  ```
  
  Anyone, please help me in develop this?
  
  **Tags:** bash, gatk, wes, mutect2
  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_1 on 2023-08-21. Score: 6 | Accepted answer
  
  There are a couple ways to approach this, one is to create output files based on the `BAM` input file name - the other approach is to make a new folder (based on the input `BAM` filename containing result files with the same name for each input. I generally use the first method, but either work.
  
  Here is a script that may run the 5 `GATK` functions. It calculates a prefix to be used for derivative output files `BASE`, and should run your pipeline on all the `BAM` files within the `INPUT_FOLDER`.
  
  ```
  #!/bin/bash
  
  set -e
  
  INPUT_FOLDER="/data/Continuum/WES/results/"
  
  for BAM_FILE in "$INPUT_FOLDER/"*.bam
  do
          BASE=$(basename "$BAM_FILE" .bam) # NG-27280_CLTSS_LTS_0017_lib506243_7661_2_MarkedDup
  
          # function 1
          ./gatk CollectSequencingArtifactMetrics \
                  -I "$BAM_FILE" \
                  -O "$BASE" \
                  --FILE_EXTENSION .txt \
                  -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
          # function 2
          ./gatk GetPileupSummaries \
                  -I "$BAM_FILE" \
                  -O "$BASE.targeted_sequencing.table" \
                  -V af-only-gnomad.hg38.vcf.gz \
                  -L resources_broad_hg38_v0_wgs_calling_regions.hg38.interval_list \
                  -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta
  
          # function 3
          ./gatk CalculateContamination \
                  -I "$BASE.targeted_sequencing.table" \
                  -O "$BASE.targeted_sequencing.contamination.table" \
  
          # function 4
          ./gatk GetSampleName \
                  -I "$BAM_FILE" \
                  -O "$BASE.targeted_sequencing.sample_name"
  
          # function 5
          ./gatk Mutect2 \
                  -R resources_broad_hg38_v0_Homo_sapiens_assembly38.fasta \
                  -I "$BAM_FILE" \
                  -O "$BASE".vcf \
                  -tumor "$BASE.targeted_sequencing.sample_name" \
                  --af-of-alleles-not-in-resource 2.5e-06 \
                  --germline-resource af-only-gnomad.hg38.vcf.gz \
                  -pon 1000g_pon.hg38.vcf.gz
  done
  
  ```
  
  It may need some adjustments - I did not test this, and used the `GATK` incantations you provided (thank you for including all of your steps!)
  
  It includes the addition of the command `set -e`, which depending on GATKs return status may stop the script if any command fails - otherwise the script will happily run through a failed command.
  
  You also may need to adjust the location of your output, and the directory the script is in upon execution, as you are running `GATK` from the current directory (`./gatk`)
  
  In order to make these long lines readable, I use the method of placing a `\` at the end of a line - this allows the command to wrap to the next line. Be sure that the `\` character is at the end of your line for this to function properly.
  
  > **Comments:**
  > **username_2**: Wow, cool, thanks.  
  > **username_0**: OMG Thanks a million the script is running  
  > **username_0**: Sorry @username_1, as the Mutect2 part of the script takes 15 hours in average for each sample, does this mean that for function 5, script will wait for the first sample to being completed by Mutect2 then goes for the second sample or Mutect2 will run of all of 18 samples in parallel? In the first case, the script will take 18 days?  
  > **username_3**: Shell scripts are run line by line. If a line takes 15 hours, the interpreter will not move to the next line before that line is done. Loops are merely instructions to run a set of lines multiple times in series, so there's no parallel action happening here. You should ideally write a snakemake script with a rule for each GATK step then run a bunch of snakemake processes in parallel using whatever parallelization means are available to you.  
  > **username_1**: Hi @Angel - yes, the current version of this script runs sequentially, and is not parallelized. If function5 `Mutect2` takes a long time, then it may be worthwhile to pull it out of the loop, and run it upon completion of all other steps, and addition, explore if it can be run in parallel.  
  > **username_4**: You really don't want to use CAPS for shell variable names. Since, [by convention](https://stackoverflow.com/q/673055/1081936), global environment variables are capitalized, it is bad practice to also use caps for your own, local variables because this can lead to naming collisions and hard to find bugs.  
  
  ---
  
  