Overview

The Integrative Genomics Viewer (IGV) from the Broad Center allows you to view several types of data files involved in any NGS analysis that employs a reference genome, including how reads from a dataset are mapped, gene annotations, and predicted genetic variants.

Learning Objectives

  1. Create a custom genome database (usually used for microbial genomes).
  2. Load a pre-existing genome assembly (usually used for the genomes of model organisms and higher Eukaryotes).
  3. Load output from mapping reads to a reference genome.
  4. Load output from calling genetic variants.
  5. Navigate the view of the genome and interpret the display of this data.

Theory

Because NGS datasets are very large, it is often impossible or inefficient to read them entirely into a computer's memory when searching for a specific piece of data. In order to more quickly retrieve the data we are interested in analyzing or viewing, most programs have a way of treating these data files as databases. Database indexes enable one to rapidly pull specific subsets of the data from them.

The Integrative Genomics Viewer is a program for reading several types of indexed database information, including mapped reads and variant calls, and displaying them on a reference genome. It is invaluable as a tool for viewing and interpreting the "raw data" of many NGS data analysis pipelines.

With all that being said, the goal of visualizing the data with IGV is not to look at every read, or even every base in the reference genome (which is actually the smaller of the 2 possibilities!). IGV works best to dig deeper on something that you already are interested in either because of the gene itself, or is something that seems confusing (ie trying to figure out why mapping quality declines for a particular variant), or are trying to familiarize yourself with the concepts of what is going on. 

Installing IGV

This is done on your local computer not on TACC. IGV can not be installed on TACC which should make some sense to you as IGV is a program designed to let you visualize information and we know TACC doesn't allow GUIs.

There are multiple ways to launch IGV on a local computer. For this course I recommend in a separate a web browser window/tab going to: http://www.broadinstitute.org/software/igv/download and selecting the appropriate operating system with java included. Mac users will need to unzip the file and launch the application. Window's users will need to download, choose an installation location, agree to some licenses, navigate to your installation location and then launch the program. Believe it or not, this is a significantly improved process compared to actions that used to be required.

Viewing E. coli data in IGV

Data files

You can start this tutorial two ways:

  1. If you have completed the Mapping tutorial and the SNV calling tutorial, then you should use those files for part 1 of this tutorial. You can proceed with either one alone, or with both.

    I strongly encourage everyone to complete these tutorials rather than relying on the canned data. I routinely work with raw fastq files, map them, and call variants on them. Outside of this course I use IGV maybe 1 time a year, and I don't think I'm an uncommon case.
  2. This block is collapsed because I want to strongly strongly discourage anyone from jumping straight to this option. The mapping and variant calling tutorials are way way more important than this visualization tutorial. 
    $BI/gva_course/mapping/IGV  # location of example files
cp -r /corral-repl/utexas/BioITeam/gva_course/mapping/IGV .  # example command to copy to current directory
scp -r <username>@stampede2.tacc.utexas.edu:/corral-repl/utexas/BioITeam/gva_course/mapping/IGV . # to copy to a local computer skipping the step of copying to a stampede2 directory and secure copying from there.

    If you are sure you this is how you want to proceed and another offer of assistance with any problems you may be having with the mapping or SNV tutorial, won't tempt you, then skip down to #Launching IGV.

Prepare a GFF feature file for the reference sequence

IGV likes its reference genome files in GFF (Gene Feature Format) rather than the fasta or gbk formats we've been working with. While you may assume this is a job for our old friend bp_seqconvert.pl, that script actually doesn't deal with GFF files. So, we're going to use another tool for sequence format conversion called Readseq. Install the readseq package from bioconda to a new conda environment to get started.

You have done this several times now, you should likely be able to do this without expanding this code block. Try to figure the command out yourself and check your work rather than relying on a copy paste. 
conda create -n GVA-readseq -c bioconda readseq
readseq --version

Unlike several of the previous programs we have installed, the --version flag actually prints the entire help file rather than just versioning information. If you scroll up, you will see the first line of the output is actually the version. In this case, 2.1.30

readseq is a java based program which means it is envoked in a very different manner than anything we have worked with thus far. Luckily for us the conda package actually includes a wrapper allowing us to envoke the command simply by typing the readseq name like all the other programs we have worked with.

Review previous year's tutorial to read about how to envoke the program using java without the readseq wrapper

This is one of only 2 java based programs that this course covers. As the readseq wrapper conda provided makes this so much easier to envoke, we will use it. It is recommended to look back at a previous years tutorial to see how this was handled without the wrapper incase you encounter a java based program in your own work that doesn't have such a helpful wrapper, and need to know where to start.


It's a bit hard to figure out how to build the command yourself as, unlike most conventions, the program requires the unnamed arguments before the optional flag arguments, there is no example command in the help. To do the conversion that we want, and get things where they need to be for the rest of the tutorial use the following:

cds
mkdir GVA_IGV
cd GVA_IGV
readseq $SCRATCH/GVA_bowtie2_mapping/NC_012967.1.gbk -f GFF -o NC_012967.1.gbk.gff

A final oddity of the readseq program is that rather than displaying any kind of status message, or being silent when executed, the program actually displays the version of the readseq program itself. This is something that initially made me assume the conversion had failed. Take a look at the contents of the original Genbank file and the new GFF file and try to get a handle on what is going on in this conversion using commands like head and tail.

File naming conventions

You may notice that the output file appended a ".gff" ending to the ".gbk" ending rather than replacing it. This is can be done to demonstrate the order of operations performed on the file (in this case taking a gbk file and converting it to a gff file) a longer list of operations, such as sequential filtering a vcf file for frequency above 90%, with mapq scores above 20, on chromosome 7,  between 10,000,000 and 190,000,000 bp might result in file contents looking like the following:

sample.vcf
sample.vcf.filtered.freq90
sample.vcf.filtered.freq90.mapq20
sample.vcf.filtered.freq90.mapq20.chr7
sample.vcf.filtered.freq90.mapq20.chr7.10MB-190MB

some programs do not like unknown or effectively nonsense file endings, in which case you may need to append ".vcf" to the names above .


Copy files to your local computer

Again, since IGV is an interactive graphical viewer program that we'll be running on our local computer, we need to get the files we want to visualize onto your desktop machine.

The files we need include:

  • Indexed reference FASTA files
  • GFF reference sequence feature files
  • Sorted and indexed mapped read BAM files
  • VCF result files
  • In practice, depending on your analysis you may need or can use additional files but those won't be discussed here


Rather than transferring each file individually, from multiple different directories, the easiest (and most common, and best practice) thing to do is:

  1. create a new directory (suggest you include a keyword like export)
  2. copy all the files you want to transfer into the new directory
  3. compress the directory to speed up the transfer

In the case of this tutorial, since many of the tutorial output files had the same names (but resided in different directories) we need to be sure to give them unique destination names when you copy them into the new directory together. Additionally, to ensure you don't overwrite files that you want (here is another reminder about there being no undo command in linux) you can (and is good practice to) use the -n or -i option with the cp command. On stampede2 the -n command will not allow you to overwrite files, while the -i command will prompt you before overwriting anything. It is worth noting that there is a slight difference in different versions of linux/bash that with these command options and they may not work on all systems so double check with the help flag before using it. 

This is probably the largest code box in the entire course, note that some of these lines likely extend beyond the right side of the window and you may need to scroll to the right to see the entire command
cds
mkdir GVA_IGV_export
cp -i $SCRATCH/GVA_IGV_Tutorial/NC_012967.1.gbk.gff GVA_IGV_export  # this is the new file you just created above
cp -i $SCRATCH/GVA_bowtie2_mapping/NC_012967.1.fasta GVA_IGV_export
cp -i $SCRATCH/GVA_samtools_tutorial/NC_012967.1.fasta.fai GVA_IGV_export
cp -i $SCRATCH/GVA_samtools_tutorial/SRR030257.vcf GVA_IGV_export
cp -i $SCRATCH/GVA_samtools_tutorial/SRR030257.sorted.bam GVA_IGV_export/bowtie2.sorted.bam
cp -i $SCRATCH/GVA_samtools_tutorial/SRR030257.sorted.bam.bai GVA_IGV_export/bowtie2.sorted.bam.bai
tar -czvf GVA_IGV_export.tar.gz GVA_IGV_export

Now, copy the entire compressed IGV directory back to your local Desktop machine. Remember there is a separate tutorial that covers scp file transfers in more detail that can be found here. In this case, you would replace README in that tutorial with GVA_IGV_export.tar.gz and need to determine the full path to that file using the pwd command.

In the Remote terminal (on the Right)

pwd

Expect this to be just the $SCRATCH directory, but remember your local computer won't know what $SCRATCH means.

In the Left terminal (on the Left). Remember to replace the text in the <> brackets with your own information:

scp <username>@stampede2.tacc.utexas.edu:<full_path_to_file>/GVA_IGV_export.tar.gz .
# enter your password/tacc token

Once you have transferred the compressed file, you will need to expand the file on your local machine so that you can access the individual files in IGV. 

You can use the tar --help option to try to figure out how to extract the contents, but I suggest the the following command:
tar -xzvf GVA_IGV_export.tar.gz


Screens may differ slightly than images show

Most of the differences will be trivial differences between windows and mac versions, or difficulties in getting the identical level of zoom in/out on the genome, but there may be some that are important and impact learning. If you are unclear if the difference is important get my attention on zoom and we'll work through it.


For the remainder of the tutorial, work on your local machine. NOT TACC!

Load genome into IGV

From the main window of IGV:

  • click on Genomes >Load Genome from File...
    •  select the fasta version of your genome (NC_012967.1.fasta)
  • click on File> Load from File...
    • select the newly made gff annotated genome (NC_012967.1.gbk.gff)

Load mapped reads into IGV

From the main window of IGV, click on File > Load from File.... Choose bowtie2.sorted.bam

After importing your reference genome and loading an alignment file, click on the + button in the upper right until reads appear! Then, your screen should look similar to the following:

Load variant calls into IGV

We're really interested in places in the genome where we think there are mutations. In the Variant calling tutorial we identified such locations but lacked a good way to visualize them. This is your opportunity to visualize them. We have already transferred the SRR030257.vcf file back to your local computer, but before we can visualize them, we need to (guess what?) index it. To answer some common questions:

To answer some common questions:

  1. Yes this index is different than the index you made for the reference genome.
  2. Yes it is different than the index you made of the sam files.
  3. Yes it is different than the index you made of the sorted bam files.
  4. Yes this is confusing.
  5. No there was nothing that required 4 different programs/groups to call this by the same name.
  6. Yes there are more than 4 different types of indexes you may encounter.

Getting back to our tutorial, you can create your new index from within IGV:

  1. Choose Tools > Run igvtools....
  2. Choose "index" from the commands drop-down menu.
  3. Select the SRR030257.vcf file  for "Input File"
  4. Click the "run" button in the middle of the window.

It will look like nothing has happened aside from the appearance of "Done" in the messages box, but you can now close the igvvtools window and choose File > Load from File. If you navigate to your IGV directory, you will now see a brand new SRR030257.vcf.idx file. You can now load the SRR030257.vcf file, and it will show up as a new track near the top of your window. Note that it is the VCF file, not the newly made index file that you are selecting. like with the other files you loaded above, IGV requires you to have an index in the same directory as the file you are loading, but it is the file, not its index that you load.

You can also index BAM and FASTA files the same way inside of IGV if you haven't already created indexes for them. But, it's usually easier and quicker to do this on the command line at TACC as indexing BAM files in particular can be a computationally hefty task.


You are now free to investigate different areas and their alignments in the genome.

Navigating in IGV

There are a lot of things you can do in IGV. Here are a few:

  • Zoom in using the slider in the upper right. Do this until you see mapped reads and finally individual bases appear.
  • Navigate by clicking and dragging in the window. This is how you move left and right along the genome.
  • Navigate more quickly. Use page-up page-down, home, end.
  • Jump to the next point of interest. Click on a track name on the left side of the window (Ex: SRR030257.vcf), to select it. You can then use control-f and control-b to jump forward and backward within that list of features. Try this on the variant calls track to jump to regions of interest.
  • Jump right to a gene. (If you have gene features loaded.) Type its name into the search box. Try "topA".
  • Load multiple BAM alignments or VCF files at once. You can use this feature to compare different mappings of bowtie2 and bowtie2 with improved quality reads, or different filtering options from the SNV tutorial.
  • 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.

See the IGV Manual for more tips and how to load other kinds of data.

Exercises

  • Right click on one of the reads to bring up a set of options. Look at the color alignment by section to see what it currently is and what you might prefer it as.

  • Interested in determining the probability that a read is not where it should be? What is a typical mapping quality (MQ) for a read?

    The estimated probability that a read is mapped incorrectly is 10^(-MQ/10). Where MQ is the mapping quality.

  • Can you find a variant where the sequenced sample differs from the reference? This would be like looking for a needle in a haystack if not for the use of variant callers and the control-f and control-b options 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.

    • Gene is pcnB, mutation is a snp

    • Gene is infB, mutation is a snp

    • Deletion of the rbsD gene

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

    • There is a 16 base deletion in the gltB gene reading frame.

    • The read pairs are discordantly mapped. There was an insertion of a new copy of a mobile genetic element (an IS150 element) that exists at other locations in the reference sequence.

    • See if you can find more interesting locations. There are ~190 mutations total in this sample MOST of which are false positives.

Viewing Human Genome Data in IGV

Now that you've familiarized yourself with IGV using a "simple" bacteria, let's look at something a "little" more complex: the human genome.

Advanced exercise: human data scavenger hunt 

Throughout this class we have used the wget command several times to download files directly to TACC which is running bash on a linux environment. As IGV only works on your local computer rather than TACC, downloading files directly to your laptop is more practical. Unfortunately wget is a linux only command and unix does not have it so instead we shall use the curl command to download some human data to IGV and look at it. The following commands should be run in a terminal window that is not logged into tacc.

get some data if working on a mac
cd ~/Desktop
mkdir GVA_human_IGV
cd GVA_human_IGV
curl -O https://rnabio.org/assets/module_2/HCC1143.normal.21.19M-20M.bam
curl -O https://rnabio.org/assets/module_2/HCC1143.normal.21.19M-20M.bam.bai
get some data if working with powershell
#Option 1 using tacc as the intermediary:
cds
mkdir GVA_human_IGV
cd GVA_human_IGV
wget https://rnabio.org/assets/module_2/HCC1143.normal.21.19M-20M.bam
wget https://rnabio.org/assets/module_2/HCC1143.normal.21.19M-20M.bam.bai
#use scp to transfer these 2 files back to your local computer
#Option 2 downloading directly from website.
enter the 2 addresses into any web browser:
https://rnabio.org/assets/module_2/HCC1143.normal.21.19M-20M.bam
https://rnabio.org/assets/module_2/HCC1143.normal.21.19M-20M.bam.bai

If you look at the file names, you may notice that this bam file and its index correspond to human chromosome 21 from position 19 million to 20 million or less than 1/4 of our bacterial genome. This limited data set is to hopefully avoid IGV crashing. 

Steps:

  1. Close IGV (if you have it open from the first tutorial with your mapping, SNV, and SV data) and reopen it. 
  2. Select "Human hg19" as the reference genome from the top left drop down (you may need to select "more" to have hg19 as an option)
  3. Load the bam files you downloaded: File > Load from File…  and select HCC1143.normal.21.19M-20M.bam
  4. Turn on dbSNP annotations File > Load from Server… >Annotation > Variation and Repeats > dbSNP 1.4.7
  5. Navigate to chr21:19,500,000-19,500,001 to be able to view reads. (normally could see reads anywhere, but we have specifically downloaded only reads that map in a 1Mb window centered on chr21:19.5M
  6. Right click on the track name on the left and select sort alignments by start location

  7. There are 2 mutations visible in the chr21:19,479,237-19,479,814 region answer the following questions:

    Yes, both forward and reverse reads (red and blue if colored by strand) contain the SNPs compared to the reference

    The one on the left does not correspond to a dbSNP entry and is therefore more likely to be related to disease state.


  8. There are 2 SNPs visible in the chr21:19,666,833-19,667,007 region. Answer the following questions:

    No, each read only has 1 mutation on it, these are 2 different alleles each with its own SNP relative to 'wt'. Both are reported in dbSNP

    Neither is likely to be related to disease (or at least not to rare disease) as both mutations have previously been identified as naturally occurring by dbSNP


  9. Homozygous deletion. In the track on the left, right click and select 'view as pairs' to see linkage between R1 and R2 to see individual reads mapping to both sides of the deletion

  10. This is an example of poor alignment to a repetitive AluY element. Notice how of the read pairs that map with numerous SNPs have 1 read that maps with lots of SNPs and the other read maps with none? This is caused by mapping reads to a limited area of the whole genome, if these reads had been allowed to map to the entire genome it is very likely that both read pairs would map without SNPs somewhere else in the genome.

    How can we identify this region as an AluY element

    There are several methods that could accomplish this. 1 as pointed out in class would be to pull the sequence from this region and blast that region, but that sounds like a lot of work. A second and much better way would be to turn on some additional tracks in IGV. In this case Alu elements are identified as 'SINE' elements in the repeat masker data base.

    To turn this on, File > Load from Server… >Annotation > Variation and Repeats > Repeat Masker > SINE

    Information about what Repeat masker is doing and where it comes from can be found here: http://genome.ucsc.edu/cgi-bin/hgTrackUi?db=hg19&g=rmsk

  11. What other interesting things can you find?

An additional tutorial from another group working with the same human data can be found here if interested.

Optional Tutorial Exercises ...

Visualize how other samples map without calling variants:

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 mapping directory, these commands will perform the necessary steps.

Submit to the TACC queue or run in an idev shell

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

Repeat the last three commands for each SAM output file that you want to visualize in IGV.

Visualize mapped reads and predicted variants from any later tutorials

Alternative genome browsers:

Return to GVA2022 Home page

  • No labels

Confluence Documentation | Web Privacy Policy | Web Accessibility