Skip to content

PeptideMapper

Jump to bottom
Harald Barsnes edited this page Jan 26, 2022 · 129 revisions

PeptideMapper - Efficient Amino Acid Sequence Mapping

PeptideMapper is a lightweight and efficient application for the mapping of amino acid sequences onto protein databases in the FASTA format. It takes advantage of the latest advances in string pattern matching to create a self-contained index that can be rapidly queried.

We have implemented queries for peptide sequences as well as sequence tags (alternating series of amino acids and mass gaps as produced by de novo sequencing algorithms).

PeptideMapper can be used standalone in command line or integrated as a dependency in other applications.

When using PeptideMapper, please cite Kopczynski et al., Bioinformatics, 2017.

For command line options and format specifications, please refer to PeptideMapperCLI. For use as a dependency, please refer to Usage As Dependency.

Methods

Database Index Creation

The first step is to prepare the protein sequences by creating one string with the pattern: /protein1/protein2/…/proteinn/ where proteini is the ith protein sequence in the given database and '/' a separation character preventing the mapping of peptides on two consecutive proteins. Accordingly, we compute the suffix array (SA) [1] of the complete proteome using divsufsort, chosen for its performance. Subsequently, the Burrows and Wheeler transform (BWT) [2] is computed and stored in a wavelet tree (WT) data structure [3] for fast rank queries.

Only the 22 amino acid characters, the four characters for amino acid combinations and '/' complete the alphabet Ε for the WT. To save memory, the WT is constructed according to a corresponding entropy-coded tree [4] built from the number of occurrences of every character. Additionally, the accession numbers extracted from the header are stored for every protein.

We implemented an efficient function querying all character occurrences within a substring in O(σ) where σ_=|Ε|. The SA is necessary to retrieve the starting positions of the mapped peptides within a protein. To save memory, we sample the SA and can retrieve the missing entries by using the last-to-front (LF)-mapping function [5] on the BWT.

For both peptide and sequence tag mapping, we utilize the backward search [5] which finds all occurrences of a pattern P in a text T in O(p) where p = |P|.

Peptide Mapping

To manage ambiguous amino acids and amino acids with indistinguishable mass, we first split peptides into their constituent amino acids and create a set of possible amino acids at each position. E.g., ATDWIRK => {A,X}{T,X}{D,B,X}{W,X}{I,L,J,X}{R,X}{K,X}. Theoretically, for this short peptide, we would need to map all possible 384 combinations, but using a dynamic programming approach we omit mapping of combinations containing prefixes we previously discarded.

The result is a set of ranges within the suffix array indicating matches that can be back calculated to the actual positions in the corresponding proteins. Having the position, we retrieve the accession number that completes all necessary data for further processing.

Sequence Tag Mapping

Sequence Tags are defined as alternating series of masses and sequences of amino acids. They can be inferred from spectra using sequence tag algorithms or de novo sequencing tools. A simple example is a tag of the following pattern m1-s-m2, where mi is a mass and s a partial peptide sequence. Having mapped s as described previously, we try to map the mass m1. Therefore, we extend s with all possible preceding sequences in the database creating a set of sequences s'.

For every sequence in s', we store a nominal mass based on the individual amino acid masses and fixed modifications provided. We compute alternative possible masses based on the variable modifications. A sequence in s' is kept when reaching m1 making a new set of sequences s", and discarded if the nominal mass exceeds m1 plus a given tolerance mass. Subsequently, we extend every sequence in s" with its possible following amino acids as done previously until reaching m2. Since this extension is conducted in the opposite direction, we also compute an index for the reversed protein sequences. Having extended s to all sequences matching the tag, all proteins it appears in as well as the corresponding starting positions are returned.

Performance

For the following benchmark tests, we are using an Intel(R) Xeon(R) 2.80 GHz quad core desktop computer with 16 GB RAM. Note that only one core is used for the performance tests, since the task of mapping independent peptides against a sequence database can be heavily parallelized, the computation time can almost be divided by the number of cores used. We compared both compomics-utilities index methods, the ProteinTree and the new adapted FM-Index. Additionally, we compared the sequence tag mapping with TagRecon [6].

Protein Sequences Databases and Benchmark Datasets

Yeast (May 2016):

  • size: 6,721 Proteins, 3,025,143 residues, 0 X residues
  • fixed modifications: none
  • variable modifications: none
  • fragment tolerance: 0.02 Da

Mouse (May 2016):

  • size: 16,813 Proteins, 9,474,758 residues, 79 X residues
  • fixed modifications: Carbamidomethylation of C
  • variable modifications: none
  • fragment tolerance: 0.02 Da

Human (July 2015):

  • size: 20,207 Proteins, 11,326,153 residues, 670 X residues
  • fixed modifications: Carbamidomethylation of C, Oxidation of M
  • variable modifications: none
  • fragment tolerance: 0.02 Da

Proteogenomics (March 2015) [7], [8]:

  • size: 83,721 Proteins, 13,851,427 residues, 0 X residues
  • fixed modifications: Carbamidomethylation of C, Oxidation of M, Acetylation of K
  • variable modifications: Acetylation of peptide N-term, Pyrolidone from Q, Pyrolidone from E
  • fragment tolerance: 0.5 Da

Metaproteomics (January 2013) [9]:

  • size 55,152 Proteins, 100,955,085 residues, 2,561,698 X residues
  • fixed modifications: Carbamidomethylation of C
  • variable modifications: Oxidation of M, Pyrolidone from Q, Acetylation of peptide N-term
  • fragment tolerance: 0.2 Da

All UniProt proteomes concatenated (July 2016):

  • size: 551,705 Proteins, 197,114,987 residues, 8,027 X residues
  • fixed modifications: Carbamidomethylation of C, Oxidation of M
  • variable modifications: Acetylation of K
  • fragment tolerance: 0.02 Da

For every database listed above, benchmark data sets D = {D1, D2, D3, D4} of different sizes (|D1| = 1,000, |D2| = 10,000, |D3| = 100,000, |D4| = 1,000,000) were created by extracting peptides and tags at random positions in the database. All benchmark data (proteomes, data sets, parameter files) are available here.

Results - Index Creation

Database time [s] size [MB]
ProteinTree FM-Index ProteinTree FM-Index
Yeast 65.73 2.55 184.3 7.02
Mouse 212.9 6.57 601.5 22.02
Human 252.0 7.49 537.9 26.33
Proteogenomics 526.2 8.81 1211 32.38
Metaproteomics 2956 58.00 5423 238.0
All Proteomes 7924 122.0 8653 458.0

Results - Peptide Sequences Mapping Time

Database D1 [s] D2 [s] D3 [s] D4 [s]
ProteinTree FM-Index ProteinTree FM-Index ProteinTree FM-Index ProteinTree FM-Index
Yeast 6.067 0.041 16.96 0.206 22.45 1.385 72.01 12.01
Mouse 10.17 0.057 22.48 0.313 35.93 2.291 90.23 19.49
Human 8.260 0.056 25.85 0.403 123.1 2.297 o.o.m. 21.35
Proteogenomics 8.910 0.058 39.26 0.262 259.2 2.094 o.o.m. 17.71
Metaproteomics 70.72 0.224 1183 1.547 9015 13.49 o.o.m. 134.5
All Proteomes o.o.m. 0.119 o.o.m. 0.643 o.o.m. 6.270 o.o.m. 62.20

Results - Sequence Tags Mapping Time

Database D1 [s] D2 [s] D3 [s] D4 [s]
TagRecon ProteinTree FM-Index TagRecon ProteinTree FM-Index TagRecon ProteinTree FM-Index TagRecon ProteinTree FM-Index
Yeast 238.0 4.670 0.096 279.0 20.64 0.590 281.0 38.12 3.863 288.0 205.9 33.59
Mouse 915.0 7.140 0.167 874.0 29.48 1.196 1002 107.3 9.910 1226 871.5 97.10
Human 1063 8.030 0.178 1151 31.76 1.352 1328 146.0 11.67 1608 o.o.m. 114.4
Proteogenomics PTMs not supported 9.850 0.297 PTMs not supported 40.11 1.546 PTMs not supported 234.7 12.16 PTMs not supported o.o.m. 119.1
Metaproteomics PTMs not supported 63.93 6.281 PTMs not supported o.o.m. 61.37 PTMs not supported o.o.m. 609.3 PTMs not supported o.o.m. 6205
All Proteomes 21981 o.o.m. 1.417 37791 o.o.m. 12.85 37353 o.o.m. 112.8 55171 o.o.m. 1145

o.o.m. = out of memory

Quick Start

To run PeptideMapper, download the packed zip archive of the latest compomics-utilities library, version 5.0.39 or newer.

wget http://genesis.ugent.be/maven2/com/compomics/utilities/5.0.39/utilities-5.0.39.zip
unzip utilities-5.0.39.zip
cd utilities-5.0.39
java -cp utilities-5.0.39.jar com.compomics.cli.peptide_mapper.PeptideMapperCLI -p exampleFiles/PeptideMapping/yeast.fasta exampleFiles/PeptideMapping/yeast-pep-1k.csv results.csv

A detailed description of the command line instructions is available here.

Troubleshooting

Should you encounter any issue with the usage of PeptideMapper, please create a new entry in the issue tracker.

References

[1] Manber and Myers, 1st Annual ACM-SIAM Symposium on Discrete Algorithms, 1990
[2] Burrows and Wheeler, Technical report, Digital Equipment Corporation, 1994
[3] Grossi et al., Proceedings of the Fourteenth Annual ACM-SIAM Symposium on Discrete Algorithms, 2003
[4] Huffman, Proceedings of the Institute of Radio Engineers, 1952
[5] Ferragina and Manzini, Proceedings of the 41st Annual Symposium on Foundations of Computer Science, 2000
[6] Dasari et al., Journal of Proteome Research, 2010
[7] Menschaert et al., Molecular & Cellular Proteomics, 2013
[8] Crappé et al., Nucleic Acids Research, 2014
[9] Tanca et al., PLoS ONE, 2013

Clone this wiki locally