In this tutorial, we're going to view the aligned reads and variants that we called in the past two lessons in the Integrated Genomics Viewer from the Broad Center. You'll need the output from Introduction to mapping (bowtie, BWA) and Introduction to variant calling (Samtools).
Getting situated
We're going to assume that you have an introduction_to_mapping directory, and that you have existing results in those subdirectories. If not, see us and we can tell you where to copy "canned" results from.
cd $SCRATCH/intro_to_mapping
Prepare a GFF feature file for the reference sequence
IGV likes its reference genome files in GFF (Gene Feature Format). Unfortunately, our old friend bp_seqconvert.pl doesn't do GFF. Fortunately, it's cousin bp_genbank2gff3.pl does.
module load bioperl cp /corral-repl/utexas/BioITeam/ngs_course/scripts/bp_genbank2gff3.pl . ./bp_genbank2gff3.pl NC_012967.1.gbk
Where's the output? Take a look at the *.gbk and *.gff files and try to get a handle on what is going on in this conversion.
(NB: I had to "fix" bp_genbank2gff3.pl to work on our GenBank file, by removing a line in the Perl script that caused an error. If someone knows of an easier way to get a GFF file from a file downloaded from Genbank, please share!)
Copy files to your desktop
IGV is an interactive graphical viewer program. You can't run it on TACC, so we need to get the relevant files back to your desktop machine.
- Indexed reference FASTA files
- GFF reference sequence feature files
- Sorted and indexed mapped read BAM files
- VCF result files
We're going to copy all of these into a new directory called IGV to make it easier to just transfer the ones that we want.
mkdir IGV cp NC_012967.1.gbk.gff IGV cp samtools_bowtie/NC_012967.1.fasta IGV cp samtools_bowtie/NC_012967.1.fasta.fai IGV cp samtools_bowtie/SRR030257.sorted.bam IGV/bowtie.sorted.bam cp samtools_bowtie/SRR030257.sorted.bam.bai IGV/bowtie.sorted.bam.bai cp samtools_bwa/SRR030257.sorted.bam IGV/bwa.sorted.bam cp samtools_bwa/SRR030257.sorted.bam.bai IGV/bwa.sorted.bam.bai cp comparison/* IGV
Now, copy this entire IGV directory back to your local Desktop machine.
Launching IGV
Navigate a web browser to this page: http://www.broadinstitute.org/software/igv/download
Go ahead and click on the "Launch with 2 GB" option. This will download a "Java Web Start" file that you can launch by locating on your Desktop and double-clicking.
Load genome into IGV
From the main window of IGV, click on File -> Import Genome and you should be presented with the following window.
Enter the ID and Name of the Genome you are working with and select the path to your *.fasta file (the index, *.fai file needs to be in the same directory), then select the path to your *.gff file for the Gene File.
Load mapped reads into IGV
From the main window of IGV, click on File -> Load File. Choose bowtie.sorted.bam
After importing your reference genome and loading an alignment file, your screen should look similar to the following:
And you are now free to investigate different areas and their alignments in the genome.
Load variant calls into IGV
We're really interested in places in the genome where we think there are mutations. You can load the VCF file to check out those spots, but first you need to (guess what?) index it.
You can do this from within IGV:
- Choose File -> Run igvtools....
- Choose "index" from the commands drop-down menu.
- Select your *.vcf file (Ex:
bowtie.vcf) for "Input File" - Click the "run" button.
It will look like nothing has happened, but you can now close the "Run" window and choose File -> Load File. If you navigate to your IGV directory you will now have a brand new bowtie.vcf.idx file. You can now load the file bowtie.vcf, and it will show up as a new track at the top of your window.
Tip: You can index BAM and FASTA files the same way if you haven't already created indexes for them.
Navigating in IGV
There are a lot of things you can do in IGV. Here are a few:
- Jump to the next point of interest. Click on a track name on the left side of the window (Ex: bowtie.vcf), then it will be selected. You can then use
control-fandcontrol-bto jump forward and backward within that list of features. Try this on the variant calls. - Jump right to a gene. Typing its name in to the search box. Try "topA".
- Navigate more quickly. The
homeandendkeys jump a page backward and forward. - Load multiple VCF and BAM alignments at once. Try this to compare a few different regions between the bowtie and BWA results.
- Change the appearance of genes. Right click on the gene track and try "expanded". Experiment with the other options.
- Change the appearance of reads. Right click on a BAM track and choose "show all bases" and "expanded". Experiment with the other options.
Exercises
- Check out the rbsA gene region? What's going on here?
- What is going on in the pykF gene region? You might see red read pairs. What does that mean?
- Go to coordinate 2,698,092. Compare the bowtie and BWA alignments. Can you explain what's going on here?
- Go to coordinate 475,263. What's the story here?

