Gene Structure Annotation and Analysis Using PASA

PASA was originally developed at The Institute for Genomic Research in 2002 as an effort to automatically improve gene structures in Arabidopsis thaliana. Since then, it has been applied to numerous Eukaryotic genome annotation projects including Rice, Aspergillus species, Plasmodium falciparum, Schistosoma mansoni, and Aedes aegypti. We have also successfully applied it to lesser extents in both Mouse and Human among others.

Functions of PASA include:

PASA is composed of a pipeline of utilities that perform the following ordered set of tasks:

PASA runs on a UNIX/LINUX-based architecture. PASA involves components written in Perl and C++. Utilities used by PASA, including GMAP, are wrapped by Perl code. Results are stored primarily within a MySQL database, and are available for analysis using the companion suite of Web-based tools and command-line utilities. Running PASA to generate alignment assemblies requires only two inputs: a multi-fasta file for the genome and a multi-fasta file for the transcripts; if FL-cDNAs are available, then an additional file containing listing the accessions of those FL-cDNAs can be provided as a third input. In order to compare the assemblies to existing gene structure annotations, preexisting gene structure annotations can be provided in GFF3 format, or imported by a user-customized data adapter (described below).

Sample data and a preconfigured complete PASA pipeline are available for demonstration purposes, all included in the software distribution.

Download the latest version of the PASA software straight from Sourceforge

or (preferred) pull the bleeding-edge code using CVS like so:

In addition to the PASA software obtained here, you will need the following:

Move the PASA distribution to a location on your filesystem that we can call PASAHOME, such as /usr/local/bin/PASA. From henceforth, we'll refer to this location as $PASAHOME.

The PASA distribution includes the following utilities that you should build and centrally install:

After installing each of the software tools above, all that is needed before running PASA is to configure it. The PASA configuration relies on the file: $PASAHOME/pasa_conf/conf.txt

Recursively copy (cp -r) the $PASAHOME area to the cgi-bin directory of your webserver. Change permissions on everything so that it is world executable (ie. % chmod -R 755 ./PASA ) Now, visit the URL for the status report page for the pasa database you created during the pasa run above.

http://yourServerName/cgi-bin/PASA/cgi-bin/status_report.cgi?db=$mysqldb

This will provide some summary statistics and links to additional web-based utilities for navigating the results from your pasa run.

Now that you have a URL for your base PASA url, update your original configuration file at: $PASAHOME/pasa_conf/conf.txt to set the value of BASE_PASA_URL=http://yourServerName/cgi-bin/PASA/cgi-bin/

For more info, visit the Tour of the PASA web portal

Have each of these files in the same working directory. Then, run the seqclean utility on you transcripts like so:

If you have a database of vector sequences (ie. UniVec), you can screen for vector as part of the cleaning process by running the following instead:

This will generate several output files including transcripts.fasta.cln and transcripts.fasta.clean Both of these can be used as inputs to PASA.

For my PASA primary configuration file, I'm using $PASAHOME/pasa_conf/sample_test.conf symlinked to conf.txt in that same directory. In it, the mysql server and user/password information our set. This configuration file need only be established once per PASA software installation. You should configure your conf.txt file immediately if you haven't done so already.

Sample inputs are provided in the $PASAHOME/sample_data directory. We'll use these inputs to demonstrate the breadth of the software application, including using sample DATA ADAPTERs to import existing gene annotations into the database, and tentative structural updates out.

The PASA pipeline requires separate configuration files for the alignment assembly and later annotation comparison steps, and these are configured separately for each run of the PASA pipeline, setting parameters to be used by the various tools and processes executed within the PASA pipeline. Configuration file templates are provided as $PASAHOME/pasa_conf/pasa.alignAssembly.Template.txt and $PASAHOME/pasa_conf/pasa.annotationCompare.Template.txt, and these will be further described when used below.

The next steps explain the current contents of the sample_data directory. You need not redo these operations:

The following steps, you must execute in order to demonstrate the software.

The output from this step includes several output files that are written:

A couple of additional files that are sometimes useful include the following: (here is where you go when you want to find out why certain alignments didn't make it into PASA assemblies or why they failed validation):

Incorporating PASA Assemblies into Existing Gene Predictions, Changing Exons, Adding UTRs and Alternatively Spliced Models

The PASA software can update any preexisting set of protein-coding gene annotations to incorporate the PASA alignment evidence, correcting exon boundaries, adding UTRs, and models for alternative splicing based on the PASA alignment assemblies generated above.

Comparing to and updating existing gene structure annotations requires that we import these annotations into the PASA database, and are able to extract the suggested updates. PASA utlizes annotation data adapters to achieve this. GFF3 data adapters are included in the PASA distribution, but you can write your own, and directly tie the PASA pipeline to your own informatics infrastructure (ie. other relational database). If you'd prefer to not use GFF3 and to write your own data adapters, visit the PASA data adapter cookbook.

A sample gff3-formatted annotation file is provided in our sample_data directory as orig_annotations_sample.gff3 and can be loaded like so:

Before loading your own GFF3-formatted annotation files, be sure to check them for PASA compatibility like so:

The above gff3-validator will report any entries in your gff3 file that it does not recognize, understand, or otherwise parse properly. It's not a general purpose gff3-validator since it cares only about your protein-coding genes. (note that you should only feed protein-coding genes to PASA using the loader above).

Now that the original annotations are loaded, we can perform a comparison of the PASA alignment assemblies to these preexisting gene annotations, to identify cases where updates can be automatically performed to gene structures in order to incorporate the transcript alignments.

I've copied the ../pasa_conf/pasa.annotationCompare.Template.txt file to our working directory as annotCompare.config. Then, I replaced the MYSQLDB=<MYSQLDB> line with MYSQLDB=pasa_sample_db as before with the alignAssembly.config file. Notice this config file contains numerous parameters that can be modified to tune the process to any genome of interest. We'll leave these values untouched for now, relying on the defaults used by PASA, and we'll revisit parameterization later. For most purposes, the defaults are well suited. Run the annotation comparison like so:

Once the annotation comparison is complete, PASA will output a new GFF3 file that contains the PASA-updated version of the genome annotation, including those gene models successfully updated by PASA, and those that remained untouched. This file will be named output.gene_structures_post_PASA_updates.$pid.gff3, where $pid is the process ID for this annotation comparison computation.

You should revisit the status_report.cgi web page as described above under Setting Up the PASA Web Portal. There, you will be able to navigate the results of the comparison and examine the classifications for annotation updates assigned to each pasa alignment assembly.

The results from running PASA on our sample data set can be examined via the PASA web portal. For example purposes, I've saved a few of the reports generated by the PASA web displays (note pages are generated on-the-fly, however these are provided as static only for example purposes).

The above are just a few examples. Install the PASA portal and navigate your PASA results.

If seqclean was used to clean the transcript sequences, and both the cleaned and original transcript databases were provided in the alignment assembly run of the PASA pipeline as described, then the polyadenylation sites as evidenced in the original transcript sequences and identified as part of the seqclean process were mapped to the genome. The termini of the polyadenylated transcripts are compared to the genome, and those transcripts that truly appear to be polyadenylated and not resulting from an artifact of internal priming to an A-rich region, are reported as candidate polyA sites. The genome coordinate reported as the polyA site is the nucleotide to which polyA is added, so it corresponds to the last non-polyA nucleotide of the polyadenylated transcript. An example of a candidate polyA site can be extracted from one of the output files (default output.polyAsite_analysis.out) like so:

An additional fasta file (default output.polyAsites.fasta) summarizes all mapped polyA sites supported by the transcripts. A 100 bp segment of the genome sequence is extracted and oriented, and the last nucleotide in uppercase corresponds to the residue to which polyA is added in the processed transcript. The site corresponding to our example above is as follows:

The accession is bundled like so:

The rest of the header indicates the number of transcripts supporting this polyA site followed by the list of those transcript accessions. The examples above were extracted from our sample data set provided. A more compelling example for Arabidopsis, using spliced transcripts only, is as follows:

PASA is a tool well suited to the identification and classification of alternative splicing isoforms as evidenced by incompatible transcript alignments. Overlapping alignments found incompatible in that they have some structural difference within their overlapping region, and due to their nature of incompatibility, they are relegated to different but overlapping alignment assemblies. PASA performs and all-vs-all comparison among the clustered overlapping alignment assemblies to identify the following categories of splicing variations:

The automated alternative splicing analysis can be run like so as exemplified from the sample_data directory:

The results are available in the default output files: output.alt_splice_label_combinations.dat:: a tab-delimited listing that contains all unique splicing labels for each pasa alignment assembly labeled with a variation. For example:

The PASA web portal provides numerous reports, graphs, and illustrations to navigate the results of the automated alternative splicing analysis.

In our current working directory, there's a file clusters_of_valid_alignments.txt that contains all the clusters of valid alignments in a simple text format like so:

The transcribed orientation is +,-, or ?. The ? orientation should be used only for single-exon transcript alignments for which the orientation of transcription is ambiguous. By default, PASA assigns all single-exon transcripts that lack evidence of polyadenylation to the ambiguous transcribed orientation. Given this input file, we can demonstrate the pasa alignment assembler like so:

Each cluster of transcript alignments is assembled separately and the results are outputted to stdout with illustrations.

Example input

Corresponding Output

The PASA alignment assemblies can be used to automatically extract gene structures for protein coding genes to be used for training gene finding software. To do this, the longest ORF is found within each PASA alignment assembly, allowing for 5' and 3' partials. Those long ORFs encoding at least 300 amino acids (configurable) are reported. After running the PASA to assemble all transcript alignments as described above, you can run the following:

for example, from our sample_data directory, we would run

which generates two files: trainingSetCandidates.fasta and trainingSetCandidates.gff3 The fasta file contains the protein sequences corresponding to the longest ORF found in each relevant PASA alignment assembly. The gff3 file provides the corresponding gene structure in the genome. The header of each entry in the fasta file provides the following, by example:

The accession corresponds to the PASA assembly. The contig and the coordinate span of the alignment assembly is provided, followed by the number of exons, and the length of the longest ORF. The type indicator can be any of the following: complete, 5prime_partial, 3prime_partial, or internal. The 5prime_partial are missing a start codon and translate to the very 5' end, 3prime_partial are missing a stop codon and translate to their very 3' end, and internal translate from the first to the last basepair in the sequence, missing a start and a stop codon. The complete category are of greatest interest for the prospects of genefinder training. Typically, we would search these complete proteins against the nonredundant protein database at GenBank and identify those ORFs that have good database matches across most of their length. Such entries can be confidently used for training, in addition to those particularly long ORFs that do not match known proteins and are sufficiently complex in sequence.

The PASA software and its original application are described in:

The use of PASA to analyze polyadenylation signals is described in:

Enhancements to PASA that automate the identifiation and classification of alternative splicing variations are described here:

The PASA pipeline was an important component of the eukaryotic genome annotation at The Institute for Genomic Research. Since that time, PASA has been successfully used by genome sequencing and annotation groups across the globe. I (Brian Haas) developed and actively maintain the software. Please contact me directly at bhaas@broad.mit.edu if you should have any questions or require assistance.