Tutorial

Job Arrays

45 Minutes

Objectives

• Convert a batch script into a job array
• Run a whole pipeline, from preparation to a summary, from one submission script.
• Compare the pipeline's resource use with the job you started with.
• Recognise when a workflow manager would serve you better than batch scripts.

Prerequisite

• You can find why a job failed. See Debugging Failed Jobs.

What is a job array?

Job arrays allow you to use a single script to submit a set of nearly identical Slurm jobs. Each job in an array is allocated the same resources and submitted to the queue simultaneously. They are the solution to the problem known as embarrassingly parallel tasks. Embarrassingly parallel tasks are those which require no communication or dependency between tasks, meaning they can be run entirely independently of each other.

One job, many samples

In Job Dependencies, qc.sl and align.sl each work through every sample with a loop, one sample after another:

for infile in untrimmed_fastq/*_1.fastq.gz
do
    base=$(basename ${infile} _1.fastq.gz)
    # ... the commands for one sample ...
done

The samples don't depend on each other, so they could run at the same time.

Why not a job for each sample?

Another way would be to submit a separate job for each sample, in a loop. This works for three samples. What problems might it cause with hundreds of samples?

Solution

• Every sample is a separate job, with its own job ID, so the jobs are harder to watch or cancel together.
• Submitting many jobs one by one puts much more load on the scheduler than one job array.
• It is easy to lose track of which jobs belong to which run.

A job array runs the same script many times as one job, with one task for each index you give --array. Each task gets its index in the environment variable SLURM_ARRAY_TASK_ID, which the script uses to pick its input.

Converting the scripts

Instead of looping over the raw read files, each task of the array picks one of them. files=(untrimmed_fastq/*_1.fastq.gz) makes a bash array: a list of the matching files, in alphabetical order. ${files[0]} is the first file, ${files[1]} the second, and so on, so each task can take the file at its own index:

files=(untrimmed_fastq/*_1.fastq.gz)
base=$(basename ${files[SLURM_ARRAY_TASK_ID]} _1.fastq.gz)

Because the list starts at 0, the task indexes start at 0 too: for three samples, --array 0-2.

Here are the quality control and alignment scripts, with the loop replaced by those two lines:

qc.sl

#!/bin/bash -e

#SBATCH --job-name      qc
#SBATCH --account       nesi99991
#SBATCH --time          00:15:00
#SBATCH --mem           2G
#SBATCH --cpus-per-task 2
#SBATCH --output        log/%x_%A_%a.out

set -u

# Quality control for one sample: the raw read file at this array task's index.
files=(untrimmed_fastq/*_1.fastq.gz)
base=$(basename ${files[SLURM_ARRAY_TASK_ID]} _1.fastq.gz)

raw_1=untrimmed_fastq/${base}_1.fastq.gz
raw_2=untrimmed_fastq/${base}_2.fastq.gz
adapters=${EBROOTTRIMMOMATIC}/adapters/NexteraPE-PE.fa
mkdir -p results/fastqc trimmed

fastqc --outdir results/fastqc ${raw_1} ${raw_2}
trimmomatic PE ${raw_1} ${raw_2} \
    trimmed/${base}_1.trim.fastq.gz trimmed/${base}_1un.trim.fastq.gz \
    trimmed/${base}_2.trim.fastq.gz trimmed/${base}_2un.trim.fastq.gz \
    SLIDINGWINDOW:4:20 MINLEN:25 ILLUMINACLIP:${adapters}:2:40:15

align.sl

#!/bin/bash -e

#SBATCH --job-name      align
#SBATCH --account       nesi99991
#SBATCH --time          00:15:00
#SBATCH --mem           1500M
#SBATCH --cpus-per-task 4
#SBATCH --output        log/%x_%A_%a.out

set -u

# Alignment and variant calling for one sample.
# Uses the same list of raw read files as qc.sl,
# so each task index is the same sample in both.
files=(untrimmed_fastq/*_1.fastq.gz)
base=$(basename ${files[SLURM_ARRAY_TASK_ID]} _1.fastq.gz)

ref_genome=ref_genome/ecoli_rel606.fasta
fq1=trimmed/${base}_1.trim.fastq
fq2=trimmed/${base}_2.trim.fastq
sam=results/sam/${base}.aligned.sam
bam=results/bam/${base}.aligned.bam
sorted_bam=results/bam/${base}.aligned.sorted.bam
bcf=results/bcf/${base}.bcf
vcf=results/vcf/${base}.vcf
final_variants=results/vcf/${base}.var.vcf

mkdir -p results/sam results/bam results/bcf results/vcf
gunzip -kf ${fq1}.gz ${fq2}.gz

bwa-mem2 mem -t ${SLURM_CPUS_PER_TASK} ${ref_genome} ${fq1} ${fq2} > ${sam}
samtools view -S -b ${sam} > ${bam}
samtools sort ${bam} -o ${sorted_bam}
bcftools mpileup -O b -o ${bcf} -f ${ref_genome} ${sorted_bam}
bcftools call --ploidy 1 -m -v -o ${vcf} ${bcf}
vcfutils.pl varFilter ${vcf} > ${final_variants}

align.sl builds the list from the raw read files too, not the trimmed ones, so each task index is the same sample in both scripts. Indexing (index.sl) runs once for all samples, so it stays as it is.

#

Output filenames

A job array script has #SBATCH --output log/sample_${SLURM_ARRAY_TASK_ID}.out. Where do the tasks write their output?

All of them write to one file, literally named log/sample_${SLURM_ARRAY_TASK_ID}.out.

Slurm reads #SBATCH lines before the job starts, so the variable is not replaced. Use %a instead.

Task 1 writes to log/sample_1.out, task 2 to log/sample_2.out, and so on.

That is what log/sample_%a.out does.

The job fails to submit.

Running the array

Solo learners

Replace nesi99991 in the --account line with your project code.

Submit a job array with --array, giving one task for each raw read file. Use %A (the job ID of the whole array) and %a (the task index) in the log filenames, so each task writes its own:

sbatch --array 0-2 --output log/%x_%A_%a.out qc.sl

squeue --me shows the tasks still waiting together on one line, such as JOB_ID_[0-2], and each running task as JOB_ID_0, JOB_ID_1 and so on. To check one task, give seff the task's ID, for example seff JOB_ID_1. Replace JOB_ID with the job ID that sbatch printed.

To limit how many tasks run at the same time, add % and a number to --array. --array 0-2%2 runs at most two tasks at once.

Don't make each task too small

Each task queues and starts like a separate job, and counts towards the queue limit in Job Limits. If each piece of work takes only a few seconds, group several pieces into each task rather than giving each its own.

Putting it together

The job array was the last piece. We can now queue the whole workflow, from quality control to a summary, with one command.

The finished pipeline

This is the final version of pipeline.sh, building on the one from Job Dependencies:

#!/bin/bash -e

set -u

# Load the software here. sbatch passes this environment on to each job.
module purge
module load FastQC/0.12.1
module load Trimmomatic/0.39-Java-1.8.0_144
module load bwa-mem2/2.3-GCC-12.3.0
module load SAMtools/1.22-GCC-12.3.0
module load BCFtools/1.22-GCC-12.3.0

# Limit Java's memory, or Trimmomatic sizes it from the whole node.
export JAVA_TOOL_OPTIONS="-Xmx1g"

# One array task per sample. Task indexes start at 0, so the last is one less
# than the number of files.
files=(untrimmed_fastq/*_1.fastq.gz)
last=$(( ${#files[@]} - 1 ))

# Quality control as a job array, and indexing alongside.
qc=$(sbatch --parsable --array 0-${last} qc.sl)
index=$(sbatch --parsable index.sl)

# Each alignment task waits for the QC task of the same sample, and for indexing.
align=$(sbatch --parsable --array 0-${last} \
    --dependency aftercorr:${qc},afterok:${index} align.sl)

# Once every alignment task has succeeded, summarise the results.
summary=$(sbatch --parsable --dependency afterok:${align} summarise.sl)

# If any alignment task fails, resubmit only the failed tasks, with more memory.
resubmit=$(sbatch --parsable --dependency afternotok:${align} --job-name resubmit \
    --time 00:05:00 --mem 512M --output log/%x_%j.out \
    --wrap "bash resubmit_failed.sh ${align}")

echo "qc: ${qc}  index: ${index}  align: ${align}"
echo "summarise: ${summary}  resubmit: ${resubmit}"

For the three samples, it queues two job arrays and three single jobs. Each alignment task waits for the quality control task with the same index and for indexing. The summary runs if every alignment task succeeds, and the resubmit job runs if any fails:

flowchart TB
    subgraph qc["qc.sl job array"]
        qc1["Task 0<br>SRR2584863"]
        qc2["Task 1<br>SRR2584866"]
        qc3["Task 2<br>SRR2589044"]
    end
    index["index.sl"]
    subgraph align["align.sl job array"]
        align1["Task 0<br>SRR2584863"]
        align2["Task 1<br>SRR2584866"]
        align3["Task 2<br>SRR2589044"]
    end
    qc1 -- aftercorr --> align1
    qc2 -- aftercorr --> align2
    qc3 -- aftercorr --> align3
    index -- afterok --> align
    align -- afterok --> summarise["summarise.sl"]
    align -. afternotok .-> resubmit["resubmit_failed.sh"]

• Quality control (qc.sl) runs as a job array, one task per sample, and indexing (index.sl) runs at the same time.
• Alignment (align.sl) is a second job array. aftercorr starts each of its tasks as soon as the quality control task with the same index has succeeded, and afterok waits for indexing. The , means both conditions must be met.

• summarise.sl waits for every alignment task to succeed, then counts the variants found in each sample. With afterok, a job that depends on a job array waits for every task:

  #!/bin/bash -e
  
  #SBATCH --job-name      summarise
  #SBATCH --account       nesi99991
  #SBATCH --time          00:05:00
  #SBATCH --mem           512M
  #SBATCH --output        log/%x_%j.out
  
  summary=results/variant_counts.txt
  
  # Count the variants found in each sample.
  echo -e "sample\tvariants" > ${summary}
  for vcf in results/vcf/*.var.vcf
  do
      sample=$(basename ${vcf} .var.vcf)
      count=$(grep -vc '^#' ${vcf})
      echo -e "${sample}\t${count}" >> ${summary}
  done
  
  cat ${summary}
  

• If any alignment task fails, the resubmit job runs instead. For a job array, afternotok waits until every task has finished and at least one has failed. It runs resubmit_failed.sh, which finds the failed tasks and submits only those again, with more memory:

  #!/bin/bash -e
  
  # align.sl uses the software loaded here.
  module purge
  module load FastQC/0.12.1
  module load Trimmomatic/0.39-Java-1.8.0_144
  module load bwa-mem2/2.3-GCC-12.3.0
  module load SAMtools/1.22-GCC-12.3.0
  module load BCFtools/1.22-GCC-12.3.0
  
  # Limit Java's memory, or Trimmomatic sizes it from the whole node.
  export JAVA_TOOL_OPTIONS="-Xmx1g"
  
  # Find the failed tasks of the job array with ID $1, as a list like "1,3".
  failed=$(sacct --jobs "$1" --array --allocations --noheader \
      --parsable2 --format JobID,State |
      grep -E 'FAILED|OUT_OF_MEMORY|TIMEOUT' |
      cut -d '|' -f 1 | cut -d '_' -f 2 | paste -sd ,)
  
  # Resubmit only those tasks, with more memory.
  if [ -n "${failed}" ]; then
      sbatch --array "${failed}" --mem 3G \
          align.sl
  fi
  

--wrap submits a single command as a job, without a batch script. If tasks are resubmitted, summarise.sl is cancelled, because its dependency can no longer be met, so run it yourself once they finish.

Save the scripts next to each other, then run bash pipeline.sh on the login node.

Compare with where you started

"CPU time booked" is the time each job ran multiplied by the CPU cores it had, which is what counts against your project.

Where did the savings come from?

Which change, from which lesson, made the biggest difference?

Solution

• Splitting the job (Job Dependencies) let each step request only what it needs, instead of holding four CPUs, 5 GB and four hours for the whole run.
• The job array (this lesson) runs the samples at the same time, so the pipeline finishes sooner.
• Passing ${SLURM_CPUS_PER_TASK} to the tools (Using Variables) lets alignment use the CPUs it asks for.

When to use a workflow manager

A submission script works well for a few steps. As a workflow grows, a workflow manager can take over: you describe each step with its inputs and outputs, and it works out the order, submits the Slurm jobs for you, and can skip the steps whose results are already up to date. Snakemake and Nextflow are both available on Mahuika.

#

Choosing a tool

Which of these suggests it is time to move to a workflow manager?

You often need to rerun only the stages whose inputs changed.

Your submission script has grown to many stages, with different dependencies between them.

You need to run the same command on 50 input files.

A single job array handles that well.

Submit

Keypoints

• A short submission script can queue a whole pipeline, with each stage sized for its own needs.
• Measure before and after: the same work can use far fewer resources.
• Workflow managers take over when pipelines outgrow scripts.

What next

• Snakemake and Nextflow on Mahuika.
• For reference: Job Dependencies and Job Arrays.

Quiz Progress

• Answered: 0 / 0 (0%)
• Correct: 0 / 0 (0%) Reset