<p>Dear BioStars Community,</p>

<p>I have a problem with an alignment to the transcriptome.</p>

<p>I have 8 RNA-Seq libraries sequenced on Illumina HiScanSQ system in one lane (2x100bp, paired-end) per sequencing run. These 8 libraries (1 pool) were put into two sequencing runs to obtain a decent number of reads. After demultiplexing (using bcl2fastq-1.8.4) the reads were trimmed using TrimGalore and aligned to the previously assembled transcriptome (because there is no reference genome for the organism - <em>Pinus sylvestris</em> - I am trying to analyze...) by Bowtie2-2.2.6. In the case of 7 libraries there were almost no difference in the alignment efficiency (~85-95%, with ~60-85% of uniquely mapped reads), but in case of one library something strange happened:</p>

<p><strong>Run1:</strong></p>

<pre>
11109859 reads; of these:
  11109859 (100.00%) were paired; of these:
    1701658 (15.32%) aligned concordantly 0 times
    6666961 (60.01%) aligned concordantly exactly 1 time
    2741240 (24.67%) aligned concordantly &gt;1 times
    ----
    1701658 pairs aligned concordantly 0 times; of these:
      11078 (0.65%) aligned discordantly 1 time
    ----
    1690580 pairs aligned 0 times concordantly or discordantly; of these:
      3381160 mates make up the pairs; of these:
        3218192 (95.18%) aligned 0 times
        86430 (2.56%) aligned exactly 1 time
        76538 (2.26%) aligned &gt;1 times
85.52% overall alignment rate</pre>

<p><strong>Run2</strong>:</p>

<pre>
14719563 reads; of these:
  14719563 (100.00%) were paired; of these:
    7641835 (51.92%) aligned concordantly 0 times
    4991995 (33.91%) aligned concordantly exactly 1 time
    2085733 (14.17%) aligned concordantly &gt;1 times
    ----
    7641835 pairs aligned concordantly 0 times; of these:
      7874 (0.10%) aligned discordantly 1 time
    ----
    7633961 pairs aligned 0 times concordantly or discordantly; of these:
      15267922 mates make up the pairs; of these:
        15039673 (98.51%) aligned 0 times
        94443 (0.62%) aligned exactly 1 time
        133806 (0.88%) aligned &gt;1 times
48.91% overall alignment rate</pre>

<p>So my question is: <strong>what should I do to find out what went wrong?</strong> I excluded (maybe too soon...) the human error because this was the same pool used for two runs (from one Eppendorf tube).</p>

<p>I also did FastQC on demultiplexed and trimmed reads - links for this library with low alignment efficiency are provided below:</p>

<p>Run1 (the good one), demultiplexed: <a href="http://twrzes.wtvk.pl/run1_R1_fastqc.html">http://twrzes.wtvk.pl/run1_R1_fastqc.html </a>and <a href="http://twrzes.wtvk.pl/run1_R2_fastqc.html">http://twrzes.wtvk.pl/run1_R2_fastqc.html</a></p>

<p>Run1, after trimming: <a href="http://twrzes.wtvk.pl/run1_R1_trimmed_fastqc.html">http://twrzes.wtvk.pl/run1_R1_trimmed_fastqc.html</a> and <a href="http://twrzes.wtvk.pl/run1_R2_trimmed_fastqc.html">http://twrzes.wtvk.pl/run1_R2_trimmed_fastqc.html</a></p>

<p>Run2 (the bad one), demultiplexed: <a href="http://twrzes.wtvk.pl/run2_R1_fastqc.html">http://twrzes.wtvk.pl/run2_R1_fastqc.html</a> and <a href="http://twrzes.wtvk.pl/run2_R2_fastqc.html">http://twrzes.wtvk.pl/run2_R2_fastqc.html</a></p>

<p>Run2, after trimming: <a href="http://twrzes.wtvk.pl/run2_R1_trimmed_fastqc.html">http://twrzes.wtvk.pl/run2_R1_trimmed_fastqc.html</a> and <a href="http://twrzes.wtvk.pl/run2_R2_trimmed_fastqc.html">http://twrzes.wtvk.pl/run2_R2_trimmed_fastqc.html</a></p>

<p><strong>Command-line commands I used for:</strong></p>

<p><strong>1) Demultiplexing:</strong></p>

<pre>
/path/to/configureBclToFastq.pl --input-dir /path/to/folder/with/BCLs/Data/Intensities/BaseCalls --output-dir /path/to/folder/with/BCLs/Unaligned --sample-sheet /path/to/folder/with/BCLs/sample-sheet.csv --fastq-cluster-count 0 --mismatches 1 --with-failed-reads</pre>

<p><strong>2) Trimming (TrimGalore-0.4.0, a wrapper for cutadapt-1.8.3):</strong></p>

<pre>
trimgalore --paired --quality 20 --illumina --stringency 1 -e 0.2 --length 40 -o /path/to/trimmed/fastq --trim1 run1_R1.fastq run1_R2.fastq</pre>

<p><strong>3) Alignment (Bowtie2-2.2.6)</strong></p>

<pre>
bowtie2 -p 12 -I 0 -X 2000 --dovetail --very-sensitive-local -N 1 -x /path/to/index/index -1 run1_R1_trimmed.fastq -2 run1_R2_trimmed.fastq -S /path/to/aligned/sam/run1.sam</pre>

<p>If you need any additional info, I would be more than happy to provide it.</p>

<p>Thank you very much for your efforts on solving this problem.</p>

<p>Kind regards,<br />
Tomasz Wrzesinski</p>

<p>--<br />
Tomasz Wrzesinski, MSc<br />
PhD Student<br />
Laboratory of High Throughput Technologies<br />
Institute of Molecular Biology and Biotechnology<br />
Faculty of Biology<br />
Adam Mickiewicz University in Poznan<br />
Umultowska 89/1.117<br />
61-614 Poznan, Poland<br />
tel. +48 61 829 5833<br />
e-mail: twrzes@amu.edu.pl</p>

<p>My first approach would be assemble the run2 unmapped reads and blast the contigs, looking for contaminants. Alternatively, you could blast a sample of the unmapped raw reads. It may be interesting to assemble and map run1 unmapped reads as well, as a control.</p>

<p>Are you using a stranded or unstranded library preparation protocol?</p>

<p>Looking at the FastQC reports, both runs seems just fine. The only suspicious thing I noticed is GC content seems to be slightly different (1%-2%) between runs, and on run2 %A is consistently higher than %T.</p>
