6
$\begingroup$

I am trying to map single-cell chromatin accessibility data from doi:10.1038/nature14590, obtained using scATAC-seq, to the reference genome. Example paired-end reads are here.

The reads contain various adapter sequences, which the authors "trimmed from FASTQs using custom python scripts" [ref]. I am having problems reproducing this trimming. They used Nextera Transposase Adapters

TCGTCGGCAGCGTCAGATGTGTATAAGAGACAG GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAG

and custom barcode adapters which are listed [here]. An excerpt from this list:

Custom Barcodes Adapter 1 (index i5): AATGATACGGCGACCACCGAGATCTACACTAGATCGCTCGTCGGCAGCGTCAGATGTGTAT AATGATACGGCGACCACCGAGATCTACACTAACCAAGTCGTCGGCAGCGTCAGATGTGTAT AATGATACGGCGACCACCGAGATCTACACAGGTTGCCTCGTCGGCAGCGTCAGATGTGTAT Custom Barcodes Adapter 2 (index i7): CAAGCAGAAGACGGCATACGAGATACCCAGCAGTCTCGTGGGCTCGGAGATGTG CAAGCAGAAGACGGCATACGAGATCCCAACCTGTCTCGTGGGCTCGGAGATGTG CAAGCAGAAGACGGCATACGAGATCACCACACGTCTCGTGGGCTCGGAGATGTG

Analyzing the FASTA-files with FastQC showed that the reverse complement of those sequences was over-represented in the data and also confirmed Nextera Transposase Adapters.

I tried to trim the sequences using cutadapt:

cutadapt --minimum-length=20 -a CACATCTCCGAGCCCACGAGAC -a CTGTCTCTTATA -A ATACACATCTGACGCTGCCGACGA -A CTGTCTCTTATA -o SRR1780164_1_trimmed.fastq.gz -p SRR1780164_2_trimmed.fastq.gz SRR1780164_1.fastq.gz SRR1780164_2.fastq.gz

with

  • CACATCTCCGAGCCCACGAGAC being the reverse complement of the 2nd fix part of the i7 adapters
  • ATACACATCTGACGCTGCCGACGA being the reverse complement of the 2nd fix part of the i5 adapters
  • CTGTCTCTTATA for the Nextera Transposase Adapters

FastQC does not indicate overrepresented sequences in the resulting FASTA-files, however still indicates Nextera Transposase Adapters in the _1 files. On the other hand, after TrimGalore removes CTGTCTCTTATA, Nextra Adapters are not indicated anymore, for some reason.

Furthermore, mapping the reads using bowtie2

bowtie2 -x GRCh38_primary_assembly_index -1 SRR1780164_1_trimmed.fastq.gz -2 SRR1780164_2_trimmed.fastq.gz > SRR1780164.sam

results in a low alignment rate

1352 (2.79%) aligned concordantly exactly 1 time
3420 (7.07%) aligned concordantly >1 times
11.25% overall alignment rate

which gives me the impression I did something wrong during adapter trimming.

Any help would be much appreciated.

--

Edit: As I am not sure about the sequences to use for trimming, I also tried -a CTGTCTCTTATACACATCTCCGAGCCCACGAGAC -A CTGTCTCTTATACACATCTGACGCTGCCGACGA, as found in a suggestion for Nextera Adapters, however it did not remove the sequence over-representation in the -1 file. -a CACATCTCCGAGCCCACGAGAC -a CTGTCTCTTATACACATCTCCGAGCCCACGAGAC -A CTGTCTCTTATACACATCTGACGCTGCCGACGA solves adapter-issues according to FastQC, but also does not give good mapping results.

$\endgroup$
4
$\begingroup$

First the Nextera adapters and the custom barcode adapters overlap each other

Nextera 1                                     TCGTCGGCAGCGTCAGATGTGTATAAGAGACAG 
Barcode 1   GATACGGCGACCACCGAGATCTACACTAGATCGCTCGTCGGCAGCGTCAGATGTGTAT 

Nextera 2                                    GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAG
Barcode 2    CAAGCAGAAGACGGCATACGAGATACCCAGCAGTCTCGTGGGCTCGGAGATGTG 

Looking at some read 1s and read 2s you were correct that the adapters appear at the 3' end of the read reverse complimented (The |> points to the adapter in the sequence). You can see that CTGTCTCTTATACACA would be the ideal sequence to use for the 3' end of the reads. Since the Nextera adapter and the barcode overlap you only need to use the one that is adjacent to the insert.

SRR1780164_2.fastq: GCCCGAGGGGACGCTGCCCCACGTGCTCCAGG |> CTGTCTCTTATACACATCTGACGCTGCCGACGATGCACGAAGTG
SRR1780164_1.fastq: GGGCTGGACAGCGGCGGGGACGCGGACCTGCAG |> CTGTCTCTTATACACATCTCCGAGCCCACGAGACACAAAGTGA

Of the 64k reads there are ~11k with the Nextera adapter at the 3' end based on a rough sequence grep

On the 5' end of the read there is some Nextera adapter present but it is not reverse complimented and there are only a handful of reads with the adapter present ~250 total between read 1 and read 2.

SRR1780164_2.fastq: TGGGCTCGGAGATGTGTATAAGAGACAG <| CTCCCCCGCCTCCAGCATCCGGGCGAGGTAGTGCATCGACGCGTCCAC

I used this command:

cutadapt --minimum-length=20 -A CTGTCTCTTATACACA -a CTGTCTCTTATACACA -G GATGTGTATAAGAGACAG 
-g GATGTGTATAAGAGACAG -o SRR1780164_1_trimmed.fastq.gz -p SRR1780164_2_trimmed.fastq.gz 
SRR1780164_1.fastq SRR1780164_2.fastq

Which gives a trimming result of

   === Summary ===

Total read pairs processed:             64,832
  Read 1 with adapter:                  16,721 (25.8%)
  Read 2 with adapter:                  16,996 (26.2%)
Pairs that were too short:                  46 (0.1%)
Pairs written (passing filters):        64,786 (99.9%)

I think this best represents the library construct and the trimming of the adapter sequences

| improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.