Unsupervised Building and Exploitation of Composite Descriptors
Generally, the present invention provides a way of determining in an unsupervised manner additional members for a family that is defined initially through exemplar sequences. The present invention is unsupervised in that it proceeds without any information related to the exemplar sequences defining the family, without aligning the sequences, without prior knowledge of any patterns in the exemplar sequences, and without knowledge of the cardinality or characteristics of any features that may be present in the exemplar sequences. In one aspect of the invention, a method is used to take a set of unaligned sequences and discover several of many patterns common to some or all of the sequences. These patterns can then be used to determine if candidate sequences are members of the family. In another aspect of the invention, a method is used to take a set of sequences and to determine a set of maximal patterns common to a number of sequences The maximal patterns are determined without any previous knowledge about any properties of features that may be present in the processed sequences.
Latest IBM Patents:
- EFFICIENT RANDOM MASKING OF VALUES WHILE MAINTAINING THEIR SIGN UNDER FULLY HOMOMORPHIC ENCRYPTION (FHE)
- MONITORING TRANSFORMER CONDITIONS IN A POWER DISTRIBUTION SYSTEM
- FUSED MULTIPLY-ADD LOGIC TO PROCESS INPUT OPERANDS INCLUDING FLOATING-POINT VALUES AND INTEGER VALUES
- Thermally activated retractable EMC protection
- Natural language to structured query generation via paraphrasing
This application is a divisional application of U.S. patent application Ser. No. 09/712,638, filed Nov. 14, 2000 incorporated by reference herein.
FIELD OF THE INVENTIONThe present invention relates to sequences of symbols and, more particularly, to unsupervised building and exploitation of composite descriptors.
BACKGROUND OF THE INVENTIONSequences of symbols are useful in a number of areas. One such area is DNA. DNA (deoxyribonucleic acid) may be described through a long sequence of symbols. DNA is commonly described through the characters A, G, C, or T. These characters may be thought of as the alphabet of DNA. Another area where sequences of symbols are important is proteins. Proteins are sequences of amino acids, where each amino acid can be described by a character or letter The “alphabet” of amino acids comprises the character's of A, C, D, E, F, G, H, I, K, L, M, N, P, Q, R, S, T, V, W, and Y Sequences of symbols are also important in encryption and coding. For example, computers commonly store character data in numeric format. For instance, the word “the” could be coded in the American Standard Code for Information Interchange (ASCII) format as decimal symbols 116, 104, and 101. Encryption schemes change these numbers to conceal the underlying information.
For amino acids, there are very large databases of knowledge that consist of sequences of proteins. Similar proteins are usually grouped into “families.” Family members should have the same properties associated with them; once the properties of one of the family members is known, it is assumed that the other family members will have similar properties. Additionally, once the family is known, the family may be used to determine which candidate proteins are members of the family Therefore, there has been tremendous research to determine how to best group proteins into families.
Generally, there are four different methods used to group proteins. One method is to determine a pattern of symbols that all of the sequences share. This is called a single descriptor approach, which looks for particular patterns of characters. The patterns are series of expected amino acids, described by alphabetic characters. In the pattern, some locations could be important and some locations might not be. An example pattern for a single descriptor might require certain amino acids to be in one particular location, then allow several “don't care” locations where any amino acid could reside, and then require only a particular amino acid in a final location. The patterns are based on observations that, in nature, specific amino acid positions seem to be preserved in a biased way. These specific amino acids positions are “conserved” even though their neighbors can undergo mutations. Thus, researchers used the concept of conservation to describe the members of the family. A very large, well known database of the single descriptor type is the Prosite database. There are about 1100 families in this database. To find the patterns contained in each family, the proteins contained there were first aligned. Then, the most conserved region of the family was located and the pattern (the single descriptor) contained in all or most of the family members was determined. However, there could be members of a family that did not share the single descriptor. This generates false negatives, as members of the family were incorrectly not discovered as such.
An improvement on the single descriptor method is the composite descriptor method. The composite descriptor method examines a candidate protein for several alphabetic patterns, as opposed to only one pattern with the single descriptor method. Again, this method generally requires aligning the proteins so that the multiple patterns, i.e., the composite descriptor, properly align within their respective blocks.
The conceptual underpinnings are the same across all the methods that rely on composite descriptors. Any differences have essentially to do with either the manner in which multiple alignments are used to construct the descriptors or whether the descriptors are explicitly (e.g., a “regular expression”) or implicitly (e.g., a “profile”) represented in the composite description. Additional characteristics common to these approaches include: (a) an iterative component; (b) the availability of a set of known (or alleged) family members (=“training” set) that provides an initial “bootstrapping” stage; (c) the computation of a multiple-sequence alignment involving members of the training set—these alignments are typically verified manually or semi-automatically and can be used to derive profiles that allow the generation of quality measures when evaluating the results; (d) a range of quality control checks that are optionally applied on the generated results; and, (e) the need to study the collection under consideration in order to identify a minimum set of components that will form the composite description.
There are several problems with these approaches. For instance, in step (c), it is implicitly assumed that there is a multiple-sequence alignment involving all of the members of the training set; the alignment may either be a global alignment of both conserved and non-conserved regions, or a local alignment of the most conserved regions. This requirement unnecessarily burdens these methods. Additionally, multiple alignment programs usually work best when the parameters are optimized for the set of sequences which are being considered.
Steps (d) and (e) presuppose the availability of biological information pertaining to the set under consideration, and this biological information may not always be present. As a matter of fact, step (e) results in the selection and use of features which are conditional on each other. Although easy to describe, an additional assumption here is that the identity, cardinality, and properties of these features are available and also agreed upon ahead of time. For example, a statement such as “G protein-coupled receptors (GPCRs) are proteins involved in signal transduction in eukaryotic organisms that consist of seven transmembrane helices composed typically of hydrophobic amino acids” represents a body of knowledge that has been used by researchers in the building of composite descriptors for GPCRs. With the supervised approaches described above, a detailed and frequently manual study of the collection under consideration is unavoidable.
In addition to descriptor approaches, there are also “windowing” approaches that build descriptors for a family. In these methods, one or more windows ale used instead of character patterns. A single window method is called the PROFILE approach. All of the sequences of each of the family members are aligned with respect to their best-conserved region. Researchers then determined a probability distribution for locations in each column of the implied window. For each such block, they determined a probability of expecting an amino acid at some location within the window and thus built a ‘profile’ of expected probabilities for each of the columns of the window. The researchers would slide this set of probabilities against an unknown protein. If this candidate protein matched the expected probabilities, they included the protein as a member of the family. This approach was more tolerant than the single descriptor approach. Subsequently, researchers began to use profiles for multiple windows. There could be two, three, four windows where the members of the family could agree on content. Sometimes, a profile was not built explicitly but rather was maintained as a collection of the instances across the known or alleged family members of the conserved legion under consideration.
The windowing methods again rely on alignment of proteins, which can be relatively complex and computationally lengthy Typically, these windowing methods are supervised and biological information pertaining to the family can facilitate the analysis. With supervised approaches, a detailed and frequently manual study of the collection under consideration is unavoidable
Therefore, there exists a need to provide a way of determining and using family members of sequences in an unsupervised manner; without knowledge of biological information related to the family, and without aligning the sequences.
SUMMARY OF THE INVENTIONGenerally, the present invention provides a way of determining in an unsupervised manner additional members for a family that is defined initially through exemplar sequences. The present invention is unsupervised in that it proceeds without any information related to the exemplar sequences defining the family, without aligning the exemplary sequences, without prior knowledge of any patterns in the exemplar sequences, and without knowledge of the cardinality or characteristics of any features that may be present in the exemplar sequences. The cardinality of a set is the number of items in a set. For instance, the cardinality of the set of letters in the English alphabet is 26. In one aspect of the invention, a method is used to take a set of unaligned sequences and discover several or many patterns common to some or all of the sequences. These patterns can then be used to determine if candidate sequences are members of the family. In another aspect of the invention, a method is used to take a set of sequences and to determine a set of maximal patterns common to a number of sequences. The maximal patterns are determined without any previous knowledge about any properties or features that may be present in the processed sequences.
A mote complete understanding of the present invention, as well as further features and advantages of the present invention, will be obtained by reference to the following detailed description and drawings.
BRIEF DESCRIPTION OF THE DRAWINGS
Generally, the present invention provides a way of determining in an unsupervised manner additional members for a family that is defined initially through exemplar sequences. The present invention is unsupervised in that it proceeds without any information related to the exemplar sequences defining the family, without aligning the sequences, without prior knowledge of any patterns in the exemplar sequences, and without knowledge of the cardinality or characteristics of any features that may be present in the exemplary sequences. The cardinality of a set is the number of items in a set. For instance, the cardinality of the set of letters in the English alphabet is 26. In one aspect of the invention, a method is used to take a set of unaligned sequences and discover several or many patterns common to some or all of the sequences. These patterns can then be used to determine if candidate sequences are members of the family. In another aspect of the invention, a method is used to take a set of sequences and to determine a set of maximal patterns common to a number of sequences. The maximal patterns are determined without any previous knowledge about any properties or features that may be present in the processed sequences.
As previously stated, the present invention provides a way of determining family members in an unsupervised manner. By “unsupervised” it is meant that no predetermined or a priori information is needed/known about the exemplar sequences or is employed by the discovery process. Additionally, there is no need for user supervision or intervention. For instance, the present invention does not require knowledge of biological information related to the family, aligned sequences, knowledge of properties of the exemplary sequences defining the family, and/or knowledge of the cardinality or characteristics of features of the exemplar sequences. It is possible to exclude one or more of the these restrictions. For instance, the present invention could be used on a set of aligned sequences. The present invention would still determine a composite descriptor suitable for examining candidate sequences and either including these sequences in or excluding them from the family. However, a great benefit of the present invention is that it does not need aligned sequences or the knowledge of predetermined properties and features that may be present in the exemplar sequences. Aligning sequences and determining properties and features in the exemplar sequences that originally define the family is time consuming, complex, and at times intractable. Instead, the present invention can determine a composite descriptor without such time intensive efforts.
Concerning features and properties of a sequence of symbols, it is not easy to define what a feature is. The definition of a feature is directly related to the representation of the items that are studied, i.e., the way each of the objects processed by the system under consideration is represented and stored in a computer. Such a representation is in turn related to the way an object can appear in the context of the sensor data, and is unavoidably application specific. For example, in the context of image processing by a computer, the following image characteristics have been used as features: linear and curvilinear segments, curvature extrema, curvature discontinuities, and identifiable conics. In the context of computational biology, an example of a feature can be a combination of amino acids with understood behavior and possibly known 3-dimensional structure. For instance, for a helix-turn-helix (HTH) motif that mediates the binding of many regulatory proteins to regulatory control sites of DNA, the two features are the two helices at the beginning (7 a.a.) and the end (9 a.a.) of the 20 a.a. stretch that corresponds to an instance of the HTH motif. A property can be thought of as an attribute of a feature: in the case of the HTH, a property would be the tact that the two features (helices) are held together through non-polar interactions of their side chains. It should be stressed that the concept of the feature is also intrinsically connected to the task at hand. For example, for some applications, individual a.a. letters can be thought of as “features.”
What is important is that previously researchers had to (a) know something about the set of sequences, or (b) align the exemplar sequences, or (c) perform both (a) and (b) before they could determine those motifs that were peculiar to the exemplar sequences and, thus, by extension specific to and characteristic of the family defined by the exemplar sequences. The researchers knew and exploited properties of sequences, knew and exploited features of the sequences, and/or aligned the sequences. The present invention is unsupervised, meaning that no information about the exemplar sequences need be known, and the present invention will still determine patterns that can subsequently be used to define the family implied by the exemplar sequences as well as analyze candidate sequences for inclusion into this family.
In an embodiment of the present invention, a training set of family members is searched in an unsupervised manner to determine statistically significant, common patterns between some or all of the family members. Each family member comprises a sequence, which itself comprises a series of characters. The present invention may be used on any sequence of symbols that can be described as a linear stream of events, e.g., DNA (deoxyribonucleic acid), proteins, languages, and numbers. Preferably, a predetermined sequence-support threshold will initially be set. This predetermined sequence-support threshold determines how many of the sequences in the family need to have a pattern for the pattern to be considered common to the training set. For instance, if there are 100 sequences in the family, the predetermined sequence-support threshold could be set to 50. This means that a pattern must be found in 50 of the sequences for the pattern to be considered common to the family members in the training set. Generally, this threshold is initially set to the number of sequences in the training set. Should no common patterns be found, the sequence threshold may be modified.
If common patterns are found, they are examined to determine if they are statistically significant. Any remaining statistically significant patterns may be used to describe the family members and, subsequently, to ascertain if candidate members are part of the family. Preferably, the statistically significant and common patterns become part of a composite descriptor. Once the statistically significant and common patterns are found for a set (which could include all) of the family members, the sequences containing the patterns are removed from the training set. This results in a smaller training set.
This modified training set is again searched for common patterns. The sequence threshold may be modified to search for fewer sequences of the modified training set or to search for all of the sequences in the training set. If any common and statistically significant patterns are found, the composite descriptor is modified to add the new patterns. This process preferably continues until either all sequences are removed from the training set or until common patterns cannot be found between the remaining sequences.
Once the composite descriptor is determined, the composite descriptor may be used to determine if a candidate sequence is past of the family. In particular, the composite descriptor may be used to search a database of sequences to determine if individual sequences in the database are members of the family described by the composite descriptor. Usually, a pattern-support threshold will be used to make this determination. The pattern-support threshold determines the number of patterns that must match between the candidate sequence and the patterns in the composite descriptor. For example, if there are 1000 patterns in the composite descriptor, the pattern-support threshold may require matches on 995 of the patterns for the candidate sequence to be considered a member of the family. Moreover, after more members of the family are found by using the current composite descriptor, these new members may be added to the original training set to create a new training set. The composite descriptor method may again be run on the new training set. This will provide even greater sensitivity and allow the composite descriptor to “learn” new patterns common to the family.
While the present invention can determine statistically significant and common patterns with aligned sequences, the present invention does not need aligned sequences. To align two sequences, one or more patterns common to both sequences are aligned in a left-to-right order. For example, assume that the pattern being aligned is ABC. The sequence of characters {DEFXYZABC} would be aligned with {ABCDEF} by either aligning the ABC patterns in a left-to-right manner or by aligning the DEF patterns. Thus, when aligning the ABC patterns, the XYZ of the first pattern would not be aligned with characters in the second pattern and the DEF of the second pattern would not align with characters in the first pattern. For this example, there is no unique alignment and it is easy to see how the situation can be complicated further as the number of sequences to process increases. Because the present invention preferably searches for patterns common to the sequences, the present invention would determine that ABC was common to the two sequences, regardless of their alignment.
The present invention also does not need the availability of biological information related to the family. While such information could be used, the present invention will determine statistically significant and common patterns within the family members without biological information. Moreover, because outliers are expected to not contribute much in the way of statistically significant patterns to the composite descriptor, outliers have less of an impact on the present invention.
Turning now to
As is known in the art, composite descriptor method 200 may be distributed as an article of manufacture that itself comprises a computer readable medium having computer readable code means embodied thereon. The computer readable program code means is operable, in conjunction with a computer system such as computer system 100, to carry out all or some of the steps to perform the composite descriptor method 200. The computer readable medium may be a recordable medium (e.g., floppy disks, hard drives, Compact Disks, or memory sticks), or may be a transmission medium (e.g., a network comprising fiber-optics, the world-wide web, cables, or a wireless channel using time-division multiple access, code-division multiple access, or other radio-frequency channel). Any medium known or developed that can store information may be used.
Composite descriptor method 200, as shown in
Method 200 begins in step 205 when a training set is provided. It should be noted that the sequence of steps are not necessarily in order. The training set, T, is preferably N unaligned sequences si for which there is reason to believe that the sequences are related. There should exist identifiable local similarities among members of T at the amino acid level, although it is assumed that no other information is available for the members of T, e.g., known or identifiable secondary structures, known or identifiable domains, functional information, physio-chemical properties, or physical properties. If no identifiable local similarities exist among members of T, method 200 will not provide a suitable composite descriptor for the family, as a composite descriptor does not exist for the family.
Each sequence is a series of symbols from an alphabet. For proteins, one can denote by Σ the alphabet of all amino acids; i.e., Σ={A, C, D, E, F, G, H, I, K, L, M, N, P, Q, R, S, T, V, W, Y}. On this alphabet, regular expressions can be defined that can range from very simple n-grams to more general ones containing wild cards and capturing strings of variable length. The ‘ ’ (referred to as the “don't care character”) is used to denote a position in a sequence or pattern that can be occupied by an arbitrary residue. A bracket is meant to denote a “one of” choice; i.e., [KR] means that the position this bracket corresponds to can be occupied by exactly one of K or R. A bracket can have a minimum of 2 alphabet characters but not more than |Σ|-1.
In step 210, the sequence threshold, K, is set. It is possible to set K=|T|, which is the number of sequences in the training set. In actuality, it has proven beneficial to assign a small starting value to K that is a fraction of the number of sequences in T. Experiments have shown that a starting value of K=|T|/b with b=4 or 5 is a good choice across many data sets. Note that the smaller the value of b, the higher the redundancy of the composite descriptor will be. The selection of K also can depend on how conserved, or similar, the family members are. If the family members are well conserved, then K can be higher; if the family members are not well conserved, then K can be lower.
In step 215, a set of maximal patterns in the K sequences is determined. In general, this step tries to determine common patterns between the K sequences. Not only should the patterns be common, but they should also be as large as possible. These large patterns may further be mathematically defined as “maximal” in a way described below. Any of the available algorithms which can guarantee that all sought patterns are discovered and that they are maximal can be used here. For the experiments related below, a Teiresias algorithm was used. This algorithm is described in floratos, et al, U.S. Pat. No. 6,108,666, “Method and Apparatus for Pattern Discovery in 1-Dimensional Systems”; Flotatos, et al, U.S. Pat. No. 6,092,065, “Method and Apparatus for Discovery, Clustering and Classification of Patterns in 1-Dimensional Event Streams”; Rigoutsos, I and A. Floratos, “Combinatorial Pattern Discovery in Biological Sequences: the Teiresias Algorithm,” Bioinformatics, 14(1):55-67, 1998; and Rigoutsos, I. and A. Floratos, “Motif Discovery Without Alignment Or Enumeration,” Proceedings 2nd Annual ACM International Conference on Computational Molecular Biology, New York, N.Y., March 1998, the disclosures of which are incorporated by reference herein.
A short introduction to this method follows. A pattern S is a regular expression on Σ that defines a language G(S). The elements of the language are all the strings that can be obtained from the regular expression that S stands for. A protein is said to match a given pattern S if and only if it contains at least one substring (i.e., a block of consecutive residues) that belongs in G(S). A pattern S′ is said to be more specific than a pattern S if G(S′)⊂G(S). Given a pattern S and a database D, an offset list of a pattern of S may be defined with respect to D (or simply the offset list of S, when the database D is unambiguously implied) to be the following set: LD(S)={(i, j)| the i-th sequence of the database D matches the pattern S at offset j}. A pattern S is called maximal with respect to a database D if there exists no pattern S′ which is more specific than S and such that |LD(S)|=|LD(S′)|. A maximal pattern cannot be made mole specific without simultaneously reducing the cardinality of its offset list. A pattern S is called an <L,W> pattern (with L≦W) if every substring of S with length W contains L or more non-don't care positions. Note that a given choice for the parameters L and W has a direct beating on the degree of remaining similarity among the instances of the domain that is captured by the regular expression: the smaller the value of the ratio L/W, the higher the degree of sought similarity.
The Teiresias algorithm is a pattern discovery algorithm that can guarantee the discovery of all <L,W> patterns that are maximal and supported by K or more input sequences The pattern discovery is carried out while allowing the symbols of Σ to be partitioned in equivalence classes. Any symbol within a given class is able to replace any other symbol of the (same) class. One such example would be the partition: {A, G}, C, {D, E}, {F, Y}, H, {I, L, M, V}, {K, R}, {N, Q}, P, {S, T}, W. In fact, the various symbol classes do not have to form a partition of Σ. In other words, a given symbol can belong to more than one class. One such set of classes can be obtained by using a distance threshold with any of the currently available scoring matrices such as the PAM and BLOSUM series. PAM is described in Dayhoff, “Atlas of Protein Sequence and Structure,” vol. 5, National Biomedical Research Foundation, 1978; and BLOSUM is described in Henikoff, “Amino Acid Substitution Matrices from Protein Blocks,” Proc. Natl. Acad. Sci. USA, 89:100915-100919, 1992, the disclosures of which are incorporated by reference.
The Teiresias algorithm permits the discovery of all <L,W> patterns that are maximal and supported by K or more input sequences, in the presence of stated equivalences involving symbols from the input alphabet. Each pattern S that the Teiresias algorithm will discover is of the form:
(Σ∪[ΣΣ*Σ])(Σ∪[ΣΣ*Σ]∪{ })*(Σ∪[ΣΣ*Σ])
Associated with each pattern S is the sensitivity of the pattern, which is directly related to the number of sequences in D that contain S. The sensitivity is a measure of how many members of the training set T do not match S (=false negatives). Also associated with S is the pattern's specificity, which is a direct measure of how many members of the database D match the pattern, but are not true members of the collection that the training set T represents (=false positives). The choice of the values for the parameters L and W is a function of the collection under consideration. Experimental work has shown that a choice supporting moderate degree of local similarities (e.g., ˜40-50%) is a good choice across a very large variety of test cases.
In step 225, it is determined if any patterns are found. In no patterns are found (step 225=NO), the sequence threshold, T, can be decreased. Preferably, this is done by setting K=|T|/b, where b is usually set to 4 or 5. It is also possible to set b to smaller values, such as 2 or 3. Setting b to smaller values increases the amount of processing time it might take to determine maximal patterns. For instance, if there are 1000 sequences in T and K=|T|=1000, and no common maximal patterns are found, it is necessarily the case that changing K to 999 will not find any common maximal patterns. Changing K from 1000 to 250, however, will make it more likely that common maximal patterns may be found. After K has been changed (step 230), it is determined if K meets a predetermined minimum limit. This limit has been set, in the example of
If maximal patterns are found in step 215 (step 225=YES), in step 235, it is determined if the maximal patterns are statistically significant. In general, in step 235, it is determined, for each maximal pattern, what the probability is that the maximal pattern occurs in a sequence. This probability should meet a predetermined threshold. This step is important because the patterns will be exploited, as part of the composite descriptor, to determine additional family members. If relatively general patterns are used, the patterns could include candidate members into a family when the candidate members are not members of the family. For instance, for the English language, the pattern “the” is much more likely to appear in a sentence than is the pattern “quit.” The pattern “the” would be much more likely to include candidate members as part of the family than would the pattern “quit.” This would be appropriate if the family was defined as any sentence having the pattern “the.” However, a much more likely occurrence is to define a sentence as any sentence having the pattern “quit,” and if the pattern “the” is used as part of a composite descriptor, it is possible that this pattern will generate too many false family members.
From the set of maximal <L,W> patterns that are discovered, the set MS is selected that contains only those that are statistically significant. With appropriate modifications, any of several published methods can be used at this step, the disclosures of which are herein incorporated by reference: Atteson, “Calculating the Exact Probability of Language-like Patterns in Biomolecular Sequences,” Proceedings of the Sixth International Conference on Intelligent Systems for Molecular Biology (ISMB '98), Menlo Park, Calif., AAAI Press, 1998; Jonassen, I Collins, and Higgins, “Finding Flexible Patterns in Unaligned Protein Sequences,” Protein Science, pp 1587-1595, 1995; Nicodeme, Salvy and Flajolet, “Motif Statistics,” INRIA Technical Report No 3606, January 1999; Pevzner, Botodovksi and Mironov, “Linguistic of Nucleotide Sequences: the Significance of Deviation from Mean Statistical Characteristics and Prediction of the Frequencies of Occurrences of Words,” Journal of Biomolecular Structure Dyn., 6:1013-1026, 1989; Regnier, “A Unified Approach to Word Statistics,” Proceedings 2nd Annual ACM International Conference on Computational Molecular Biology, New York, N.Y., March 1998; Sagot, and Viari, “A Double Combinatorial Approach to Discovering Patterns in Biological Sequences,” Proceedings of the Seventh Symposium on Combinatorial Pattern Matching, pp. 186-208, 1996; Sewell and Durbin, “Method for Calculation of Probability of Matching a Bounded Regular Expression in a Random Data String,” Journal of Computational Biology, 2(1):25-31, 1995; and Wooton, “Evaluating the Effectiveness of Sequence Analysis Algorithms Using Measures of Relevant Information,” Computers Chem., 21(4):191-202, 1997.
For simplicity, the probabilities of the discovered patterns, as disclosed in the Examples section below, were determined with the help of a 2nd order Markov chain method, as described in Salzberg, Delcher, Kasif, and White, “Microbial gene identification using interpolated Markov models,” Nucleic Acids Res., 26(2):544-8, 1998, which is incorporated herein by reference. The natural logarithm of the estimated probability was used as the measure of a pattern's significance. This threshold can be estimated as a function of the size of the database to be searched with the composite descriptor.
The cardinality of the sub-selected set Ms of patterns ought to be high because of the redundancy of sequence segments from T that ate captured by the patterns. This will guarantee a strong signal-to-noise ratio when the composite descriptor is used as a predicate. It is worth pointing out at this point that even if the training set has just a few members, the cardinality of Ms (and thus the redundancy) can be high since there is a multitude of patterns that one can generate even from a few sequences.
Once the statistically significant patterns are found, these patterns are removed from the training set, T, of sequences. This occurs in step 240. It should be noted that steps 240, 245 and 250 do not have to occur in this order and could even occur in parallel. Preferably, each sequence of the training set is examined to determine whether it matches any of the significant patterns of Ms. After all patterns of Ms have been exhausted, all sequences that matched one or more patterns are added to a temporary set A. Upon completion of the iteration, one or mote sequences from T will have been entered into the set A; these are essentially the sequences that have been accounted for. What remains of T after the removal of these sequences, i.e., T\A, is used as the training set for the next iteration. Thus, the training set T is modified (step 245), which could include marking which sequences in an array of sequences are no longer valid, or copying the remaining sequences into a new array.
In step 250, the composite descriptor is modified. Preferably, the composite descriptor is a union of the composite descriptor and the set Ms. The set of significant patterns Ms which was discovered during this last iteration is added to the composite descriptor by adding those patterns in Ms that are currently not in the composite descriptor.
In step 255, it is determined if the training set, T, is empty. If the training set is not empty (step 255=NO), the method continues in step 215 and repeats. If the training set is empty (step 255=YES), and after step 220=YES, the method ends in step 280. Optionally, Steps 260, 270 and 275 may be performed at this point.
At the end of this stage, the composite descriptor contains a set of patterns that by design ate specific and sensitive for the collection that the training set T represents. Several properties distinguish this composite descriptor from previous collections of patterns, such as the Prints database of patterns. For example, the building of the composite descriptor is automatic, it does not requite manual intervention and does not necessitate the computation of multiple alignments. Additionally, there is no need for biological knowledge specific to the training set T that will impose helpful constraints during generation of the composite descriptor. Also, highly similar sequences need not be removed from the training set prior to the building of the composite descriptor. Additionally, as discussed below in reference to step 260, the training set can safely contain a small percentage of potential outliers, i.e, sequences that have questionable membership in the collection that the training set represents. Because of the redundant, iterative nature of the building phase, the resulting composite descriptor is not expected to contain any statistically significant patterns that are shared by both the outliers and the rest of the sequences in T. Through the initial selection of the support value (small K) the composite descriptor can be made sensitive and contain patterns that are specific for the set T (i.e., large probability threshold, Thrprob). Finally, the fact that the composite descriptor contains all those patterns which are specific, significant, and which by design account for every member of the training set, guarantees a strong signal-to-noise ratio when using composite descriptor as a multi-valued predicate (which takes place in step 260) Steps 205 through 255 may be expressed in pseudo-code as follows
In step 260, the composite descriptor is exploited to determine if candidate sequences are members of the family described by the composite descriptor. Generally, a database of sequences (such as database 150 of
Because method 200 relies on searching through a family and determining the common, statistically significant and maximal patterns that compose the composite descriptor, outliers tend not to matter as much for the present invention. An outlier is a sequence that has been erroneously included within the family. Some simple examples will help to explain why outliers are not a hindrance to the present invention.
Assume that there are 100 members of the family; assume also that 93 members of the family are accounted for but there are 7 outliers that were erroneously included as members of the family. Since, by definition, the latter set comprises the outliers, it is generally true that the number of patterns that will be shared among them and the remaining 93 sequences should be very small (if not 0) when compared to the number of patterns that will be shared by the 93 truly related sequences. This will thus generate very small (if any) support for sequences that are not true members of the family being studied. Moreover, these erroneous patterns will be further filtered out through the statistical significance filtering stage. Finally, when the composite descriptor, which contains patterns common to all 100 sequences, is used to determine if a new sequence is part of the family, the composite descriptor will be used with a pattern-support threshold. In other words, there will be some minimum number of patterns that the new sequence must have in order to be considered part of the family. This threshold will usually be high enough such that outliers, even if they contribute patterns, will not cause non-family members to be included within the family.
In step 260, the composite descriptor can be used as a multi-valued predicate that can determine the membership of a query sequence in the collection that the original training set T defines. The composite descriptor can be used to examine a candidate-for-membership-in-T sequence Scand for instances of the permitted patterns. Given Scand, as many local counters as the length of the sequence may be allocated and initialized to 0. A global counter for the sequence may also be allocated and also initialized to 0. If it is determined that a segment of the query matches a pattern m, the local counters at the sequence positions matching the pattern are incremented by an amount equal to d. The possible choices for d include among others “the number of occurrences om of m in T” and “the number 1.” The former choice favors segments that match patterns supported by a lot of sequences in T whereas the latter gives comparatively increased support to segments that are only moderately conserved. The choice for the amount d by which to increment the local counters modifies the semantics of the predicate's output value.
If the value of d is set to ‘1’ then the predicate is a measure of how many distinct patterns generated from T are matched by the query sequence. In this case, large values indicate that the result is corroborated by multiple patterns which are specific for the collection T. Smaller values are at the very minimum indicative of the existence of local similarities that are shared by the query and one or more members of the training set T. Such similarities can imply one of two things: either the query is a true but distant member of the collection under consideration or it is not a true member but it nonetheless shares one more regions of similarity with members of the collection.
If the value of d is set to ‘the number of occurrences’ of the respective pattern in the training set T, the predicate is a measure of how many distinct sequence fragments in T are similar to the respective query fragment. Large values indicate regions that are shared by a large number of sequences in T and can be indicative of a conserved active site, for example. Both choices of d have merit and the one to use depends on the task at hand.
Independent of what the choice for d is, every time a segment of Scand matches a pattern m, the global counter associated with Scand is incremented by d. After all of Scand have been examined, the values of the global counter are inspected for Scand; if they exceed Thresrand, Scand is reported as a candidate for membership in the collection defined by T.
The value of Threscand depends on the actual contents of the composite descriptor and can be determined as follows: beginning with the composite descriptor that was built from the training set T, one can scan as outlined above a randomized version of a very large database such as GenPept or Swiss-Prot. Essentially, each sequence of such a database is treated as a potential query. Upon completion of the scanning process, one can accumulate support for all the sequences that matched one or more patterns of the composite descriptor and histogram the support values to obtain their distribution. The value of Thresrand may be determined by identifying the q-th percentile of this last distribution. Typically, q is set to 95 or higher.
After step 260 has been performed, it is possible to take the new members found and add them to a new training set that comprises the old training set and the new members. Then steps 205 through 280 may be run again (step 275) to further refine the composite descriptor for this fully. Thus, the present invention allows learning to be performed, if this is desired.
The present method does not suffer from drawbacks related to (a) the need for good multiple sequence alignments, (b) the inclusion of outliers, (c) the inherent dependence of the results on the selection of the scoring matrix that is used, and (d) overtraining. Indeed, building of the composite descriptor does not require the computation of any multiple sequence alignments, whereas the redundancy of representation that is inherent in composite descriptor is expected to more than counterbalance the inclusion of any small number of outliers. Additionally, this will prevent the system from including even more outliers during the following iteration. Moreover, after each iteration, only the sequence fragments whose support exceeds threshold ate considered thus allowing the process to remain ‘focused’ on what has been deemed important and relevant for the dataset under consideration.
Finally, it should be noted that the training set T, which is given at the very beginning of this iterative process, impacts on the quality of the results (i.e., sensitivity and specificity) that the method will produce. For example, if the original training set is not sufficiently representative of all instances of a family's members (e.g. GPCRs), or of the construct of interest (e.g. the helix-turn-helix DNA binding motif), the generated composite descriptor should not be expected to discover all instances relating to the training set. This last observation holds true for all methods that try to build single or composite descriptors by starting with a training set T. Since the augmented training sets at the beginning of the i+1-st iteration preferably only comprise the sequence fragments which exceeded threshold during the i-th iteration, the composite descriptor will maintain its ‘focus’ on what is essentially dictated by the original training set. That is not to say that that the composite descriptor will not be sensitive; on the contrary, the composite descriptor will be sensitive to the extent that the processed data permit while at the same time remaining in lock-step, so to speak, with the originally provided training input. As a matter of fact, the experimental results discussed below on three specific datasets demonstrate that even starting with small training sets allows discovery of a large number of representatives of the same group.
EXAMPLESNow that the method and apparatus have been described, some exemplary results are shown in this section. In this section, results are described from the building and use of composite descriptors for three distinct collections of data. The collections were chosen in such a way so as to showcase the ability of the present invention to handle input sets across a variety of contexts.
The first collection comprises sequences from PROSITE entry PS50040 of elongation factor T gamma chain sequences; in Release 15.0 of the PROSTIE database, only a matrix profile is available for this collection.
The second collection comprises complete sequences as well as fragments of G protein-coupled receptors, a very important and diverse family of proteins that has traditionally been used as a benchmark test for gauging the quality of pattern-based approaches.
Finally, the third collection comprises sequence fragments that are known to contain an instance of the helix-turn-helix DNA binding motif, a structural motif of great importance.
First, the composite descriptors were built for each of the three collections and evaluated by treating the sequences in Swiss-Prot Release 38.0 as candidates for membership in each of the respective three collections.
Once the behavior of the descriptors is characterized in the context of Swiss-Prot the 19,099 ORFs were searched in the complete genome of Caenorhabditis elegans and these results reported below.
Before proceeding, here are some methodological details and parameter choices that are common in all three cases. In particular, the value d, by which the counters are incremented, is set to 1, essentially favoring those sequences that contain more instances of distinct patterns over others. The value of Thrprob is determined by assuming that the patterns ought to be able to discriminate among sequences in a database as large as GenPept; although for a database of this size an estimated log-probability of −25 or less ought to suffice. Thus, the more stringent threshold of Thrprob=−30 was used with the understanding that this will result in a sacrifice in sensitivity. But as the results will demonstrate, even with this stringent threshold, the redundancy of each composite descriptor leads to a sensitivity that is satisfactory. Also, in all three cases the following a.a. equivalences are assumed: {A, G}, C, {D, E}, {F, Y}, H, {I, L, M, V}, {K, R}, {N, Q}, P, {S, T}, W.
The First Example EF1G/PS50040An application of the above described methodology is in the context of the PROSITE database. Although numerous entries in PROSITE contain succinct and specific patterns capturing most or all of the members of the corresponding collection, there exist entries for which only a profile/matrix is available: PS50040, the family of elongation factor 1 gamma chain proteins is one such example thus making it an ideal candidate for processing with the described method.
PS50040 comprises 10 full sequences (EF1G_ARTSA, EF1G_CAEEL, EF1G_HUMAN, EF1G_RABIT, EF1G_SCHPO, EF1G_TRYCR, EF1G_XENLA, EF1G_YEAST, EF1H_XENLA, EF1H_YEAST) and 1 fragment (EF1G_PIG). The reported profile matrix captures all 10 full sequences, misses the one fragment and generates no false positives when the target database is Swiss-Prot Rel 38.0.
It should be noted here that if one relaxes the constraints imposed by the chemical equivalence classes shown above, it is possible to discover a specific pattern that belongs to all 11 members of PS50040 and generates no false positives when used in conjunction with Swiss-Prot Rel 38.0. In fact, this pattern is
[ILMV]. . . [NW][ILMV]. . . [AG]. . . [RI][ILMV]. . . [KT]. . . F . . . [ILMV]. . . [GH]. . . [AG]
and can be used to describe and capture elongation factor 1 gamma chain proteins; the deviations from the above chemical equivalence classes awe shown in boldface.
The composite descriptor was built for this collection by setting the Teiresias parameters to L=5 and W=10; since the dataset is small there was only a single iteration over the dataset with a threshold choice of K=6. In other words, the composite descriptor was built by discovering patterns that involved a minimum of 5 non-wild cards in any rolling window that spans 10 positions and begins/ends with a literal, a relatively high-degree of local similarity (i.e. 50% or higher). Those patterns whose estimated log-probability was equal to −30.0 or less were selected and this generated a composite descriptor that comprised 2,260 patterns.
First, a corresponding DFA (deterministic finite automaton, which will only recognize instances of the composite descriptor patterns in a query sequence and which performs method step 260) was used to search a randomized version RAND-Swiss-Prot of Swiss-Prot (Release 38.0) that was obtained by applying a randomly chosen permutation to the amino acids of each of the valid sequences. Both the composition and lengths of individual sequences were maintained by this operation. The global counter for each randomized sequence was derived by summing up the local counters from each sequence region that received non-zero support. The sequences were then sorted in order of decreasing global-counter value. Twenty seven (27) randomized sequences received non-zero support with global counter values that ranged between 1 and 2 inclusive. Thresrand was thus set to 3, and the DFA was subsequently used to search the actual Swiss-Prot database. Of the 69 sequences that received non-zero support, only 16 exceeded the predefined threshold. The support values for the 16 sequences were: EF1G_HUMAN 861, EF1G_RABBIT 846, EF1G_XENLA 791, EF1H_XENLA 765, EF1G_ARTSA 349, EF1G_CAEEL 228, EF1G_YEAST 110, EF1G_SCHPO 110, EF1H_PIG 96, EF1H_YEAST 94, EF1G_TRYCR 88, SYV_FUGRU 7, GTT1_RAT 5, GTT1_MOUSE 5, SYEP_HUMAN 3 and GTH4_MAIZE 3.
Note that the 5 hits SYV_FUGRU, GTT1_RAT, GTT1_MOUSE, SYEP_HUMAN and GTH4_MAIZE are clearly separated from the 11 top scoring sequences. They however obtained scores which were above threshold and thus are studied in more detail. In all 5 cases, one or more sizeable regions that were shared with one or more members of the PS50040 collection were discovered The Clustal-W alignment of EF1G_XENLA and the N-terminus of SYV_FUGRU, a valyl-trna synthetase from Fugu rubripes, are shown in Table 1 below. Table 1 shows a Clustal-W alignment of EF1G_XENLA and the N-terminus of SYV_FUGRU, and this shows a strong similarity. As can be seen, the similarity among these two sequences is pretty extended and the Clustal-W score for the shown alignment equaled 462.
Similar shared legions are present in GTT1_RAT & GTT1_MOUSE (a glutathione s-transferees 5 from Rates norvegicus and a glutathione s-transferase theta 1 from Mus musculus respectively), SYEP_HUMAN (a multi-functional aminoacyl trna-synthetase frrom Homo sapiens) and GTH4_MAIZE (a glutathione s-transferase IV from Zea mays). The Clustal-W alignments for these cases are shown in Tables 2 through 4 below. Table 2 shows a Clustal-W alignment showing a substantial similarity between GTT1_RAT, GTT1_MOUSE and EF1G_ARTSA. The Clustal-W score is 1577. Table 3 shows a Clustal-W alignment between a fragment from EF1G_CAEEL (a.a 100 through 243) and a fragment from SYEP_HUMAN (a.a 1 through 180) showing a shared region The Clustal-W score for this alignment is 74 Table 4 shows a Clustal-W alignment showing a strong similarity between EF1G_RABIT and GTH4_MAIZE. The Clustal-W score is 215.
It should be noted that a search of MEDLINE has indicated that with the exception of the similarity between the EF1G family and the valyl-tRNA from Fugu rubripes, none of the other similarities shown here has been reported in the literature.
In summary, the composite descriptor has correctly picked out the members of PS50040 from the contents of Swiss-Prot as well as has identified several substantial similarities with other sequences in the database.
The family of G protein-coupled receptors has a long evolutionary history and is of particular importance for signal transduction in all eukaryotes. Spanning the lipid bilayer of the plasma membrane with seven helices, they bind and form signal transducing couples that are at the center of many key processes such as visual excitation, olfaction, histamine secretion in allergic reactions, and chemotaxis. G protein-coupled receptors form a very diverse family and extensive studies have shown that single descriptor approaches do not suffice to characterize the family's members.
Despite considerable efforts, very few membrane proteins have yielded high-resolution X-ray crystallographic data; this led to increased use of electron microscope approaches. The first such data were in fact obtained for bacteriorhodopsin, the bacterial analogue of rhodopsin, where a 3 Å electron-microscopy reconstruction of it has established directly the presence of the seven transmembrane helices. The significant sequence similarity that the members of this family exhibit indicates that they ought to have the same topology.
In order to demonstrate the power of the present invention and its ability to generalize, the experiment began with the contents of the GPCRDB as they existed in May 1998. Note that from this collection the hypothetical proteins from Caenorhabditis elegans are excluded since it was intended to carry out GPCR-discovery in this genome. The bacterial analogues of rhodopsin as well as all listed G-proteins were also excluded. What was left was a total of 1,019 GPCR entries, of which 862 were complete sequences and 157 were fragments. This set was intersected with an older release of Swiss-Prot (Release 35.0 from November 1997) and determined that the intersection of the two databases comprised a total of 804 sequences and fragments. Starting with data that were almost two years old was intentional since it was important that the ability of the composite descriptors to generalize and identify additional candidate sequences in the much larger databases of today would be shown.
The collection of 804 GPCR sequences and fragments contained several classes (e.g. thodopsin-like, secretin-like, pheromone, etc.) of proteins. In turn, each of these classes comprised several representatives. Instead of selecting representatives from each of the identified classes, the order of the sequences in this set of 804 members were randomized. Note that the contents of the sequences themselves remained unchanged, only their order of appearance was modified. For example, the 613-th sequence was now listed 4-th, the 11-th sequence now appeared in the 45-th position, and so on. Subsequently, a training set T was formed by collecting the sequences and fragments listed in the first 80 positions, arguably a very small set if one considers the diversity of the GPCR family. Essentially, slightly less than 1/10-th of the available dataset were randomly sub-selected for the purposes of building the composite descriptor. Table 5 below contains a listing of the labels of the 80 sequences in this training set. Table 5 shows the Swiss-Prot labels of the 80 sequences in the training set for the G protein-coupled receptor experiment. The labels are listed in the order they were selected and they correspond to both sequences and sequence fragments.
As in the previous example, the patterns were discovered assuming the equivalence classes {A, G}, C, {D, E}, {F, Y}, H, {I, L, M, V}, {K, R}, {N, Q}, P, {S, T}, W. The Teiresias parameters were set to L=5, W=10, whereas the successive threshold choices were K=80, K=16 and K=3. It was set out to discover patterns that involved at least 5 non-wild cards in any rolling window that spans 10 positions and begins/ends with a literal, which is a relatively high-degree of local similarity (i e., 50% or higher). Those patterns whose estimated log-probability was equal to −30.0 or less were selected and this generated a composite descriptor that comprised 1,703 patterns.
First, the corresponding DFA (deterministic finite automaton, which will only recognize instances of the composite descriptor patterns in a query sequence and which performs method step 260) was used to search a randomized version RAND-Swiss-Prot of Swiss-Prot (Release 38.0) (see also relevant discussion in the PS50040 example). The sequence regions with non-zero local counters were identified and the maximum counter values from each such legion were summed up; the sum-total was attached to the sequence label and the sequences were sorted in order of decreasing sum value. A total of 1,564 sequence fragments from RAND-Swiss-Prot received non-zero support and the actual histogram of these values is shown in
Subsequently, the same DFA was used to search the actual Swiss-Prot database testing each of its 80,236 sequences for membership in the G protein-coupled receptor family. Sum values were attached to each sequence as above and only 947 sequences from Swiss-Prot that received support greater than or equal to Thresrand=10 were kept.
In order to determine the quality of the composite descriptor and determine the number of true and false positives that the descriptor gives rise to, the Swiss-Prot annotation (keyword “KW” lines) was used for each of these 947 sequences. Of these retrieved sequences, 928 are actually listed as ‘G protein-coupled receptor's, 10 are eukaryotic transmembrane proteins (SUR7_YEAST, C561_HUMAN, YIPC_YEAST, NU4M_APIME, SCG2_XENLA, GTR2_LEIDO, GARP_HUMAN, CIN6_HUMAN, CIN3_RAT, PLSC_COCNU), 2 are hypothetical eukaryotic transmembrane proteins (YJZ3_YEAST,YMJC_CAEEL), 2 are hypothetical proteins (YKY4_YEAST, YCX7_YEAST), and finally 5 are bacterial false positives (PIP_BACCO, VIRR_AGRT6, YQGP_BACSU, HBD_CLOTS PROA_HAEIN).
This is a very notable result, given the comparatively small amount of information that is captured by the 80-sequence input set and the diversity of the G protein-coupled receptor family. Table 6 below contains a listing of the labels of the 947 Swiss-Prot sequences whose support exceeded threshold; the labels are listed in order of decreasing value of the global counter that was associated with the corresponding sequence, and the 5 false positives are shown in boldface Table 6 shows the labels of the 947 sequences from Swiss-Prot Release 38.0 that received support above threshold in the G protein-couple receptor example. The 5 false positives are shown with an “(FP).”
The third example that showcases the present invention corresponds to the helix-turn-helix motif that mediates the binding of many regulatory proteins to regulatory control sites of DNA. This 20 amino-acid long structural motif consists of two helices (7 and 9 a.a. respectively) that are separated by a 4 amino acid turn that are held together through non-polar interactions of their side chains. It has been argued that sequence-based analysis using traditional approaches cannot unambiguously identify helix-turn-helix motifs unless it is combined with the use of stereo-chemical constraints. More recently, a pattern-based approach started with 91 carefully-selected, aligned sequence fragments that corresponded to known helix-turn-helix instances and produced significant results by essentially estimating a pattern-based profile for the helix-turn-helix binding motif. This set of 91 fragments is particularly interesting because it is a very diverse collection of helix-turn-helix motif instances that share very little at the sequence level.
In the experiment carried out, a subset of 70 fragments from the set of 91 were selected (excluding those of the helix-turn-helix instances that corresponded to pieces of homeoboxes) and no alignment information was assumed. Additionally, each of the fragments was extended to the left and to the right by including an additional 10 amino acids, thus producing fragments that were 40 amino acids long. Again, the patterns were discovered assuming the equivalence classes {A, G}, C, {D, E}, {F, Y}, H, {I, L, M, V}, {K, R}, {N, Q}, P, {S, 1}, W. The Teiresias parameters were set to L=5, W=10 whereas the successive threshold choices were K=70/5-14, K=3 and K=2. It was set out to discover patterns that involved at least 5 non-wild cards in any rolling window spanning 10 positions that begins/ends with a literal, a relatively high-degree of local similarity (i.e. 50% or higher). From the discovered set, those patterns whose estimated log-probability was equal to −30 0 or less were selected, thus giving rise to a composite descriptor with 517 patterns. Table 7 below lists the labels of the 70 fragments in this training set. Table 7 shows Swiss-Prot labels of the 70 sequence fragments with length 40 a.a. in the training set for the helix-turn-helix experiment.
The resulting DFA (deterministic finite automaton, which will only recognize instances of the composite descriptor patterns in a query sequence and which performs method step 260) was used to search the randomized version RAND-Swiss-Prot of Swiss-Prot (Release 38.0) and therein were discovered a total of 277 randomized sequences that received non-zero support. Of the 2.77 randomized sequences, 275 received a support value that was less than or equal to 6. Thus, Thresrand was set equal to 7. This threshold choice corresponded to the 99.2-th percentile.
Subsequent search of the actual Swiss-Prot database gave rise to 193 sequences that received support greater than or equal to Thresrand=−7. The support values ranged from the minimum allowed value of 7 to a maximum value of 66.
Next, the Swiss-Prot annotation (feature table “FT” lines and description “DE” lines) was used for each of these 193 sequences. Of these, 169 are actually listed in Swiss-Prot as containing a helix-turn-helix motif, 2 are listed as belonging to an H-T-H group from PFAM (Y4WC_RHISN, Y4AM_RHISN) and 3 are listed as having dna-binding properties (VR2B_BPT4) or being putative DNA replication proteins (Y4CK_RHISN) or being a cytosine-specific methyltransferase (MIE8_ECOLI). Of the remaining proteins, 1 is listed as hypothetical protein (YP60_METTM), 1 is listed as a hypothetical transcription factor containing a helix-turn-helix motif (Y558_METJA), 1 is listed as being involved in DNA packaging (XIMA_BACSU), 1 is listed as having strong similarity to MJ1545 which is a putative transcription repressor protein containing a helix-turn-helix motif (YO14_ARCFU), 3 have very good blastp P-values with all the similarities confined in the helix-turn-helix region of the input fragments (PRPD_SALTY, PRPD_ECOLI, Y0FO_MYCTU), and finally, 2 are likely to be false positives (YOAE_ECOLI, CIPE_MYCTU). Table 8 below contains a listing of the labels of these 193 hits in order of decreasing value of accumulated support. Table 8 shows the Swiss-Prot labels of the 193 sequence fragments that are discovered using the composite descriptor derived from the original set of 70 fragments.
Starting now with the set of all 193 discovered sequence fragments, one more iteration of the described method was carried out using this set as the new training set, T. The training set for this iteration was formed by collecting the individual sequence fragments whose support exceeded threshold. As before, the Teiresias parameters were set to L=5 and W=10 whereas the successive threshold choices were K=193/5=38, K=7 and K=2. Sub-selecting those patterns whose estimated log-probability was equal to −30.0 or less produced 1,061 patterns which were added to the previous set of 517 to form a new augmented composite descriptor. The DFA resulting from the latter descriptor was applied to RAND-Swiss-Prot. Of the 537 sequence fragments that received non-zero support, 534 received support 9 or less thus establishing the value 10 as the new Thresrand (=99.2-th percentile). Processing Swiss-Prot with this last DFA, an additional 96 sequence fragments were discovered that exceeded threshold for a grand total of 289 fragments. Table 9 here lists the labels for this additional set of fragments. Table 9 shows the Swiss-Prot labels of the additional 96 sequence fragments that are discovered after augmenting the original composite descriptor with the patterns that are discovered from treating the first set of 193 discovered fragments as a training set.
An analysis of the additional hits using the feature tables in Swiss-Prot showed that 81 of those are true positives, 4 are listed as DNA binding (TRA4_BACFR,V57A_BPT4,NADR_KLEPN) or transcription regulation proteins (EBSC_ENTFA), and 2 are listed as hypothetical proteins (YFED_ECOLI, YD16_SCHPO). Finally, 8 hits probably correspond to false positives (Y4BA_RHISN, EMRD_ECOLI, VG36_BPMD2, VG36_BPML5, TETP_CLOPE, YSCC_YEREN, MERB_SERMA, MLID_ECOLI).
A Fourth Example Searching the C. elegans Genome for EFIG, GPCR and HTH CandidatesThe three composite descriptors were used to search the collection of 19,099 ORFs that were reported for the C. elegans genome, by the Washington University in St. Louis, School of Medicine, Genome Sequence Center, as of Jun. 13, 1999. In all three cases, the corresponding values of Thresrand that were established by searching RAND-Swiss-Prot were used.
Elongation Factor 1 Gamma Chain
First, this ORF collection was searched using the 2,260 pattern composite descriptor that was built for the elongation factor gamma chain (PS50040 above). Of the 13 ORFs that received non-zero support only one, F17C 11.9, exceeded threshold. This ORF is the one listed in Swiss-Prot (and in PS50040) as EF1G_CAEEL.
G-Protein Coupled Receptor's
Next, the C. elegans genome was searched using the composite descriptor for the G protein-coupled receptor that comprised 1,703 patterns. Note that for this particular experiment, it was not set out to discover and enumerate all putative G-protein coupled receptors in C. elegans but rather to show that even when starting with a small knowledge base that contains no GPCR sequences from the genome under consideration it can be effective to mine a complete genome such as C. elegans.
In Table 10 below, the labels of the 101 C. elegans ORFs whose support exceeded threshold are shown. For each of those ORs, the Score and the P and N values are shown for the top scoring sequence obtained from running a BLASIP search against the set of 804 Swiss-Prot Rel. 35.0 sequences that are known to be true GPCRs (see also discussion above). Table 10 shows the 101 ORFs from C. elegans that were discovered using a composite descriptor for the GPCR family and whose support exceeds threshold. For each of the reported ORFs, also listed are the top scoring sequence from running blastp against the set of 804 Swiss-Prot Rel 35 sequences that are known to be true GPCRs.
In addition to the above 101 C. elegans ORF s that exceeded threshold and as testimony to the stringent thresholds use, there is also listed in Table 11 below an additional 19 ORFs whose scores were just below threshold and which generated blast-search P values that were significant. As before, the blast searches were carried out against the set of 804 Swiss-Prot Rel. 35.0 known true GPCRs. Table 11 shows an additional 19 ORFs from C. elegans that receive scores just below threshold but show significant blast-search P values when compared against the set of 804 true GPCRs from Rel 35.0 of Swiss-Prot.
Several comments are in order here. First, it should be stressed that the above analysis is not implying that there is only 120 G-protein coupled receptors in C. elegans. Instead, what is attempted to be demonstrated is that even if one begins with a small knowledge base of only 80 known GPCRs that have been selected randomly, one can still build a pretty useful composite descriptor for the family and use it to explore a largely-unexpected genome such as C. elegans. In order to have a complete enumeration of the GPCRs that are present in C. elegans, the composite descriptor should be built by using all of the GPCRs that are present in GPCRDB and not only 80 of them. Second, it was opted the BLAST searches against the set of 804 sequences in order to show the ability of the proposed method to extrapolate. As such, blast-search results with P values that are relatively high (e.g. E-02) should not be surprising since the target database of 804 true GPCRs is but a small fraction of the current contents of GPCRDB. Indeed the November 1999 release of GPCRDB contained 1,704 GPCR sequences and 431 GPCR sequence fragments for a grand total of 2,135 entries
Helix-Turn-Helix
Finally, the 19,099 ORFs of C. elegans was searched for instances of the helix-turn-helix binding motif using the corresponding 2,288 (=1,896+392) pattern composite descriptor. Of the 169 sequences that received non-zero support, only 5 exceeded threshold : Y94H6A—142.g (in the region delineated by a.a. 65 through 95), C16C2.1 (in the region delineated by a.a 59 through 89), F18C5.2 (in the region delineated by a.a. 850 through 880), Y39F10A a (in the region delineated by a.a 125 through 155), Y48C3A.s (in the region delineated by a a. 113 through 143), and Y48C3A.s (in the region delineated by a.a 113 through 143),
The fragments were:
and their respective top-scoring sequences from the training set of 70 helix-turn-helix segments, blast scores, P and N values are:
It is to be understood that the embodiments and variations shown and described herein are merely illustrative of the principles of this invention and that various modifications may be implemented by those skilled in the art without departing from the scope and spirit of the invention.
Claims
1. A method for unsupervised building and exploitation of composite descriptors, the method comprising the steps of;
- i. providing a training set of sequences, each sequence comprising a plurality of symbols;
- ii determining a set of maximal patterns, each of the maximal patterns being common to a predetermined number of the sequences, wherein the step of determining a set of maximal patterns is performed without any knowledge about properties or features of sequences in the set of unaligned sequences;
- iii determining which, if any, of the maximal patterns are statistically significant; and
- iv creating a composite descriptor from the statistically significant maximal patterns.
2. The method of claim 1, wherein the sequences in the training set are unaligned.
3. The method of claim 1, wherein the step of creating a composite description from the statistically significant maximal patterns further comprises the steps determining which of the statistically significant maximal patterns are currently not part of the composite descriptor, adding those statistically significant maximal patterns that are currently not part of the composite descriptor to the composite descriptor, and removing the added statistically significant maximal patterns from the training set of sequences.
4. The method of claim 3, wherein each symbol comes from an alphabet that describes DNA (deoxyribonucleic acid) or proteins.
5. The method of claim 1, wherein the symbols are numerical.
6. The method of claim 3, further comprising the steps of iterating steps (ii) through (iv) until either the training set contains no sequences or there are no statistically significant maximal patterns common to the sequences in the training set.
7. The method of claim 3, further comprises the step of determining if a candidate sequence comprises a predetermined number of the statistically significant maximal patterns.
8. The method of claim 7, comprising the steps of if the candidate sequence comprises the predetermined number of the statistically significant maximal patterns, adding the candidate sequence to the set of sequences to create a new training set of sequences and performing the steps (ii) through (iv) on the new training set of sequences.
9. The method of claim 1, wherein the step of determining which, if any, of the maximal patterns ate statistically significant comprises the step of determining for each of the maximal patterns if a probability that this maximal pattern occurs in a sequence meets a predetermined threshold.
10. The method of claim 1, wherein the set of maximal patterns is empty and wherein the step of determining a set of maximal patterns further comprises the steps of reducing the predetermined number of sequences and performing step (ii) again.
11. A system for unsupervised building and exploitation of composite descriptors, comprising:
- a memory that stores computer-readable code; and
- a processor operatively coupled to said memory, said processor configured to implement said computer-readable code, said computer-readable code configured to; i. provide a training set of sequences, each sequence comprising a plurality of alphabetic symbols; ii. determine a set of maximal patterns, each of the maximal patterns being common to a predetermined number of the sequences, wherein the maximal patterns are determined without any knowledge about properties or features of sequences in the set of unaligned sequences; iii. determine which, if any, of the maximal patterns are statistically significant; and iv. create a composite descriptor from the statistically significant maximal patterns
12. An article of manufacture for unsupervised building and exploitation of composite descriptors, comprising:
- a computer readable medium having computer readable code means embodied thereon, said computer readable program code means comprising; a step to provide a training set of sequences, each sequence comprising a plurality of alphabetic symbols; a step to determine a set of maximal patterns, each of the maximal patterns being common to a predetermined number of the sequences, wherein the maximal patterns are determined without any knowledge about properties or features of sequences in the set of unaligned sequences; a step to determine which, if any, of the maximal patterns are statistically significant; and a step to create a composite descriptor from the statistically significant maximal patterns.
Type: Application
Filed: Jun 21, 2007
Publication Date: Jan 17, 2008
Applicant: International Business Machines Corporation (Armonk, NY)
Inventors: Isidore Rigoutsos (Astoria, NY), Yuan Gao (Croton On Hudson, NY), Aristidis Floratos (Long Island City, NY)
Application Number: 11/766,256
International Classification: G01N 33/48 (20060101);
