Content Comparison

...

You will need to index your reference FASTA and convert your SAM output files into sorted and indexed BAM files. The "why?" behind these steps is described more fully in the Variant calling tutorial. If you are in your intro_to_mapping directory, these commands will perform the necessary steps.

samtools faidx NC_012967.1.fasta
samtools view -b -S -o bowtie/SRR030257.bam bowtie/SRR030257.sam
samtools sort bowtie/SRR030257.bam bowtie/SRR030257.sorted
samtools index bowtie/SRR030257.sorted.bam

Copy files to your desktop

...

• Why are some reads different colors? Hint: Try changing the display options to show read pairs and editing some of the distance constraints.
• What is a typical mapping quality (MQ) for a read? Convert this to the probability that it is mismapped.

  Expand
  title Remember the formula for a Phred quality score?
  The estimated probability that a read is mapped incorrectly is 10^(-MQ/10).

• Can you find a variant where the sequenced sample differs from the reference? This is going to be like looking for a needle in a haystack. Fortunately, we are going to learn how to use variant callers tomorrow and then we'll be able to zoom right to areas where there are discrepancies between reads and the reference genome that might indicate there were mutations in the sequenced E. coli.

  Expand
  Some interesting locations to look at for the time being... Some interesting locations to look at for the time being...
  Coordinate 161,041. What gene is this in and what is the effect on the protein sequence? Coordinate 3,248,957. What gene is this in and what is the effect on the protein sequence? Coordinate 4,015,892. WHat What is different about the reads mapped to this location? Coordinate 3,894,997. What type of mutation is this? Coordinate 1,733,647. What's going on here? See if you can find more interesting locations. There are ~40 mutations total in this sample.

...

• Check out the rbsA gene region? What's going on here? Go
  What is going on in the pykF gene region? You might see red read pairs. What does that mean?

  Expand
  Answer Answer
  There was a large deletion. Can you figure out the exact coordinates of the endpoints?

• Navigate to coordinate 475,288. Compare the bowtie, BWA, and bowtie2 alignments. Can you explain what's going on here?

  Expand
  Answer Answer
  There is a 16 base deletion in the gltA gene reading frame.

  Go to coordinate 1,733,647. Compare the bowtie,
• BWA, and bowtie2 alignments. What's the story What is going on in the pykF gene region? You might see red read pairs. What does that mean? Can you guess what type of mutation occurred here?

  Expand
  Answer Answer
  The read pairs are discordantly mapped. There was an insertion of a new copy of a mobile genetic element (an IS150 element) . Note the reads with discordant pair mapping.

Workflow 2: Viewing Human Genome Data in IGV

...

Data from the CEU trio from the 1000 Genomes Project can be found directly from the Broad's server for IGV.

Find the dbSNP accession number for the SNP apparent in the two 1000 genomes project trios in the intron between exons 8 and 9 of the GABBR1 gene.

Steps:

1. Download and install the Integrative Genome Viewer from the Broad Institute.
2. Select "Human hg18" as the reference genome
3. Get some data: File -> “Load from Server…” -> 1000 genomes -> CEU and YRI trios
4. Navigate to the GABBR1 gene and the rightmost exons
5. Zoom in until you find the SNP
6. Load and look at the SNP track: File -> Load from server -> Annotations -> Variants -> dbSNP

This is whole genome coverage data; later we'll look at exome data.

...

rs29220

From here...

You can also use IGV to visualize RNA-seq data in later tutorials.

Check out alternative genome browsers:

• Tablet - a lightweight NGS data browser

  that exists at other locations in the reference sequence.