An Introduction to NCBI BLAST

The Basic Local Alignment Search Tool (BLAST) is a program that can detect sequence similarity between a query sequence and sequences within a database. The ability to detect sequence similarity allows us to identify putative orthologous genes in a sequence from another species. It also allows us to determine if a gene or a protein is related to other known genes or proteins.

BLAST is popular because it can quickly identify regions of local similarity between two sequences, which we can then use to infer homologous sequences. More importantly, BLAST uses a robust statistical framework that can determine if the alignment between two sequences is statistically significant. In this exercise, we will use the National Center for Biotechnology Information (NCBI) BLAST service to help us annotate a sequence from the Drosophila yakuba genome (unknown.fna in the exercise package).

Before we begin the analysis, we should first familiarize ourselves with the NCBI BLAST web interface. Open a new web browser window and navigate to the NCBI BLAST home page. In this exercise, we will only use a few of the tools available on the NCBI BLAST website. To learn about the more advanced options available (such as setting up My NCBI accounts), click on the "Help" link on the main navigation bar to access the documentations for NCBI BLAST (Figure 1).

All the NCBI BLAST pages have the same header with four links:

Besides the main toolbar, there are two other sections of the NCBI BLAST web interface that are of interests: the "Web BLAST" section contains links to the common BLAST programs and the "Specialized searches" section contains links to additional tools for performing sequence searches (e.g., use CD-search to identify conserved domains within a query sequence). The type of BLAST search you need to use will depend primarily on the type of query sequence and the database you would like to search.

Four of the five common BLAST programs are available through the "Web BLAST" section of the NCBI BLAST home page (Figure 2, top). The program tblastx, which translates the nucleotide query and nucleotide database when it performs the sequence comparisons, is not listed under the "Web BLAST" section. However, you can access this program by clicking on any of the BLAST programs in the "Web BLAST" section and then click on the "tblastx" tab in the NCBI BLAST search form (Figure 2, bottom).

The basic BLAST programs are summarized below:

In this exercise, we will characterize an unknown genomic sequence (unknown.fna) and determine if it has sequence similarity to any known genes. One strategy we can use is to search for sequence similarity to mRNA sequences in the NCBI Reference Sequence (RefSeq) database.

When we set up a BLAST search, there are three basic decisions we must make: the BLAST program we want to use, the query sequence we want to annotate, and the database we want to search. In addition, we can change several optional parameters (such as the Expect threshold and low complexity filters) in order to modify the behavior of BLAST. In this case, we will set up our BLAST search using mostly default parameters. We will use the blastn program to search our sequence (query) against the NCBI Reference Sequence (RefSeq) RNA database (Figure 3).

Once the search is complete, a new web page will appear with the BLAST report. For teaching purposes, the BLAST output (blastnInitial.txt) is available in the package for this exercise.

The top left panel (Figure 5) of the BLAST results page shows the parameters used in the BLAST search (e.g., database name, query ID, query length). The controls in the top right panel (Figure 5) can be used to filter the BLAST hits by organism, percent identity, and Expect value (E-value).

The details of the BLAST results are organized into the four tabs below these two panels: "Descriptions", "Graphic Summary", "Alignments", and "Taxonomy". We will go through each of these sections in order to interpret our blastn output.

Now that we have a better understanding of how the BLAST report is organized, we are ready to interpret the blastn results. The "Descriptions" and the "Graphic Summary" tabs (Figure 6 and Figure 8) show that many of the top hits are much more significant (with E-values of 0.0) than the rest of the blastn hits. Most of these top hits contain regions of sequence similarity that span the entire length of the query sequence (Figure 8). Looking at the descriptions and the corresponding GenBank records, it appears that these blastn hits correspond to the gene legless (also known as BCL9) in different Drosophila species.

From the blastn hit list, click on the description that corresponds to the best hit to an experimentally verified gene to navigate to the alignment between this hit and the unknown sequence.

While it is useful to examine the mRNA sequence, often we are specifically interested in the protein products of the gene. Most RefSeq mRNA records include both translated and untranslated regions (i.e., the 5' and 3' UTRs). If we are interested in the protein, we would need to determine the location of the Coding DNA Sequences (CDS). Because every mRNA in the RefSeq RNA database has a corresponding sequence in the RefSeq Protein database, we will search our D. yakuba nucleotide sequence against the RefSeq Protein (refseq_protein) database with blastx. We now have all the information we need to setup the blastx search.

For teaching purposes, the blastx search result (blastxRefSeqProtein.txt) is available in the package for this exercise. (Note that this blastx search can take several minutes to complete.)

The blastx report is similar to the blastn report. It consists of the "Descriptions", "Graphic Summary", "Alignments", and "Taxonomy" tabs. The "Graphic Summary" tab shows the highly significant hits to the legless protein in D. melanogaster and similar proteins in other Drosophila species. It also shows a few significant hits to transposases in the region between 6000-8000 bp of our sequence (Figure 21).

These hits suggest our sequence contains a type of repetitious element called a transposable element. In future exercises, we will learn how we can reduce the number of spurious hits in our BLAST reports by masking these elements prior to performing the BLAST search. For now, we will ignore these additional matches and focus on the best manually curated RefSeq hit — the D. melanogaster legless protein (NP_651922.1; Figure 22).

We can analyze the blastx alignments in the same way that we have previously analyzed the blastn report. However, because blastx translates the input sequence in all 6 reading frames before comparing our sequence with the protein database, there is an additional "Frame" field in each alignment block. The frame begins with either + or -, which corresponds to the relative orientation of our sequence compared to the protein. The number following the relative orientation in the frame field ranges from 1 to 3, which reflects the reading frame that produces the translated peptide sequence. Collectively, the relative orientation and the number can be used to represent all 6 reading frames. A frame shift between two alignment blocks in the blastx match often indicates that the two alignment blocks correspond to different coding exons. The "Positives" field corresponds to the number of amino acids that are either identical or have similar chemical properties between the translated query and the subject sequences (Figure 23).

Similar to the blastn alignment, each alignment block in our blastx report also consists of three lines: the query sequence, the matching sequence, and the subject sequence. Note that the query sequence has been translated into the corresponding amino acid sequence in the reading frame specified by the "Frame" field. However, the coordinates of the query sequence are still relative to the original nucleotide sequence. Like our blastn alignment, the grey lowercase residues in the query sequence correspond to low complexity sequences that were masked by BLAST.

There are some minor differences in the matching sequence of the blastn and blastx outputs. Residues in the matching sequence represent amino acids that are identical between the query and subject sequences. The "+" character denotes amino acids that are different between the query and subject, but these different amino acids have similar chemical properties. A space indicates that the two aligned amino acid in the query and subject are different, and they have different chemical properties.

The blastn and blastx search results suggest that we have identified the putative ortholog of the legless gene in our D. yakuba sequence. However, in order to construct a complete gene model, we must resolve the discrepancies in the alignments of our blastn and blastx output. Because coding sequences are under strong selective pressure and are more likely to be conserved than other regions of the genome, our first step is to identify the coding sequences of our putative gene. To begin the more detailed analysis, we will perform a series of BLAST searches using the amino acid sequence of each exon in the D. melanogaster version of the legless gene. It will be helpful to our annotation efforts if we can obtain the amino acid sequence that corresponds to each exon individually. Fortunately, we can easily obtain the individual exon sequences using the Gene Record Finder.

In the "mRNA Details" section of the gene report, we notice that there is only one isoform of the legless gene in D. melanogaster (lgs-RA, the A isoform of lgs). We can access all the transcribed exons through the "Transcript Details" tab and all the coding sequences (CDS) through the "Polypeptide Details" tab (Figure 25).

To retrieve the amino acid sequence for each coding sequence, click on the row that corresponds to the coding sequence in the CDS table (Figure 26).

To determine the locations of the coding sequences, we will perform BLAST searches to compare the individual exons from D. melanogaster with our D. yakuba sequence. Because we are comparing a protein sequence against a nucleotide sequence, we will use the tblastn program for our search. In order to prevent BLAST from masking low complexity regions in our protein, we will turn off the low complexity filter. In addition, because we are only comparing two sequences, we will also turn off compositional adjustments under scoring parameters.

For teaching purposes, the BLAST output (bl2lgsExon1_tblastn.txt) is available in the package for this exercise.

The "Alignments" tab of the tblastn output shows that the first coding exon has a length of 60 amino acids, and it aligns to 9344-9168 of the query sequence when it is translated in frame -2 (Figure 28).

We can use the same strategy to map the rest of the coding exons.

We could apply the same strategy to map the locations of the putative untranslated regions (UTRs) using bl2seq with the blastn program. Go back to the Gene Record Finder web browser window. The genome browser image in the "mRNA Details" panel shows that the first CDS of legless (1_9474_0) is part of a larger exon (lgs:1) that includes the 5' UTR (Figure 29). To estimate the location of the 5' UTR in our unknown sequence, we will perform a blastn search of the exon lgs:1 against the unknown sequence, and then compare the alignment with the tblastn search result for CDS 1_9474_0.

To retrieve the exon sequence for the first transcribed exon, select the "Transcript Details" tab in the Gene Record Finder and then click on the first row in the exon table (Figure 30).

We can then use the following steps to perform the blastn search:

For teaching purposes, the BLAST output (bl2lgsExon1_blastn.txt) is available in the package for this exercise.

In this exercise, we have used multiple BLAST programs to identify and characterize a putative gene in a genomic sequence from D. yakuba. You are now ready to tackle some of the more challenging BLAST exercises on the GEP website: