Skip to main content

STEP 3/4 Alignments

This page covers the alignment and post-processing steps for ATAC-seq and ChIP-seq data - mapping sequencing reads to a reference genome.

STEP 3: Alignment

The most commonly used aligners for ATAC-seq and ChIP-seq are:

Bowtie2

(for ChIP-seq, ATAC-seq, and other DNA-seq data)

Bowtie2 is a fast and memory-efficient aligner suitable for short reads. It requires an indexed reference genome to perform the alignment.

Note: The alignment steps often take 2-8 hours depending on the size of the FASTQ files and the computational resources available.

Loading Bowtie2

Bowtie2 is loaded on the CCB with:

module load bowtie2

Indexed Reference Genomes

The CCB provides pre-built Bowtie2 indices for common reference genomes. You can find the available indices in the directory /databank/igenomes/.

The reference genome index for human (hg38): /databank/igenomes/Homo_sapiens/UCSC/hg38/Sequence/Bowtie2Index/genome

The reference genome index for human (mm39): /databank/igenomes/Mus_musculus/UCSC/mm39/Sequence/Bowtie2Index/genome

If you have a custom reference, a Bowtie2 indexed genome can be generated from a genome-wide fasta file:

bowtie2-build -f custom_genome.fa custom_genome_index

Running Bowtie2

Common Bowtie2 Parameters

  • -x - Index basename
  • -U - Unpaired reads
  • -1 / -2 - Paired-end reads (forward/reverse)
  • -S - Output SAM file
  • -p - Number of threads (for parallel processing)
  • --very-sensitive - Slower but more accurate alignment mode

We add special cases to Bowtie2:

  • -k 2 = allow reads which multi-map to max 2 regions to continue into downstream analysis, very important for Hba-a1/2. Change -k to alter this.
  • -N 1 = allow maximum of one mismatch in seed. Change -N to alter this.
    • Allowing fewer mismatches will increase accuracy;
    • Allowing more will increase sensitivity.
  • --maxins 1000 = maximum fragment length allowed is 1000bp. Change –maxins if you want to alter this
    • perhaps change to 150bp if you want only mononucleosomes

Paired-end reads:

To do a paired-read alignment to mm39:

bowtie2 -k 2 -N 1 -p 4 \
-1 sample_R1.fastq.gz -2 sample_R2.fastq.gz \
--maxins 1000 \
-x /databank/igenomes/Mus_musculus/UCSC/mm39/Sequence/Bowtie2Index/genome \
-S sample_alignment.sam

(Rare) Single-end reads:

There may be occasions when you are analysing single-end ATAC-seq or ChIP-seq data:

bowtie2 -k 2 -N 1 -p 4 \
-x /databank/igenomes/Mus_musculus/UCSC/mm39/Sequence/Bowtie2Index/genome \
-U sample.fastq.gz -S sample.sam

Alterntive to Bowtie: BWA

BWA (Burrows-Wheeler Aligner) is another widely used aligner for short reads. It is particularly effective for mapping reads to large reference genomes.

BWA is loaded on the CCB with:

module load bwa

STEP 4 - Post-alignment Processing

After alignment, the SAM file needs to be converted to BAM format and processed.

Samtools

Samtools is a suite of programs for interacting with high-throughput sequencing data in SAM and BAM formats.

Loading Samtools

On the CCB:

module load samtools

Keep only Mapped Reads

Convert SAM to BAM, keep only mapped reads:

samtools view -@8 -bS -F4 -o mapped.bam alignment.sam ; # Convert SAM to BAM, keep only mapped reads

Sort BAM File

Sort the BAM file by genomic coordinates:

samtools sort mapped.bam -@8 -o mapped.srtC

Remove Duplicates

Remove PCR duplicates from the sorted BAM file:

samtools rmdup mapped.srtC filtered.bam ; # Remove PCR duplicates

Index BAM File

Index the final BAM file for efficient access:

samtools index filtered.bam ; # Index the final BAM File