Tutorial

Job Dependencies

45 Minutes

Objectives

• Split a job into parts that each request only the resources they need.
• Make a job wait for another job with --dependency.
• Queue a whole workflow from one submission script.

Prerequisite

• You will need pipeline.sl as set up in Using Variables.

One script per stage

The steps of pipeline.sl have different resource needs:

1. Indexing the reference genome
2. Quality control
3. Alignment and variant calling

Even so, pipeline.sl runs every step, for every sample, in one job that holds 4 CPUs, 5 GB and 4 hours the whole time. Before we can give each step its own resources, each step needs to be a script of its own.

Move quality control into its own script

Move the quality control part of pipeline.sl, with its loop and the variables it uses, into a new script, qc.sh.

Solution

qc.sh is the quality control section of pipeline.sl, with adapters moved in from the top:

#!/bin/bash -e

set -u

# Quality control for every sample.
adapters=${EBROOTTRIMMOMATIC}/adapters/NexteraPE-PE.fa
mkdir -p results/fastqc trimmed

for infile in untrimmed_fastq/*_1.fastq.gz
do
    base=$(basename ${infile} _1.fastq.gz)
    raw_2=untrimmed_fastq/${base}_2.fastq.gz
    fastqc --outdir results/fastqc ${infile} ${raw_2}
    trimmomatic PE ${infile} ${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
done

Do the same for indexing and alignment

Move reference genome indexing into index.sh, and alignment and variant calling into align.sh.

Solution

index.sh indexes the reference genome:

#!/bin/bash -e

set -u

# Reference genome indexing, once for all samples.
ref_genome=ref_genome/ecoli_rel606.fasta

gunzip -kf ${ref_genome}.gz
bwa-mem2 index ${ref_genome}

align.sh also creates the results folders, which pipeline.sl did between the steps. It unzips each sample's trimmed reads inside the loop, so it loops over the zipped files:

#!/bin/bash -e

set -u

# Alignment and variant calling for every sample.
ref_genome=ref_genome/ecoli_rel606.fasta
mkdir -p results/sam results/bam results/bcf results/vcf

for infile in trimmed/*_1.trim.fastq.gz
do
    base=$(basename ${infile} _1.trim.fastq.gz)

    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

    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}
done

The step scripts don't load any modules: run with bash from pipeline.sl, they use the modules it has already loaded.

Running the steps in one job

pipeline.sl still runs the same steps in the same order, in one job, but now calls a script for each:

#!/bin/bash -e

#SBATCH --job-name      pipeline
#SBATCH --account       nesi99991
#SBATCH --time          04:00:00
#SBATCH --mem           5G
#SBATCH --cpus-per-task 4
#SBATCH --profile       task
#SBATCH --output        log/%x_%j.out

set -u

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"

# The step scripts run with the modules loaded here.
bash index.sh
bash qc.sh
bash align.sh

It runs index.sh, then qc.sh, then align.sh. Quality control and alignment each work through every sample:

flowchart TD
    subgraph job["pipeline.sl"]
        index["index.sh"] --> qc["qc.sh"]
        qc --> align["align.sh"]
    end

Split by resource needs

Each step is now a script of its own, but pipeline.sl still runs them all in one job, holding 4 CPUs, 5 GB and 4 hours for every step.

Splitting tasks by resource needs

Let's give each step script (qc.sh, index.sh and align.sh) an appropriate resource request, based on the profile from Evaluating Job Efficiency, so each can run as a separate job.

Once a script has an #SBATCH header, it is a batch script, so rename it to end in .sl:

mv qc.sh qc.sl
mv index.sh index.sl
mv align.sh align.sl

Solution: updated scripts

There is no one correct answer here! But here is one option. The profile plot in Evaluating Job Efficiency shows the steps of the job: quality control first, a short blip that is probably indexing the reference genome, then alignment and variant calling.

So we can give the three scripts better resource requests. Each still works through every sample, and writes its output to log/, like pipeline.sl since Using Variables.

Quality control: qc.sl requests fewer CPUs, less memory and much less time.

#!/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_%j.out

set -u

# Quality control for every sample.
adapters=${EBROOTTRIMMOMATIC}/adapters/NexteraPE-PE.fa
mkdir -p results/fastqc trimmed

for infile in untrimmed_fastq/*_1.fastq.gz
do
    base=$(basename ${infile} _1.fastq.gz)
    raw_2=untrimmed_fastq/${base}_2.fastq.gz
    fastqc --outdir results/fastqc ${infile} ${raw_2}
    trimmomatic PE ${infile} ${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
done

Reference genome indexing: index.sl requests less CPUs, memory, and time.

#!/bin/bash -e

#SBATCH --job-name      index
#SBATCH --account       nesi99991
#SBATCH --time          00:15:00
#SBATCH --mem           1G
#SBATCH --cpus-per-task 1
#SBATCH --output        log/%x_%j.out

set -u

# Reference genome indexing, once for all samples.
ref_genome=ref_genome/ecoli_rel606.fasta

gunzip -kf ${ref_genome}.gz
bwa-mem2 index ${ref_genome}

Alignment and variant calling: align.sl requests the same CPUs but less memory and time.

#!/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_%j.out

set -u

# Alignment and variant calling for every sample.
ref_genome=ref_genome/ecoli_rel606.fasta
mkdir -p results/sam results/bam results/bcf results/vcf

for infile in trimmed/*_1.trim.fastq.gz
do
    base=$(basename ${infile} _1.trim.fastq.gz)

    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

    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}
done

Solo learners

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

Quality control requests two CPUs, indexing one, and alignment four. Indexing doesn't need the results of quality control either, so those two jobs can run at the same time.

Once each script has run, check it with seff and adjust its request, leaving 20 to 30% spare as in the Introduction to HPC.

Make a job wait for another

Alignment can only start once quality control and indexing have finished. Rather than waiting to submit align.sl by hand, we can submit it straight away and tell Slurm to wait, with --dependency:

sbatch --dependency afterok:JOB_ID align.sl

Replace JOB_ID with the job ID of the job to wait for. To wait for more than one job, separate the job IDs with :, for example afterok:JOB_ID_1:JOB_ID_2.

The condition decides what the job waits for:

Condition	The job can start when the listed jobs have
afterok	finished successfully.
afternotok	finished in a failed state.
afterany	finished, in any state.
after	started.

If a job you depend on fails, a job waiting with afterok can never run, and Mahuika cancels it.

#

Dependencies

Job 1001 is running and job 1002 has failed. Which of these jobs can start once their dependency is checked?

--dependency after:1001

after only needs job 1001 to have started.

--dependency afterok:1002

Job 1002 failed, so this job can never start. Mahuika cancels it.

--dependency afternotok:1002

afternotok waits for the job to fail, which it has.

--dependency afterok:1001

Job 1001 is still running. This job waits until it finishes successfully.

Submit

Queue a whole workflow

pipeline.sl already runs a script for each step. If it submits each step with sbatch instead of running it with bash, every step becomes a job of its own. sbatch --parsable prints only the job ID, so the script can save it and use it in the next submission:

#!/bin/bash -e

set -u

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"

# sbatch passes the modules loaded here on to each job,
# and each step script's #SBATCH header sets its resources.
qc=$(sbatch --parsable qc.sl)
index=$(sbatch --parsable index.sl)

# Alignment waits for quality control and indexing to succeed.
align=$(sbatch --parsable --dependency afterok:${qc}:${index} align.sl)

echo "qc: ${qc}  index: ${index}  align: ${align}"

It queues three jobs, each with its own resources. Quality control and indexing can start straight away. Alignment waits, with afterok, for both to succeed:

flowchart TB
    qc["qc.sl<br>all samples"] --> align["align.sl<br>all samples"]
    index["index.sl"] --> align

The script is no longer a batch job itself: it has no #SBATCH header, and only submits the other jobs, so it finishes straight away. As it is no longer a batch script, save it as pipeline.sh, not pipeline.sl, next to the three step scripts. Run it on the login node with bash pipeline.sh. squeue --me then shows the alignment job waiting, with the reason (Dependency), until the jobs they depend on have finished.

Why queue every step at once?

What is gained by submitting the later steps straight away, rather than waiting to submit them by hand?

Solution

Each step can start as soon as the steps it depends on have succeeded, without you being there to submit it. The whole workflow runs from one command.

Quality control and alignment still work through the samples one after another, inside one job each. The samples don't depend on each other, so the next lesson runs them at the same time, as a job array.

Keypoints

• A script for each step lets each step run as its own job.
• Split a job where its parts need different resources.
• --dependency afterok:JOB_ID makes a job wait for another to finish successfully.
• sbatch --parsable gives the job ID, so a script can queue a whole workflow.

What next

• Next lesson: Job Arrays.
• For reference: Job Dependencies.

Quiz Progress

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