I work with fastq files containing NGS reads for some human DNA regions. The reference genome is hg19. I had two fastq files (pair-ended). I generated alignment BAM files. I used "bwa" and samtools to find a possible target gene region (chr7:55,242,376-55,242,574,). This corresponds to a region of the EGFR gene.

Here is a screenshot of the gene region


And here is a screenshot of the primer region primer

I also have a list of primers and reverse primers:



And I have some adaptors



I have to answer to the following questions:

1) What is the constitution of each read (adaptor+primer+amplified region)? The forward primer is obvious and corresponds to the first primer of the list


So why do we have the three other forward primers? I don't understand. And what about the correct reverse primer?

It seems to be the complementary sequence of the end region


Finally which adaptor is used?

Is [ADAPTOR1 - PRIMER1 - AMPLIFIED REGION] the correct answer? Are there other possibilities?

2) % of reads that map the human genome? what about the unmapped?

Can you please give me some clues about how to answer this question.

I just need some clues I want to do it by myself.

Thank you very much for your help.


You have asked quite a few questions here. I will try to answer a subset of them.

To answer your question, I started by going to the UCSC genome browser, and selected the BLAT from the tools dropdown menu. Then passing your primer sequences as queries, we can clearly see where they map to the reference genome.

enter image description here

Looking at the figure above, I can guess that your data comes from a targeted sequencing experiment (exome, or gene panel, etc); where you have primers around each exon of EGFR on two strands, on 5' end of exon sequence on each strand (think about DNA polymerase acting in PCR to amplify the number of DNA molecules before sequencing).

Regarding the adapter, I am not sure which read you are referring to. Can you clarify please?

About percentage of mapped vs unmapped reads, you can use bamtools.

Example Usage:

/user/me/src/bamtools/bin/bamtools-2.3.0 stats -insert -in my-sequence-file.bam

Example output:

Total reads: 103277668

Mapped reads: 90088436 (87.2293%)

Forward strand: 58136735 (56.2917%)

Reverse strand: 45140933 (43.7083%)

Failed QC: 6529806 (6.32257%)

Duplicates: 0 (0%)

Paired-end reads: 103277668 (100%)

'Proper-pairs': 87439672 (84.6646%)

Both pairs mapped: 87910438 (85.1205%)

Read 1: 51638834

Read 2: 51638834

Singletons: 2177998 (2.10888%)

Average insert size (absolute value): 6317.39

Median insert size (absolute value): 301

| improve this answer | |
  • $\begingroup$ Thank you very much. They all use the same 2 adaptors. Is that possible? It's seems weird for me. $\endgroup$ – Newbe Mar 7 '16 at 13:06
  • $\begingroup$ just one last question How do you pass your primer sequences as queries in BLAT? And obtain such a nice figure $\endgroup$ – Newbe Mar 7 '16 at 14:01
  • $\begingroup$ Glad I was able to help. I generated a fasta format input with your eight sequences. Then after performing BLAT, selected the browser link for one one your sequences, zoomed out by a factor of 10 a few times so that all the query sequences are in the range. By default, one gets a crowded display in genome browser. You can press 'hide all' button below the genome view figure, and follow by setting BLAT and refGene/refSeq options to full. $\endgroup$ – Noushin Mar 7 '16 at 17:40

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.