<p>If you want to keep track of the unmapped reads, I recommend you to use <a href="http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml" rel="nofollow noreferrer">Bowtie2</a> instead of bwa.
Start by building an index:</p>
<pre><code>bowtie2-build genome.fa genome.fa
</code></pre>
<p>Then map your reads, using arguments to save unmapped reads in a specific file like this:</p>
<pre><code>bowtie2 -x genome.fa -1 reads_1.fq -2 reads_2.fq -S bowtie2_mapping.sam --un-conc-gz unmappedreads_%.fq.gz
</code></pre>
<p>It will write the alignment as well as unmappedreads_1.fq.gz and unmappedreads_2.fq.gz that you may use in DI-tector.</p>

CC BY-SA 4.0
2020-06-19T11:59:10.800Z
My thoughts exactly. Personally I don't know why a paper in 2018 would want to use BWA over Bowtie2. Any ideas?
M__
CC BY-SA 4.0
2020-06-19T12:30:13.193Z
after putting a glance at the article, it is no really clear how they switch from bwa mem to bwa samse to align reads. Maybe it has something to see with the flags in the output .sam file, required for DI-tector proper running. But it is only a supposition...
thomas duge de bernonville
CC BY-SA 4.0
2020-06-19T16:16:53.087Z
Thank you. But what about `samtools` step? I mean, my read has to be aligned with bwa (or bowtie2 - great I will try that) *FIRSTLY* with the human genome, then *SECONDLY* with virus genome and only after these two steps it has to be filtered through `SAMTOOLS` and only after that I have to have the unmapped reads that I can put into a DI-tector.
Dmitrii Trubetskoy
CC BY-SA 4.0
2020-06-19T19:59:10.137Z
I think it depends on the format of input data required for DI-tector. Bowtie2 will return .fastq files. In the article, they used samtools with -Bs -f 4 options probably to return unaligned reads into a .sam with specific CIGAR labels.
thomas duge de bernonville
CC BY-SA 4.0
2020-06-19T20:25:22.247Z
@thomas duge de bernonville OK, I did `bowtie2-build GCF_000001405.39_GRCh38.p13_genomic.fna human_genome.fa` first, and the output is 4 files `human_genome.fa.1.bt2, human_genome.fa.2.bt2, human_genome.fa.3.bt2, human_genome.fa.4.bt2`. I wonder how can I map my reads using these **4** files while in the command you have proposed you have only genome.fa Anyways I have already tried this `bowtie2 -x human_genome.fa /mnt/e/nastya/SLX066-02/B-dVMV-RIG-1/B-dVMV-RIG-1_ACAGTG_L008_R1_001.rc.fastq.gz -S bowtie2_mapping.sam --un-conc-gz unmappedreads_%.fq.gz` and it gave me an (ERR) bowtie2 exited
Dmitrii Trubetskoy
CC BY-SA 4.0
2020-06-20T12:17:34.827Z
try with `bowtie2 -x human_genome.fa -U /mnt/e/nastya/SLX066-02/B-dVMV-RIG-1/B-dVMV-RIG-1_ACAGTG_L008_R1_001.rc.fastq.gz -S bowtie2_mapping.sam --un-conc-gz unmappedreads_%.fq.gz`, adding the -U argument to indicate the single end read file
thomas duge de bernonville
CC BY-SA 4.0
2020-06-20T18:29:46.683Z
`18035544 reads; of these:
  18035544 (100.00%) were unpaired; of these:
    17979950 (99.69%) aligned 0 times
    55594 (0.31%) aligned exactly 1 time
    0 (0.00%) aligned >1 times
0.31% overall alignment rate` this is what I have now. I also have unmappedreads_1.fq and unmappedreads_2.fq archives of 1Kb each.. Is that ok? Why are there 2 files? Which one can I use with DI-tector? I also got a `bowtie2_mapping.sam` file which is about 3Gb. Thanks!
Dmitrii Trubetskoy
CC BY-SA 4.0
2020-06-21T11:45:02.903Z
This is because you asked for it `--un-conc-gz unmappedreads_%.fq.gz`. But I don't know if your input data is indeed paired. What is the size of each unmapped*fq?  What input is required for DI-tector?
thomas duge de bernonville
CC BY-SA 4.0
2020-06-21T20:14:13.357Z
Thanks again, and sorry for taking your time. There are two unmapped reads and the size is 1Kb each. In fact, they are empty. DI-tector requires unmatched reads in .fastq format + Virus_ref genome (I believe it's an indexed .fasta viral genome)
Dmitrii Trubetskoy
CC BY-SA 4.0
2020-06-22T06:50:15.260Z
ok I think it is my fault, you do not have paired-end reads but single end. So the correct option should be `--un-gz`, instead of `--un-conc-gz`, as described in the [bowtie2](http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml) manual.
thomas duge de bernonville
CC BY-SA 4.0
2020-06-19T04:24:31.920Z
2020-06-19T04:24:31.920Z
thomas duge de bernonville
<p>I am trying to follow the Beauclair et al. paper (free version here <a href="https://rnajournal.cshlp.org/content/24/10/1285.long" rel="nofollow noreferrer">https://rnajournal.cshlp.org/content/24/10/1285.long</a>) for identifying defective genomes using their DI-tector program.</p>
<p>According to figure 2 <a href="https://i.sstatic.net/8e1SZ.jpg" rel="nofollow noreferrer"><img src="https://i.sstatic.net/8e1SZ.jpg" alt="enter image description here" /></a></p>
<p>Before using the program itself, it is highly recommended to filter my raw read using <strong>bwa</strong> and then <strong>samtools</strong></p>
<p>Unfortunately I have no idea where to start - I have browsed official bwa and samtools sites, I have installed them using the following commands:</p>
<pre><code>sudo apt-get update -y
sudo apt-get install -y bwa
</code></pre>
<p>I have indexed my reference genome with the following command (as I read it is required before the analysis)</p>
<pre><code>bwa index sequence.fasta
</code></pre>
<p>Then, I have installed samtools</p>
<pre><code>cd ..
wget https://github.com/samtools/samtools/releases/download/1.9/samtools-1.9.tar.bz2
tar -vxjf samtools-1.9.tar.bz2
cd samtools-1.9
make
sudo apt install samtools
</code></pre>
<p>I wonder what will be the further steps to get these &quot;UNMATCHED READS at the end of step i&quot;, since Di-tector program takes virus_reference genome (I have it) and these unmatched reads.</p>
<p>First of all I want to align my .fastq files vs <code>GCF_000001405.39_GRCh38.p13_genomic.fna</code> - what will be the output file? If I got it right - I need this output file to filter with samtools, am I correct?
If I am not mistaken, the output after samtools will be the .sam file that I will have to put into a DI-tector - correct?</p>

<samtools><bwa>
How can I get unmatched reads for defective genomes analysis using bwa and samtools?
# How can I get unmatched reads for defective genomes analysis using bwa and samtools?

(Asked by: username_0 on 2020-06-18)

I am trying to follow the Beauclair et al. paper (free version here https://rnajournal.cshlp.org/content/24/10/1285.long) for identifying defective genomes using their DI-tector program.
According to figure 2 
Before using the program itself, it is highly recommended to filter my raw read using bwa and then samtools
Unfortunately I have no idea where to start - I have browsed official bwa and samtools sites, I have installed them using the following commands:
sudo apt-get update -y
sudo apt-get install -y bwa

I have indexed my reference genome with the following command (as I read it is required before the analysis)
bwa index sequence.fasta

Then, I have installed samtools
cd ..
wget https://github.com/samtools/samtools/releases/download/1.9/samtools-1.9.tar.bz2
tar -vxjf samtools-1.9.tar.bz2
cd samtools-1.9
make
sudo apt install samtools

I wonder what will be the further steps to get these "UNMATCHED READS at the end of step i", since Di-tector program takes virus_reference genome (I have it) and these unmatched reads.
First of all I want to align my .fastq files vs GCF_000001405.39_GRCh38.p13_genomic.fna - what will be the output file? If I got it right - I need this output file to filter with samtools, am I correct?
If I am not mistaken, the output after samtools will be the .sam file that I will have to put into a DI-tector - correct?

**Tags:** samtools, bwa

---

## 1 answers

### Answer by username_1 on 2020-06-19. Score: 1

If you want to keep track of the unmapped reads, I recommend you to use Bowtie2 instead of bwa.
Start by building an index:
bowtie2-build genome.fa genome.fa

Then map your reads, using arguments to save unmapped reads in a specific file like this:
bowtie2 -x genome.fa -1 reads_1.fq -2 reads_2.fq -S bowtie2_mapping.sam --un-conc-gz unmappedreads_%.fq.gz

It will write the alignment as well as unmappedreads_1.fq.gz and unmappedreads_2.fq.gz that you may use in DI-tector.

> **Comments:**
> **username_2**: My thoughts exactly. Personally I don't know why a paper in 2018 would want to use BWA over Bowtie2. Any ideas?  
> **username_1**: after putting a glance at the article, it is no really clear how they switch from bwa mem to bwa samse to align reads. Maybe it has something to see with the flags in the output .sam file, required for DI-tector proper running. But it is only a supposition...  
> **username_0**: Thank you. But what about `samtools` step? I mean, my read has to be aligned with bwa (or bowtie2 - great I will try that) *FIRSTLY* with the human genome, then *SECONDLY* with virus genome and only after these two steps it has to be filtered through `SAMTOOLS` and only after that I have to have the unmapped reads that I can put into a DI-tector.  
> **username_1**: I think it depends on the format of input data required for DI-tector. Bowtie2 will return .fastq files. In the article, they used samtools with -Bs -f 4 options probably to return unaligned reads into a .sam with specific CIGAR labels.  
> **username_0**: @thomas duge de bernonville OK, I did `bowtie2-build GCF_000001405.39_GRCh38.p13_genomic.fna human_genome.fa` first, and the output is 4 files `human_genome.fa.1.bt2, human_genome.fa.2.bt2, human_genome.fa.3.bt2, human_genome.fa.4.bt2`. I wonder how can I map my reads using these **4** files while in the command you have proposed you have only genome.fa Anyways I have already tried this `bowtie2 -x human_genome.fa /mnt/e/nastya/SLX066-02/B-dVMV-RIG-1/B-dVMV-RIG-1_ACAGTG_L008_R1_001.rc.fastq.gz -S bowtie2_mapping.sam --un-conc-gz unmappedreads_%.fq.gz` and it gave me an (ERR) bowtie2 exited  
> **username_1**: try with `bowtie2 -x human_genome.fa -U /mnt/e/nastya/SLX066-02/B-dVMV-RIG-1/B-dVMV-RIG-1_ACAGTG_L008_R1_001.rc.fastq.gz -S bowtie2_mapping.sam --un-conc-gz unmappedreads_%.fq.gz`, adding the -U argument to indicate the single end read file  
> **username_0**: `18035544 reads; of these:   18035544 (100.00%) were unpaired; of these:     17979950 (99.69%) aligned 0 times     55594 (0.31%) aligned exactly 1 time     0 (0.00%) aligned >1 times 0.31% overall alignment rate` this is what I have now. I also have unmappedreads_1.fq and unmappedreads_2.fq archives of 1Kb each.. Is that ok? Why are there 2 files? Which one can I use with DI-tector? I also got a `bowtie2_mapping.sam` file which is about 3Gb. Thanks!  
> **username_1**: This is because you asked for it `--un-conc-gz unmappedreads_%.fq.gz`. But I don't know if your input data is indeed paired. What is the size of each unmapped*fq?  What input is required for DI-tector?  
> **username_0**: Thanks again, and sorry for taking your time. There are two unmapped reads and the size is 1Kb each. In fact, they are empty. DI-tector requires unmatched reads in .fastq format + Virus_ref genome (I believe it's an indexed .fasta viral genome)  
> **username_1**: ok I think it is my fault, you do not have paired-end reads but single end. So the correct option should be `--un-gz`, instead of `--un-conc-gz`, as described in the [bowtie2](http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml) manual.  

---

