US 5,867,402 AGrant

Computational Analysis of Nucleic Acid Information Defines Binding Sites

Issue Date:1999-02-02

•57 Claims

•22 Drawing Sheets

Abstract

In accordance with the present invention, binding sites are defined based upon the individual information content of a particular site of interest. Substitutions within the binding site sequences can be analyzed to determine whether the substitution will cause a deleterious mutation or a benign polymorphism. In addition, new binding sites can be identified using individual information content. Further a computer system is described for determining and displaying individual information content of a binding site sequence.

Metadata

Assignee

The United States of America as represented by the Department of Health and Human Services

Inventors

Thomas D. Schneider
Peter K. Rogan

Application Information

Application Number:US 4941152

Filing Date:1995-06-23

Priority Date:1995-06-23

Art Unit:244

Classifications

IPC:

G06F 1716

Field of Search:

364395434496-500;578144279

Patent Drawings (22 sheets)

Description

Field of the Invention

The present invention relates to information computational methods of defining binding sites.

Background of the Invention

When studying molecular binding sites in DNA or RNA, it is conventional practice to align the sequences of several sites recognized by the same macromolecular recognizer and then to choose the most common bases at each position to create a consensus sequence (see Davidson et al., 1983. Nature (London), 301, 468-470). Consensus sequences are difficult to work with and are not reliable when searching for new sites (Sadler et al., 1983b. Nucl. Acids Res. 11:2221-2231; 26 Hawley & McClure, (1983) Nuc. Acids Res.; 11:2237-2255).

This is partly because information is lost when the relative frequency of specific bases at each position is ignored. For example, the first position of Escherichia coli translational initiation codons has 94% Adenine ("A"), 5% Guanine ("G"), 1% Uracil ("U") and 0% Cytosine ("C"), which is not represented precisely by the consensus "A". To avoid this problem, four histograms can be made that record the frequencies of each base at each position of the aligned sequences. Such histograms can be compressed into a single curve by the use of a .chi..sup.2 function (Gold et al., 1981. Annu. Rev. Microbiol. 35, 365-403; Stormo et al., 1982. Nuc. Acids Res. 10, 2971-2996). Although these curves show where information lies in the site, they have several disadvantages: the .chi..sup.2 scale is not easily understood in simple terms; it is difficult to compare the overall information content of two different kinds of sites, such as ribosome binding sites and restriction enzyme sites; and .chi..sup.2 histograms are not directly useful in searching for new site sites (Stormo et al. 1982 Nuc. Acids Res. 10, 2997-3011).

Many general methods exists for identifying sequence changes which are deleterious. However, these methods require experimentation in the laboratory. The most common method is the identification of a disease state and a corresponding genetic mutation in a particular sequence element. This method is quite labor intensive and requires that the mutation produce an identifiable phenotype. Another method uses restriction fragment length polymorphisms to identify alterations within the genome. This method is also experimental, but can only detect alterations in the genome at restriction sites, whether or not a phenotype results.

The average information contained in a set of nucleic-acid binding sites can be calculated by using the methods of information theory, and this has been useful for understanding a number of genetic control systems (Schneider et al., 1986. J. Mol. Biol., 188, 415-431; Schneider & Stormo, 1989. Nuc. Acids. Res., 17, 659-674; Eiglmeier et al. 1989. Mol. Microb., 3, 869-878; Penotti, 1990. J. Mol. Biol., 213, 37-52; Penotti, 1991 J. Theor. Biol. 150, 385-420; Schneider & Stephens, 1990. Nuc. Acids. Res., 18, 6097-6100; Herman & Schneider, 1992. J. Bact., 174, 3558-3560; Gutell et al., 1992. Nuc. Acids. Res., 20, 5785-5795; Papp et al., 1993. J. Mol. Biol., 233, 219-230). However, thus far an effective method does not exist for working with information content of single sequences or for predicting the effect of changes in information content due to sequence alterations--be it through biological evolution or by genetic manipulation.

Information analysis of normal splice junctions reveals partially conserved nucleotide sequences that are not always reflected in the corresponding consensus sequence (Stephens & Schneider, 1992. J. Mol. Biol. 228:1124-1136). Information content may be represented by a sequence logo, which depicts the relative contribution of each position of the splice site and the relative frequencies of each nucleotide at every position (Schneider & Stephens, 1990. Nucl. Acids Res. 18:6097-6100). The logo illustrates the full range of normal variants in the splice junction.

The present invention is principally directed to binding sites on a sequence. In particular, the present "Walker" program enables a scientist or clinician to identify mutations within a nucleic acid binding site which are deleterious, without extensive experimentation. This method generates a model of the binding site which is called the R.sub.i (b,l) weight matrix, which can then be used to evaluate other individual sites for their information content. The present invention allows one to analyze the effect on the binding site of changing a base at a particular position within the site.

The weight matrices of the present invention are not found in the prior art in several respects. R.sub.i values, which represent the sum of all weights at each position within a site, are on an absolute scale, rather than the relative scale found in the prior art. R.sub.i =0 is a cutoff point for functional sites within the present invention. This feature is lacking in both Staden's method (1984 Nuc. Acids Res., 12:505-519) and Berg & von Hippel's method (1987 J. Mol. Biol., 193:723-750; 1988 J. Mol. Biol., 200:709-723; 1988 Nuc. Acids Res., 16(11):5089-5105). Hence, these methods draw no distinction between polymorphisms and mutations.

Moreover, the Berg & von Hippel's method relies upon the consensus sequence as the ideal, i.e. the best binding sequence. Therefore, Berg & von Hippel had no way of distinguishing a polymorphism from a deleterious mutation.

In addition, unlike the prior art (Berg & von Hippel's statistical-mechanical theory, in particular), no assumption about the relationship between energy and information is required to obtain R.sub.i in the present invention. The statistical-mechanical approach assumes that the energy of binding, "discrimination energy", is equal to the information contained within a recognition sequence. This assumption does not allow for a situation where more than one protein could bind to a particular site and thus increase the apparent information contained within that site.

Further, the R.sub.i method described in the present invention is much more sensitive to sequence changes than the widely and almost universally used consensus sequence method. The consensus sequence destroys data by taking the most frequent base at every position as the base used in the consensus model, whereas the R.sub.i method does not alter the frequency data and so can be used to detect subtle effects.

One object of the present invention relates to the use of individual information content of the site and its comparison with the overall distribution of individual information in a set of binding sites, to determine whether a substitution is a polymorphism or a mutation.

Another object of the present invention relates to designing binding sites to adjust the activity of the site. The present invention further relates to a computer system capable of determining the individual information content of a binding sequence and identifying new binding sequences.

Yet another object of the present invention relates to the use of individual information content to determine the effect of a particular position change in a sequence acting as a binding site.

Another object of the invention is to use the "Ri" and "Walker" computer program to display the reaction of a binding macromolecule at every position in a sequence and to determine the change in information content when a particular position within a binding site is altered.

Objects and advantages of the invention set forth herein and will also be readily appreciated here from, or may be learned by practice with the invention. These objects and advantages are realized and obtained by means of instrumentalities and combinations pointed out in the specification and claims.

Summary of the Invention

The present invention relates to identifying mutations and polymorphisms within a nucleic acid region acting as a macromolecule binding site. The invention further relates to analyzing protein regions acting as binding sites for macromolecules to identify mutations and polymorphisms within the site. In either case, the instant method relates to the identification of mutations/alterations in a sequence, either nucleic acid or amino acid, which will be deleterious to the system which it affects.

In accordance with the present invention, a computer system and computation method are described for processing sequence signals by a transformation into an information content weight matrix, as represented by R.sub.i (b,l). A second transformation follows which applies a particular sequence signal to the information content weight matrix, R.sub.i (b,l) thereby producing a value, Ri, which comprises the individual information content of said particular sequence signal. An alteration of a particular position within a binding sequence provides a third signal, transforming the individual information content of the binding sequence by the amount of information either lost or gained by the position change. The third transformation produces an output record, for example a graphical representation (an X-Y graph or a numerical value) of the information content of the sequence after the alteration and defines whether the alteration will be deleterious to the cell. Such a deleterious alteration is referred to as a mutation, whereas a non-deleterious alteration is a polymorphism. The invention also relates to computer programs embodied on a computer-readable medium.

The present invention also relates to the display of the product of the transformations of the present method in the form of a graphical image.

The present invention further relates to a method for identifying and manipulating the binding affinity of a particular position within and surrounding a binding site. The instant method allows comparison of the information on particular binding sites to the individual information content of other binding sites, to distances between features of the sequence, and to their measured binding energies. The present invention further allows adjustment of the binding affinity of a binding site by manipulating positions within the site to alter its individual information content.

The present invention further relates to a method of designing sequence elements which function as binding sites.

The invention also relates to a method of diagnosing a genetically-determined disease based upon the identification of a deleterious mutation, based upon a change in individual information content of the binding sequence. In addition, the invention relates to identification and use of cryptic binding sites on a particular sequence.

Brief Descriptions of the Drawings

The file of this patent contains at least one drawing executed in color. Copies of this patent with color drawing(s) will be provided by the Patent and Trademark Office upon request and payment of the necessary fee.

FIG. 1(a) is the sequence 5' CAGGTCTGCA 3' represented in matrix format.

FIG. 1(b) is the individual weight matrix for human donor splice junctions derived from data given in (Stephens & Schneider, 1992 J. Mol. Biol., 228: 1124-1136). The weights of the matrix in (b) which are selected by the sequence in (a) are enclosed by boxes.

FIG. 2 is a histogram of individual information for 1055 E. coli ribosome binding sites. The mean and standard deviation of the R.sub.i values were fitted by a Gaussian distribution.

FIG. 3 is a histogram of individual information for 1799 human donor binding sites. Donor sites which lacked a complete sequence in the region 0 to +6 were not included.

FIG. 4 is a histogram of individual information for 1744 human acceptor binding sites. Acceptor sites which lacked a complete sequence in the region -25 to +2 were not included.

FIG. 5 is a graph illustrating the correlation between GCN4 sites and the log of their relative binding affinities.

FIG. 6 is a schematic diagram showing important landmarks on the individual information, R.sub.i scale. The consensus is the highest possible evaluation of the R.sub.i (b,l) matrix; the anti-consensus is the lowest. Sequences with R.sub.i =0 separate sites (R.sub.i >0) from nonsites (R.sub.i <0). By definition, the mean of the distribution is R.sub.sequence. The standard deviation of the distribution is .sigma..sub.R.sbsb.i. The standard deviation of R.sub.sequence is the standard error of the mean, SEM.

FIG. 7 is a flow diagram illustrating computer programs for individual information analysis in accordance with the present invention.

FIG. 8 is a graphic plot of the individual information of the Fis Promoter produced by the program Xyplo. The position of the zero base of the Fis weight matrix on the sequence is given on the abscissa, while the individual information for the sequence surrounding each position -10 to +10 is given on the ordinate. The 6 previously identified Fis sites are marked with a plus (+). Predicted sites are represented as squares above the zero line. Transcription begins at base 375 and proceeds to the right (arrow). The sequence is from GenBank accession X62399 (Ninnemann et al., 1992 EMBO J., 11:1075-1083) (see also accession M95784 (Ball et al., 1992 J. Bact., 174:8043-8056). Bottom: a larger region of sequence graphed by DNA plot shows clustering of potential fis sites around the promoter but not further downstream. The dashed line indicates the corresponding parts of the figure.

FIG. 9 is an example of a Walk Display.

FIG. 10 is an example showing the effect of mutations in a "Walk" display.

FIG. 11 is a sequence logo showing the location of the hMSH2 polymorphism in the human splice acceptor site. This sequence logo was created from 1744 wild-type acceptor sites. The height of each nucleotide is proportional to its frequency at that position, while the height of each entire stack of nucleotides corresponds to the information measure (in bits) or, equivalently, the sequence conservation at that position. When sequence conservation is measured in bits, the relative heights of the stacks can be compared to one another and the total sequence conservation in a region can be found by adding the heights of the stacks together (Shannon & Weaver, 1949, The Mathematical Theory of Communication, University of Illinois Press, Urbana, Ill.). Coordinates in the splice site are defined along the abscissa. RNA strand cleavage during splicing occurs at the vertical line between positions 0 and 1. All positions except -3 in this logo are significantly above background (p<8.times.10.sup.-8). The arrow shows the position of the T.fwdarw.C substitution of the hMSH2 gene.

FIG. 12 is a set of graphs illustrating individual information scans of inversion regions. Symbols are the same as in FIG. 5. Previously identified Fis sites are marked with a plus inside a square and named as in (Finkel & Johnson, 1992 Molec. Microb., 6:3257-3265; Finkel & Johnson, 1992 Molec. Microb., 6:1023). The proposed Fis sites are marked with a circle inside a square. Spacing between sites is indicated by numerals surrounded by dashes. Note that the spacing between proximal and distal sites is always 48 bases.

FIG. 13 is the sequence for the S. typhimurium hin mutants. The wild-type sequence containing the proximal Fis site from the S. typhimurium hin region (HW) is given on the top, flanked by EcoRI and HindIII restriction sites. The known proximal site is indicated next to the predicted medial site. In the next sequence, the right anticonsensus (HR) was used to destroy the medial site, leaving the proximal site intact. In the third sequence, the left anticonsensus (HL) Fis site sequence was used to destroy the proximal site while leaving the medial site intact. In the fourth sequence both (HB) sites were destroyed.

FIG. 14 is a matrix table for the n(b,l) and the R.sub.i (b,l) weight matrix for 76 Fis binding sites. Column l is the position relative to the center of the Fis site. Columns n(a,l), n(c,l), n(g,l) and n(t,l) give the number of bases b at positions l (the n(b,l) table). The frequency table is defined as f(b,l)=n(b,l)/.SIGMA..sub.b.sup.T =.sub.A n(b,l). The 4 columns for the R.sub.i (b,l) table give the individual information weights (in bits) for bases b at position l. This distribution of Fis sites has a mean of 8.24 bits and a standard deviation of 2.69 bits.

FIG. 15 is a sequence logo of Fis binding sites and DNA base pair structure with 38 experimentally defined Fis binding sequences and their complements. The total sequence conservation, found by adding the stack heights together, is R.sub.sequence =8.2.+-.0.6 bits per site. (this standard error of the mean=0.6 bits was calculated according to (Schneider et al., 1986 J. Mol. Biol. 188:415-431). See text for further description. Methylated guanines which interfere with Fis binding are indicated by filled circles (.circle-solid.) and methylated adenines which interfere with Fis binding are indicated by open circles (.smallcircle.) (Bruist et al., 1987 Jones Dev., 1:762-772)

FIG. 16 are mobility shift experiments for hin and cin. Top: Gel shifts of DNA containing the hin proximal and medial Fis binding sites. Each lane contains increasing concentrations of Fis protein added, beginning with no Fis protein, Fis diluted 1 to 8, etc. The 1:1 ratio is 1000 nM Fis. Letter designations refer to the sequences given in FIG. 13. Bottom: Gel shifts of DNA containing the cin proximal and external Fis binding sites with the same conditions as above.

FIG. 17 is a Scattergram showing the relationship between R.sub.i and phenotype of mutations altering splice donor sequences. The clinical presentations of each inherited abnormality studied were categorized as mild, moderate or severely affected based on the descriptions of these patients. On the ordinate axis, individuals with a mild disorder are coded as 1, moderate as 2, and severely affected individuals as 3. The individual information content of the corresponding mutations is plotted on the abscissa.

FIG. 18 is a graph showing the relationship between mutant R.sub.i and splicing efficiency for mutant donor splice sites. The relationship between the logarithm (base 2) of the mRNA splicing efficiency with the change in R.sub.i for 21 splice donor site mutations (I) or 10 acceptor site mutations (II) associated with different inherited conditions. The change in R.sub.i due to the mutation (.DELTA.ri) is expressed as the normalized absolute value of the difference between the mutant and cognate individual information content values (in bits). According to this definition, a non-functional mutation will have .DELTA.R.sub.i =1; .DELTA.R.sub.i for a polymorphic substitution will be zero. The logarithm of the splice efficiency ranges from 10 (for 100% efficiency) to -27 (for negligible levels of splicing, this has been set at 1.times.10.sup.-8 %, since the logarithm of 0 cannot be computed). Regression of the best linear fit of the data is shown as a line. the correlation coefficients of .DELTA.R.sub.i of donor and acceptor splice site mutations, respectively, are 0.45 and 0.68.

FIG. 19 is a graph showing the relationship between mutant R.sub.i and splicing efficiency for mutant acceptor splice sites. (See FIG. 18 above for details.)

FIG. 20 is a Scan plot of nrd binding sites.

FIG. 21 is a graph showing measured splice product for variations in the polypyrimidine tract of the adenovirus 2 intron of the major late promoter Leader 1 and Leader 2 splicing unit versus individual information R.sub.i of the same sequences.

Detailed Description of the Invention

The present invention relates to a method of identifying and manipulating the affinity of macromolecule binding sites.

The present invention further provides a method for identifying mutations and polymorphisms within a nucleic acid region acting as a protein or--macromolecule--binding site. The method is also used on other binding sites, such as protein--protein sites or protein binding sites for other small molecules. In particular, these sites are analyzed to determine whether a particular amino acid substitution is deleterious or not. In either case, the method relates to the identification of alterations in a sequence, either nucleic acid or amino acid, which will be deleterious to the system.

Information as the term is used herein is defined as the number of choices made by a machine, given on a logarithmic scale in bits. Information content as the term is used herein is defined as the number of choices needed to describe a sequence pattern given on a logarithmic scale in bits. R.sub.sequence as the term used herein is defined as the information content of a nucleic-acid binding site or of a protein. A "binding site" as the term is used herein is defined as the region of macromolecule which binds to another molecule. A "cryptic site" as the term is used herein is defined as a weak nucleic-acid binding site that may be revealed by mutation of the sequence or by the destruction of a neighboring strong site. Splicing efficiency is defined as the proportion of normal mRNA produced by the mutant allele relative to the normal allele.

The individual information theory methods of the present invention can be applied to genetic engineering.

For example, polymorphisms are extremely useful tools in genetic engineering, one can use individual information analysis to introduce a polymorphism into splice sites or other types of motifs. For splice sites, one might introduce a substitution that does not impair splicing (at both the donor and acceptor sites flanking an exon) and which produces cleavable restriction sites at either end of the exon. This permits the investigator to "shuffle" the exon(s) in vitro to create a novel protein with additional functions. The only prerequisite for such an application of the instant invention is that the reading frame is preserved. A real benefit to this embodiment of the present invention is that it eliminates the necessity of flanking intron sequences to be carried along with the desired intron sequence. This is an important consideration, since not all introns consist of "junk" DNA. For example there are transcription factor binding sites, i.e. Oct 1 in the immunoglobulin V region introns, internal promoters, i.e. in the murine major histocompatibility complex genes, and even more important, cryptic splices sites in the flanking intron sequences, as shown below.

The same strategy is used to move promoters from one gene to another.

There is another way human geneticists may use the instant invention. A polymorphism may be introduced, without a loss of information or change in function, in order to track a transgene or a transfected gene in a cell type where other similar sequences may be present. An example of such a system is when the introduced gene is a member of a multigene family. Using Ri analysis to insure that the polymorphism does not have an effect on splicing or other aspects of gene or protein expression is an important consideration. Introducing transgenes in this way permits distinguishing maternally and paternally derived chromosomes, thus providing another tool for identification of imprinted genes.

Another embodiment of the present invention utilizes individual information techniques to allow design of binding sites. As genetic engineering advances, it is useful to have the capability to create more complex genetic structures. The strongest binding sites are not always the most desirable. For example, since a strong bacteriophage T7 promoter will kill a bacterial cell or tax the resources of the cell by using up the free ribonucleotides, it is at times not practical to have the strongest possible promoter. The tools of the present invention allow the design a promoter of the strength required for a particular application. In the case of T7 promoters, one may find an optimum at which a promoter strength is chosen which maximizes production of a gene product because cells are still healthy. These same tools allow not only the creation of "designer promoters" and "designer genetic control systems", but also the design of the active site in an enzyme, other motifs in proteins and drug binding sites.

Sometimes it is necessary to design two binding sites that overlap each other. To insure that each site has the required strength is impossible with the consensus methods, but the individual information technique of the present invention easily allows this. The present invention allows a user to select from many weight matrices which may be stored in a library. A computer may automatically evaluate the effect on binding for each recognizer and of any changes to the sequence that the user contemplates. This allows the user to modulate the strengths of the binding sites individually so that these binding sites work together for the desired genetic effect. This embodiment allows the fine-tuning of gene expression.

The term "recognizer" as it is used herein refers to a molecule which recognizes and binds to the binding site of interest. Recognizer is further defined to mean a macromolecule that locates specific sites on nucleic acids. These may include repressors, activators, polymerases, ribosomes and spliceosome.

The method uses individual information content to determine the effect of a particular mutation at a specific binding site. The information content for a particular binding site is derived from an analysis of nucleic acid sequence information from various data bases available for that information, which are used to determine the frequency of a particular nucleic acid base being present at a particular position within the sequence of interest. This analysis results in the development of a sequence logo, which is a graphical representation of the probability that a particular nucleic acid base will be present at a particular position within the sequence of interest. The height of each nucleotide within the sequence logo is proportional to its frequency at that position, while the height of each stack of nucleotides corresponds to the information measure (in bits) or, equivalently, the sequence conservation at that position. The area under the logo represents the information content in bits (referred to as R.sub.sequence) of the binding site. The logo illustrates the full range of normal variants in the protein binding site of interest.

Preliminarily, a sequence logo may be produced. Each binding site, as represented in a sequence logo contains a specific amount of information, which is expressed as R.sub.sequence, in "bits" of information. A bit, as the term is used herein, is the amount of information needed to choose one of two equally likely possible outcomes.

In accordance with the present method the gathered sequence information regarding a binding site is converted into a weight matrix, referred to as R.sub.i (b,l) which provides a model of the recognizer which binds to the binding sites. The information weight matrix is then applied to a particular sequence, generating an individual information content, Ri for that sequence. This sequence can be further analyzed for the effect of a specific mutation at any position within the sequence and the resulting change in Ri can be measured. A nucleic acid substitution may be analyzed for a change in the individual information content, which can be displayed in the sequence logo-type image (see FIGS. 9 and 10). True mutations are expected to reside in positions where the sequence conservation in bits significantly exceeds the background variation and where the base frequency decreases significantly.

It is noted that R.sub.i (b,l) is related to R.sub.sequence, in that R.sub.sequence is the mean value generated from the R.sub.i (b,l) matrix, when that matrix is applied to the original site of binding sequences used to create the R.sub.i (b,l) matrix itself.

Practically, the invention is effectuated through the use a series of computer programs, which sequentially, retrieve selected nucleic acid information,and analyze the information content for the sequences retrieved by development of a weight matrix (in the Ri program). The weight matrix is applied to a specific sequence thereby producing an individual information content, Ri which is then loaded into a program called "Walker", which is capable of displaying the reaction of the binding protein to every base in a sequence and of determining the effect of a nucleic acid substitution at any position within the binding site, based on the information weight matrix. Where there is a change in the individual information content which deviates from the determined information content value by more than three standard deviations, or which makes the individual information content go below zero, these changes are considered mutations rather than polymorphisms. The Ri program, the Walker program and the related Scan program allow the user to investigate the effects of sequence changes in the regions around the binding site, so that the creation or destruction of binding sites nearby can be detected.

A preferred embodiment of the present invention relates to a method for assigning a sequence conservation to individual nucleic-acid binding site sequences based on a large collection of sample sites. In this method, the sample sequences bound by a particular protein or molecular complex (such as a ribosome or spliceosome) are aligned and the frequencies of bases at each position are determined. The base 2 logarithm of each frequency at every position is added to 2 and a sample size correction factor to obtain a weight matrix, R.sub.i (b,l), where b is one of the 4 bases and l is a position along the sequences. This "individual information" matrix represents the sequence conservation of the sites measured in bits of information and it can be used to rank-order the sites, to search for new sites, to compare binding sites of the same or of different kinds to one another, to compare binding sites to other quantitative data such as binding energy or distance between binding sites, and to detect errors in databases.

In accordance with the present invention the individual information matrix is: ##EQU1## where f(b,l) is the frequency of each base b at position l in the aligned binding site sequences and e(n(l)) is a sample size correction factor for the n sequences used to create f(b,l) (Schneider et al., 1986 J. Mol. Biol., 188:415-31; Penotti, 1990 J. Mol. Biol., 213:37-52). To simplify the notation, the factor e(n(l)) was separated from log.sub.2 f(b,l) and joined to "2" to create E(H.sub.n).

In a set of sequences the j.sup.th sequence by a matrix s(b,l,j) contains only 0's and 1's. For example, the sequence 5' CAGGTCTGCA 3' is represented as shown in FIG. 1a. Likewise, the R.sub.i (b,l) matrix for human donor splice junctions is shown in FIG. 1b.

The individual information of a sequence is the dot product between the sequence and the weight matrix: ##EQU2## For the donor splicing weight matrix given in the figure, the sequence 5' CAGGTCTGCA 3' is assigned 0.58+1.25+1.64+1.99+1.98+(-3.68)+(-1.59)+1.71+(-0.51)+0.05=3.42 bits. Essentially, each base of the sequence "picks out" a particular entry from a column of the R.sub.i (b,l) matrix, and these weights are added together to produce the total R.sub.i.

The average information of the n individual sequences which were used to create the frequency matrix f(b,l) is the expectation (i.e. mean) of R.sub.i : ##EQU3## We now substitute equation (1) into (2) and then substitute equation (2) into (3). By using the definition of the frequency matrix: ##EQU4## and the fact that the frequencies sum to 1: ##EQU5## we find with some manipulation that: ##EQU6## The right hand side is exactly the definition of R.sub.sequence (Schneider et al., 1986 J. Mol. Biol., 188:415-431), so we have demonstrated that the average of individual information contents is the average information content of the sites:

By expressing the formula (6) as a subtraction, we emphasize that information is a state function defined as a difference of uncertainties (Schneider, 1994 Nanotechnology, 5(1):1-18). The R.sub.i (b,l) function is unique because it can be proven that R.sub.i (b,l) is the only function whose average is R.sub.sequence, as described above.

Roots of information theory: surprisal of bases

The individual information method is consistent with early work on information theory. Selecting one symbol from a set of M symbols, requires log.sub.2 M binary decisions. For example, any corner of a cube may be specified by the answer to three yes-no questions of the form:

1. Is it on top?

2. Is it on the left side?

3. Is it in front?

That is, log.sub.2 8=3 bits. The next step is to rearrange the formula:

where P is the probability of the equally likely symbols. What can we do if the symbols are not equally likely, as is the case for frequencies of bases in binding sites? To handle this, Tribus (Tribus, 1961 Thermostatics and Thermodynamics, D. van Nostrand G., Inc., Princeton, N.J.) proposed the concept of "surprisal", h.sub.i as the negative logarithm of a symbol's probability in the midst of a stream of symbols:

where p.sub.i is the i.sup.th symbol's probability so that (9) is an extension of the form given in equations (8). For example, the less likely the ringing of a telephone is, the more startled we are to hear it. The advantage of using this definition becomes clear when we consider the average surprisal for the entire stream of symbols. To find this, we take the individual surprisals and weight them by their occurrence, p.sub.i, and find the total: ##EQU7## This is the Shannon uncertainty measure, so we have demonstrated that H is an average of surprisals.

What change does an individual "finger" of a recognizer see when the recognizer goes from non-specific binding (the before state) to specific binding (the after state)? In the before state, the average surprisal may be 2 bits, since the recognizer is not making contact with the nucleic acid bases in that state, so the composition of the genome should not matter. It is noted that 2 in equation (1) represents 2 bits of information, that is the uncertainty before a recognizer binds to a binding site. However, it may alternatively be represented by a value H.sub.g which represents the uncertainty associated with binding anywhere in a particular genome. This value will vary from one genome to the next but will be a constant for all binding sites within one genome. Thus for the difference in surprisal we write:

This is equation (1) except for the sampling correction.

We are now in a position to understand the individual information, which is the sum of R.sub.i (b,l) across a binding site, as the total surprisal decrease from the viewpoint of a particular recognizer binding to a particular sequence. This model allows a recognizer to have different responses to different sequences. Different recognizers have different surprisals for the same sequence because they have different molecular recognition surfaces.

A word of caution is in order. If the set of sequences contains gaps (as when sequence data are missing on one or both sides of a site) then the average of the individual information contents generally will not equal the R.sub.sequence as calculated from the frequencies of bases at each position. This is because the individual sequences can be strongly affected by missing data, but R.sub.sequence is not. For this reason calculation of R.sub.sequence should still be done by the original frequencies method, and individual information values taken from partial sequence data should be treated with care.

The model described above assumes that positions along the site are independent from one another. It should be possible to extend the method to cases where each base is correlated to the next, or even longer relationships. However, to do this requires many more sequences to avoid the severe effects of small sample size.

Individual Information Distribution

The R.sub.i (b,l) matrix can be applied to each sequence used to generate the R.sub.i (b,l) itself. This produces n numbers. A histogram of the number of sites with a given information versus the information displays the R.sub.i distribution (see FIG. 2 for an example). The expectation of this distribution is, by definition, R.sub.sequence.

Variance of R.sub.i

Analogous to the mean of the R.sub.i distribution is the spread or variance of the R.sub.i distribution, given by ##EQU8## For ease of calculation, this may be reexpressed as: ##EQU9## The standard deviation of the distribution is: ##EQU10## This number measures how variable the binding sites are.

Standard Error of the Mean.

By definition, R.sub.sequence is the mean of the individual information distribution. By using the R.sub.i distribution, the standard deviation of this mean can be determined, and is known as the standard error of the mean (SEM). The SEM can be determined directly from the standard deviation of the R.sub.i distribution (.sigma..sub.R.sbsb.i) by ##EQU11## where n is the number of examples (Taylor, 1982). The variation of R.sub.sequence can also be determined by a Monte Carlo method (program Rsim, as described in detail in Stephens & Schneider, 1992 J. Mol. Biol., 228:1124-1136).

Individual information at each position in a binding site.

R.sub.i (b,l) may also be used to determine the variance at each position l in the binding site. First we define the individual information at each position l of each sequence j: ##EQU12## Since the mean at each position is: ##EQU13## we have for the variance ##EQU14## The standard deviation is: ##EQU15## Finally, the standard deviation of the mean is the variation of R.sub.sequence (l) at each position in the site: ##EQU16## This measure may have practical application for producing error bars in the sequence logo display (Schneider & Stephens, 1990 Nuc. Acid. Res., 18:6097-6100).

Searches using individual information

By applying the R.sub.i (b,l) matrix to sequences other than the sites from which it was derived, we create a search tool. Since the numerical value assigned to each position in a sequence by an R.sub.i (b,l) matrix is in bits per site, the evaluations can be directly compared to the average measures R.sub.sequence and R.sub.frequency. Because information is the only measure which allows one to add together "scores" from each position in a binding site (Shannon, 1948 Bell Systems Tech. I, 27:379-423- 623-656), other proposed search methods (Mulligan et al., 1984 Nuc. Acids Res., 12:789-800; Shapiro & Senapathy, 1987 Nuc. Acids Res., 15:7155-7174; Goodrich et al., 1990 Nuc. Acids Res., 18:4993-5000) cannot be justified.

When the R.sub.i (b,l) matrix is used for sequence searches, one must be aware that if a particular base does not appear in the data set used to create f(b,l), then f(b,l)=0 and so R.sub.i (b,l)=-.infin. at that position (see equation (1)). This expresses the fact that there are no known examples of a functioning site containing the base b at position l. That is, the simple-minded mathematics reacts as if it were very "surprised" that this is a site. This cannot happen if the matrix is only used to analyze the sequences that were used to make up the matrix itself because the infinite positions are never selected. Also, when using the dot product method, the fact that ##EQU17## f=0 assures that the infinite quantities are suppressed. Search programs can handle this situation by replacing -.infin. with a large negative value. Alternatively, the search may be relaxed by using a less severe penalty. Staden suggested replacing every f(b,l)=0 with f(b,l)=1/n (Staden, 1984 Nuc. Acids Res., 12:505-519), which allows for the possibility that the base at the position is as rare as the number of sequences used to generate the matrix. Unfortunately both -.infin. and this proposed substitution will be erroneous in most cases because the true value of the frequency will usually lie somewhere between these two extremes. The computer program of the present invention therefore allows substitution with 1/(n+t), with t.gtoreq.0. For example, using t=1 suggests that the missing base would be found if just one more binding site sequence were obtained.

The individual information method was applied to a series of situations.

Single binding site conservation distributions. The individual conservation distribution for several binding sites are shown in FIGS. 2, 3, and 4. For the splice junctions, R.sub.i (b,l) was created from the data described in (Stephens & Schneider, 1992. J. Mol. Biol., 228:1124-1136). Partially sequenced sites, which tend to make negative R.sub.i evaluations, were eliminated from the distributions shown.

Correlation of a binding site conservations with a binding energy. As an example of the use of individual information to relate sequences to binding energy, the GCN4 affinity data of Arndt and Fink was chosen. (Arndt & Fink, 1986 Proc. Natl. Acad. Sci. USA, 83:8516-8520). 28 GCN4 sites were used to create the R.sub.i (b,l) matrix. When one plots affinity directly against the R.sub.i (j), the correlation coefficient is only 0.65. Although there is still a wide scatter, the GCN4 binding sites correlate better to the logarithm of the relative affinities, having a correlation coefficient of 0.78 (FIG. 5).

The present invention provides a method for evaluating the sequences of individual binding sites. It is important to realize that the method is performed in several steps. The first step is to gather a number of example sites. These are used to generate a model of the binding sites which is called the R.sub.i (b,l) weight matrix. Because this matrix can be created from a large numbers of sequences, it can give statistically significant evaluations of individual sequences. Thus there is no contradiction: the individual sites are always evaluated in the light of a model created from a large collection of sequences.

The R.sub.i evaluation is always relative to a particular nucleic-acid recognizer. For example, each position of a given nucleic-acid sequence can be searched with an R.sub.i matrix for donor splice sites and with a different R.sub.i matrix or acceptor splice sites. Each matrix provides a different evaluation as to what its respective recognizer's response should be at every position of the sequence.

The Scan program reports the evaluation of each position in three ways: the individual information (R.sub.i), the standard deviation from the wild type distribution (Z) and the one tailed probability (p). The values of p are particularly curious because sequences with evaluations significantly higher than the mean (i.e. R.sup.sequence) have low probabilities of being real sites. There is no denying this, as it is clear from the distributions (FIG. 2, FIG. 3, FIG. 4) but it is odd because we have been socially conditioned to think that stronger binding sites are always better. They may indeed be stronger, but they are less likely to appear in the set of natural sites. Evidently the sites evolve to what is required for their function (Schneider et al., 1986 J. Mol. Biol., 188:415-431; Schneider, 1988 Maximum-Entiopy and Bayesian Methods in Science & Engineering, (Erickson, G. J. & Smith, C. R. eds) vol. 2, p.147-154, Kluver Academic Publishers Dordrecht, The Netherlands).

The computer system of the present invention comprises a processor and a memory storage device. In general, the computer system may be any IBM personal computer or compatible with operating system such a MS-DOS, PC-DOS, Windows, OS2, Unix, MacIntosh (i.e., system 7). A particularly preferred computer system is SPARCstation 20/61 with a Unix System 5 operating system (Sun Microsystems, Inc., Mountainview, Calif.). Additionally, as is readily apparent to those skilled in the art, the binding site defining system of the present invention can run effectively on currently available portable computers.

RAM: The walk program (produced from walker version 3.09) currently requires 4.2 megabytes of random access memory for a 1149 base sequence and a 21 base wide R.sub.i (b,l) weight matrix. This is within the range of many small modern computers.

DISK: The program source code sizes currently are:

walker.p 118463 bytes

scan.p 40796 bytes

ri.p 66021 bytes

Of the Delila programs, Delila is largest:

delila.p 164261 bytes

Thus a 1 gigabyte disk drive is sufficient to store the files.

Various types of database software can be used with the present invention. If it is preferred that output be produced as a printout, software exists for allowing many printers to print PostScript graphics. Any standard PostScript printer will suffice for printing the graphics from Walker.

The computer system of the present invention preferably is capable of reading a Postscript program from a file (the walk) and then switching to reading user-typed PostScript commands. One such program is the Ghostscript program, which is currently freely available from two sources. Ghostscript and Ghostview are freely available from "http://www.cs.wisc.edu/.about.ghost/index.html" and "http://155.198.1.40/gnu".

The programs are preferably compiled by a Pascal compiler such aspc, the Sun Microsystems Pascal Compiler. (See Jensen & Wirth, Pascal User Manual and Report, Springer-Verlag, New York, 1975). The source code in Appendix A through J is written in the Pascal and Postscript program languages to be portable and to avoid system dependent features. Other programming languages may be used as would be known to those skilled in the art, for example Fortran or C++.

If a Pascal compiler is not available, the Pascal code can be automatically converted to C using the p2c program. The p2c translator and library is freely available from David Gillespie ([email protected]). It can be obtained by anonymous ftp to csvax.caltech.edu in the pub directory.

The computer programs of the present invention may be stored on any computer-readable medium. Preferred types of computer-readable mediums include but are not limited to floppy diskettes, laser disks, tapes and cassettes.

Computer Program Description

One embodiment of the present invention is a method for analyzing the binding sites of macromolecules on DNA or RNA. The way data flows through various programs is shown in FIG. 7. Rectangles surround the names of programs that have been described previously. Ellipses surround the names of programs of the present invention.

The Delila Program

The data flow begins with a set of DNA sequences to be analyzed. These sequences may be obtained from GenBank or from private sources and are called a "library". They are then analyzed by programs of the Delila system (Schneider, T. D., 1982 Nuc. Acids Res., 10:3013-24; Schneider, T. D., 1984 Nuc. Acids Res., 12:129-140 and in the Delila Library System: ftp://ftp.ncifcrf.gov/pub/delila/libdef ftp://ftp.ncifcrt.gov/pub/delila/delman.ps and http://www-lmmb.ncifcrf.gov/.about.toms/delila.html)). Four files are created. The inst (instruction) file, which can be created automatically or generated by band, defines the DNA fragments and coordinates on those fragments of the binding sites to be analyzed. This set of instructions is used by the Delila program to generate a subset of the library called a book. The book contains the sequences to be analyzed. Together the inst and book files define the binding sites. These files are used by several other programs (Encode and Rseq) to create the rsdata file, which contains the initial information analysis. The information analysis at this stage is for the average of the data set, not the individuals.

The Ri program

Analysis of the individual binding sites is accomplished with the Ri program. The program is controlled by a parameter file rip and it can be given quantitative experimental data about each binding site in the values file. The output of the Ri program is given in three files. The xyin file lists the individual information content values for every sequence in the inst and book files, and these data are joined to the data from the values file. The joined data can be plotted by the xyplo program (not shown in the diagram). The raw sequences of the sites are listed in the sequ file. The ribl file contains the individual information weight matrix. This is defined in equation (1) as:

where f(b,l) is the frequency of each base b at position l in the aligned binding site sequences and e(n(l)) is a sample size correction factor for the n sequences used to create f(b,l) at position l (Schneider, et al., 1986 J. Mol. Biol., 188:415-431; Penotti, 1990 J. Mol. Biol., 213:37-52). The mathematical reasoning behind this equation is given below. R.sub.i (b,l) defines how every "finger" (l) of a protein should react to every possible base (b).

The Scan program

The ribl file is used by the Scan program to search any sequences the user is interested in (book to search). The program is controlled by parameters in the scanp file, and the output is given as a data table. The table contains a list of coordinates evaluated and the evaluation of each position (in bits of information), the number of standard deviations of each evaluation from the mean R.sub.sequence (Z score) and one-tailed probability of that Z score.

The Xyplo and DNAPlot programs

The data table from Scan may be used as the input to many programs (not shown) or it may be graphed either by the general purpose Xyplo program (which is controlled with parameters in the xyplop file) or by the specific purpose DNAPlot program which is controlled with parameters in the dnaplotp file, a positions file that can define the ends of the graph and a dnasymbols file that defines symbols to put on the graph. The advantage of DNAPlot over Xyplo is that DNAPlot can handle many pages of graphs for many sequences, but Xyplo can only make one page and use one sequence. An example of graphs generated by Xyplo and DNAPlot output is given in FIG. 8. An R.sub.i (b,l) model for the E. coli Fis protein was created as described above. The graph on the top of the figure was created by Xyplo. It shows the scan of the Fis model across the promoter region for the fis gene itself. At each step of the scan, the responses by each part of the weight matrix are added together to get the total response. This response is plotted against the position in the sequence. The plus symbols (+) indicate previously known Fis sites and the arrow shows the start and direction of transcription. The graph shows that there are several other Fis sites in this region. The lower graph, created by DNAPlot shows the scan for the entire fis gene, demonstrating that the newly predicted Fis sites cluster at the promoter.

The Walker program

The Walker program collects data from several sources. The individual information weight matrix model is read from the ribl file; colors to be used in the display are read from the colors file; parameters that define the initial display are read from the walkerp file; and the sequences to study are given in the book to search. The program manipulates these data and creates a PostScript graphics program called a walk. The walk can be shown on any PostScript device, but by using the public-domain GhostScript program it can be displayed on almost any computer system. The Walk program is carefully created so that a user can type commands (user input) in a window and receive results and help in the same window (output to user). At the same time, GhostScript displays the graphics in a second window.

An example of this display is given in FIG. 9. There are 5 horizontal rows of characters. Each row represents the placement of the individual information weight matrix for the Fis protein at a particular position on the S. typhimurium hin sequence. The DNA sequence is the same in each row. As one proceeds down the figure, the walker is stepped one position to the right on the DNA sequence so that the figure shows the frames of a "movie". Normally this would be displayed on a computer screen and only one row would be needed since the user completely controls the display in real-time. The heights of the grey letters indicate the orientation of the DNA helix, with the high points of the sine wave representing the major groove facing the protein. Horizontal grey bars are used in the region of the Walker. Note that the DNA "turns" as the movie proceeds. A pink or light green vertical bar represents the 0 coordinate of the information weight matrix. This characterizes the position of the Fis protein on the DNA. The bar is a scale, with its lowest point at -4 bits and its upper point at +2 bits. The Walker itself is shown by colored letters. Letters that extend upwards represent energetically favorable DNA contacts, while those which are upside down and extend downward represent unfavorable contacts. If a contact is more unfavorable than -4 bits, the letter is surrounded by a purple box (an example is shown in the 4throw). If a contact has never been observed at a position in the weight matrix, it is given a black box.

Three numbers are reported in the vertical bar shown in FIG. 9. The first number is the position of the bar on the sequence. The second number is the R.sub.i evaluation of the entire binding site, given in bits. This is obtained by adding together the heights of all the letters in the Walker. The third number is the Z score for this evaluation. A Z score is calculated by subtracting the mean and dividing by the standard deviation of the individual information distribution. If the Z score is below a given threshold (that can be set by the user) and the R.sub.i evaluation is positive (or greater than some value set by the user) then the bar is green to indicate that a binding site has been located. Otherwise the bar is pink. Position 180 is a known Fis binding site.

FIG. 10 demonstrates the use of the mutation feature of the Walker program to distinguish mutations from polymorphic changes (see also, Example 1). The weight matrix in this case was created from human splice acceptor sites. (See Stephens & Schneider. 1992. J. Mol. Biol., 228:1124-1136, for the details regarding how this data set was constructed). Three rows of sequence are given, but unlike the previous figure, these represent modifications of one sequence. The top sequence in FIG. 10 is the human splice acceptor site given in Fishel et al. (1993. Cell, 75:1027-1038). This is the DNA found in normal colon tissue. The middle sequence is an altered sequence found in a sporadic colorectal tumor. Fishel et al. (1993. Cell, 75:1027-1038) proposed that this T.fwdarw.C change at position -5 was the cause of the cancer, but inspection of the Walker immediately shows that this change is not significant since the R.sub.i only changes from 6.5 to 6.3 bits and the absolute value of the Z score is still below 1. Thus this change represents a polymorphism and not a mutation. The true mutation lies elsewhere or this mutation represents a change in the binding site for some molecule other than the spliceosome. The bottom row shows the effect of altering the sequence in the top row: when position -1 is changed to a cytosine ("C"), the R.sub.i becomes negative and the Z score approaches significance (p<0.02). Such an alteration would probably lead to colon cancer.

Overview of Program Actions

The specific actions of each of the programs are set forth in the Appendices A, C, E, G, and H. However, a brief overview of the activity flow is helpful for further understanding of the program's operation.

"Initialize"--gather information on a number of experimentally demonstrated example binding sites.

Align the binding sites to maximize their information content:

chose an alignment of the sequences relative to a "zero" base. (Delila programs dbbk.p, catal.p, delila.p, alist.p).

Tabulate the number of bases b at each position l, n(b,l). (Delila programs encode.p, rseq.p)

sum the n(b,l) to find the number of bases at each position, n(l). (Delila program rseq.p)

calculate from n(l) the small sampling correction factor e(n(l)) for each position. (Delila program rseq.p)

calculate a frequency matrix, f(b,l) from n(b,l)/n(l). (Delila program rseq.p)

calculate Rsequence from f(b,l) and e(n(l)). (Delila program rseq.p)

repeat the previous steps with different alignments until Rsequence is maximized. (Delila program malign.p)

generate the sequence logo. (Delila programs dalvec.p, makelogo.p)

Generate the Ri(b,l) matrix from f(b,l) and e(n(l)) (program ri.p, Appendix A; file rip, Appendix B).

if f(b,l)>0, use Ri(b,l)=2-(-log 2(f(b,l))+e(n(l))).

if f(b,l)=0, use Ri(b,l)=2-(-log 2(F(l))+e(n(l))).

where F(l)=1/(t+n(l)), with t>=0. Larger values of t are more stringent. Alternatively, the program can record "negative infinity" for the Ri(b,l) rather than stopping execution.

Evaluation of a sequence:

(program ri.p, Appendix A; file rip, Appendix B;

program scan.p, Appendix C; file scanp, Appendix D;

program walker.p, Appendix H; file walkerp, Appendix I;

file walk, Appendix J)

Obtain a sequence to be analyzed. (Delila programs: dbbk.p, catal.p, delila.p).

set the zero of the Ri(b,l) matrix at a position on the sequence.

Select the values of Ri(b,l) that correspond to the sequence.

Add these values together to obtain the individual information, Ri.

Evaluation of mean:

The mean is Rsequence determined above. (This is more reliable than the average of the Ri values unless there are no gaps in the sequence data.) (Delila program rseq.p)

Evaluation of Standard Deviation: (program ri.p, Appendix A; file rip, Appendix B)

Set the Ri(b,l) matrix at the position of each sequence used to generate the n(b,l).

Evaluate each sequence by the global Ri(b,l) matrix.

Collect the distribution in a file and calculate the standard deviation for the distribution. Scan: (program scan.p, Appendix C; file scanp, Appendix D)

Step base by base across a sequence to be analyzed. (program scan.p, Appendix C, procedure scansequence).

For a particular step, evaluate the Ri at that position (program scan.p, Appendix C, procedure scansequence).

Determine the Z score for the Ri by subtracting the mean and dividing by the standard deviation (program scan.p, Appendix C, procedure writeitout).

Determine the probability of this or a higher Z score (program scan.p, Appendix C, procedure simpson).

Record the coordinate, Ri evaluation, Z and probabilities in a data file (program scan.p, Appendix C, procedure writeitout).

Plot the data file information in a graph. (programs xyplo.p or dnaplot.p, Appendix E; file dnaplotp, Appendix F; file dnasymbols, Appendix G)

Walker: (program walker.p, Appendix H; file walkerp, Appendix I; file walk, Appendix J)

collect together information:

the sequence to analyze (program walker.p, Appendix H procedure make sequence array).

the Ri(b,l) matrix (program walker.p, Appendix H procedure makeribl).

the color scheme to use (program walker.p, Appendix H procedure varchardefs).

the overall form of the walker display (program walker.p, Appendix H procedure read parameters).

specific instructions for generation of the display (program walker.p, Appendix H procedure themain).

generate the walk graphic program described below (program walker.p, Appendix H procedure themain; file walk, Appendix J)

Running the walk program: (file walk, Appendix J) Note: commands in the walk program are implemented directly as GhostScript procedures. For example, "goto" is a procedure that the user knows about from the documentation, while "movesequence" is a procedure that the user generally does not know about.

Draw the sequence using grey in one or more lines on a graphics device. The vertical scale is in bits running from some defined lower bound in bits to zero and to 2 bits. For DNA, the letters of sequence vary in height according to a cosine wave between 1 and 2 bits high with a periodicity of 10.6 letters to indicate the helical twist of the DNA. (file walk, Appendix J, procedure movesequence).

Draw the walker either inside the sequence or next to it. When the walker is inside the DNA cosine wave is given by dashes. (file walk, Appendix J, procedure movesequence)

Evaluate each base of the sequence within the range of the walker by the Ri(b,l) matrix. These letters are colored, usually by the scheme A=green, C=blue, G=orange, T=red. When the walker is next to the sequence, the letters being evaluated are colored blue. (file walk, Appendix J, procedures evaluate, sumribl)

Draw the letters of the sequence upwards for positive Ri(b,l) evaluations. These are proportional to the evaluation and are between 0 and 2 bits (file walk, Appendix J, procedure anycolorletter)

Draw the letters for the sequence downwards for negative Ri(b,l) evaluations. These letters are drawn upside down, and range from 0 to the lower bound. Letters that extend below the lower bound are placed on a purple background. Letters for positions that have negative infinity for their evaluation are placed on a black background (file walk, Appendix J, procedure anycolorletter)

The aligning base is printed on top of a colored bar that extends from the lower to the upper bound. The bar is light green if the

program finds a binding site by the current criteria. The bar is light red (pink) if not.

The use of lighter colors is important because otherwise the letter on top of the bar would sometimes be invisible (file walk, Appendix J, procedure anycolorletter).

In the space of the colored bar opposite to the base (up or down) the coordinate, the Ri evaluation, the Z score and conceivably the probability are printed. (Evaluation of probability is currently too expensive.) (file walk, Appendix J, procedure display data)

Once the basic drawing has been made, relinquish control of the graph to the user who may then type commands. At every command the walker is redrawn as appropriate. At each step the evaluation is given not only on the walker itself but also in the window that the user uses to control the walker (file walk, Appendix J.

After all procedures have been read by the PostScript interpreter, the display is generated once by a call to toggleprinting. The user may call any procedure after that point.)

The user may move the walker or the sequence to the left or to the right by one base, by direct jumps or by a series of steps as in a movie (file walk, Appendix J, procedures h, l, jump, goto)

The user may move the walker complete lines up and down (file walk, Appendix J, procedures k, j)

The user may have the walker stay still while the sequence moves instead. (file walk, Appendix J, procedure w)

The user may move the walker in and out of the sequence. (file walk, Appendix J, procedures in, out)

The user may restructure the number of lines and bases per line on the page, the position of the entire graph on the page and

the size of the entire graph on the page. (file walk, Appendix J, procedures lines, bases, left, right, up, down, height, width)

The user may turn on and off the wave that represents DNA twist. (file walk, Appendix J, procedures waveon, waveoff)

The user may redefine the criteria for locating a binding site. (file walk, Appendix J, procedures setri, setz)

The user may instruct the program to run a search for the next or previous binding site. (file walk, Appendix J, procedures f, b)

The user may reverse the direction of the weight matrix or sequence (not yet but soon-to-be implemented).

The user may change the sequence either at an absolute

coordinate or at a coordinate relative to the current position

of the walker Ri(b,l) matrix, and immediately see the effect.

(file walk, Appendix J, procedures a, c, g, t, A, C, G, T)

The user may define delays in the display (in seconds) so that

the individual steps of the walker motion can be observed on a

fast computer.

(file walk, Appendix J, procedures setwait, isasecond)

The user may turn on and off printing and erasing of the display so that several displays can be shown on one page. (file walk, Appendix J, procedures toggleprinting, tp, toggleerase, te)

User commands may be stored in the file that defines the initial graph configuration so that figures can be generated on a printer. (file walkerp, Appendix I)

The user may ask for help, refresh the current display, restart

GhostScript on the current walk file and quit the program.

(file walk, Appendix J, procedures help, ?, r, R, q, quit)

Evaluating the Effect of sequence changes

By scan:

(program scan.p, Appendix C; file scanp, Appendix D;

xyplo.p or dnaplot.p, Appendix E;

file dnaplotp, Appendix F; file dnasymbols, Appendix G)

Scan the sequence and obtain the evaluation graph.

Modify the sequence.

Re-scan the sequence and generate another graph for the changes.

Compare the graphs to determine the effects of the changes.

By walker:

Set up the walker on the sequence of interest. (program walker.p, Appendix H; file walkerp, Appendix I; file walk, Appendix J)

Move the walker to the binding site. (file walk, Appendix J, procedures w, h, j, k, l, jump, goto)

Instruct the program to make the changes that generate the mutation. (file walk, Appendix J, procedures a, c, g, t, A, C, G, T)

Observe the change in the walker at the point of the mutation and observe the change in the evaluations that the mutations engender.

ANALYSIS

For mutation/polymorphism analysis, there are two preferred methods of analysis. With either method, a database is created containing the normal and the mutant sequence, each as a component of the same book (or separate books). In one preferred method, one may use Delila instructions to select sequences around the site of interest. The size of the region selected must be at least the size of the site defined by the R.sub.i (b,l) with the R.sub.i program but is generally larger. One can then run the Scan program on both sequences, and then may plot the normal (WT) and the mutant (MT) sites with the Dnaplot program. This will display the changes in information, including the appearance of novel binding sites or cryptic binding sites (which can be particularly important in splicing for example). This approach may be more intuitive than Walker for identification of novel or cryptic sites.

Alternatively, one can make mutations in a Postscript capable software program, such as Ghostscript, using the Walk file directly. This has the advantage of being faster particularly when there are several mutations at the same site that can be studied. A disadvantage is that it is not simple to examine mutations that result from deletions, insertions, or inversions with Walker unless a user changes many bases in the starting sequence or evaluates a book with this sequence in it. User-error is more likely when multiple sequence changes are introduced.

Polymorphic substitutions in splice recognition sites would be expected to have little or no effect on mRNA splicing, whereas true mutations reduce splicing efficiency or produce aberrant messages. Ri analysis can be used to distinguish between polymorphisms and mutations. The mean information content of 26 mutant donor splice junctions responsible for a wide variety of genetic disorders is significantly lower than the cognate wild type junctions (1.9.+-.2.2 bits versus 7.0.+-.2.4 bits; p=0.0001 by 2-tailed Student's t test). Similarly, the mean Ri for 10 mutant acceptor sites is also significantly lower (2.8.+-.2.3 bits versus 9.4.+-.3.4 bits; p=0.0001). More severe mutations involving either donor or acceptor sites tend towards lower Ri values, whereas those with a mild or moderate phenotype are likely to have information contents greater than zero, but these are still significantly less than normal sites (see FIG. 17) Mutations at normal sites with high Ri values (>12 bits) may produce non-functional sites with borderline Ri values (between 4 and 5). This observation supports the notion that while there is a minimum quantity of information needed to recognize a splice site, some sites have evolved specific requirements for nominal splicing that depend on the genic context in which they reside. For example, selection for particularly strong recognizer at the IVS2 acceptor in the human beta globin gene has been imposed by the presence of a potential cryptic acceptor sequence in the intron upstream of the normal site. A mutation at a strong splice recognition site in one context may splice appropriately in another context. Conversely, even subtle mutations at a weak splice site could make it exquisitely susceptible to loss of function regardless of genic context. In accordance with this hypothesis, it is possible to predict which genes will be affected by mutations in splice sites. Clinically, this may be useful in developing a strategy for efficient screening of various classes of mutations in particular genes, since it may permit diagnostic laboratories to determine which inherited conditions should be screened for substitutions in splice sites prior to examining other types of mutations.

Of the 49 nucleotide substitutions examined in this study, 5 polymorphic changes in splice acceptor sites were identified that were presumed in the original reports to be mutations that alter splice efficiency or the sequence of the mature mRNA. These included nucleotide changes in the familial non-polyposis colon cancer gene MSH2, the p53 gene which has been associated with some instances of bladder carcinoma, the gene encoding ornithine-transcarbamylase, and the gene encoding steroid 21-hydroxylase causing adrenal hyperplasia. To show that the change in Ri in these instances was not significantly different form the wild type sequence, splicing efficiency was categorized as either normal or severely impaired and analysis of variance on Ri was performed. Splicing was assumed to be normal if either mRNA studies demonstrated nominal splicing or levels of correct, mature message or protein were observed or the true mutation was demonstrated elsewhere in this or another gene. The Ri values for individuals with normal splicing were significantly different from those with a severe splicing defect (F test=8.85, p=0.01). This indicates that the change in Ri in the normal individuals is inconsequential, and therefore, these substitutions are genetic polymorphisms.

Measuring the Ri of mutant splice sites may permit prediction of the severity of the splicing defect. According to level 2 information theory (Schneider, 1994 Nanotechnology 5(1):1-18), the Gibbs' free energy between bound and unbound recognisers is related to information at the binding site. We therefore can compare the Gibbs energy to the R.sub.i values. We substitute the logarithm of the splicing efficiency for the energy. This is plotted for those donor and acceptor sites where quantitative studies of mRNA splicing were available (FIG. 18 and 19). The relationship is approximately linear (Correlation coefficients: for 14 donor mutations, R squared=0.60; for 9 acceptor mutations, R squared=0.40). These results provide a consistent, quantifiable approach to measuring splice efficiency.

The following examples illustrate various aspects of the present invention and in no way are intended to limit the scope thereof. All books, articles, and patents referenced herein are incorporated herein, in toto, by reference. Other similar embodiments will be clear to the skilled artisan and are encompassed within the spirit and purview of the present invention.

Example 1

ANALYZING A BINDING SITE.

As an example of this method, a T.fwdarw.C transition found at position -5 of the intervening sequence of the hMSH2 gene from multiple, independent sporadic colon carcinomas and patients with Lynch syndrome (Fishel, et al., 1993. Cell 75:1027-1038) has been analyzed by the method of the present invention. Other mutations in the coding domain of this gene cause hereditary nonpolyposis colon cancer by disrupting the repair of somatic lesions that accumulate in genomic DNA (Leach et al., 1993. Cell 75:1215-1225). Although the substitution at position -5 of the splice site was proposed to cause aberrant splicing of hMSH2 mRNA (Fishel et al., 1993. Cell 75:1027-1038), our analysis using the method of the present invention indicated that such alteration was probably not deleterious to maturation of the hMSH2 message. First, upon inspection of the sequence logo, there is a nearly equal probability of observing C or T at position -5 in this set of splice acceptor sequences (FIG. 11; this corresponds to position -6 in Fishel et al. (1993. Cell 75:1027-1038)). Second, cytosine at this position does not impede the normal splicing of 691 of 1712 acceptor sites derived from numerous human genes (Stephens & Schneider, 1992. J. Mol. Biol. 228:1124-1136). Third, we find that the common allele contains 6.5 bits of information, and the substitution weakens it to 6.3 bits. The average of the distribution of sites is 9.3 bits, and the distribution has a standard deviation of 4.6 bits. Non-functional sites are predicted to be below zero on this scale. Indeed, 2 of 20 unrelated normal individuals displayed this variant, consistent with the suggestion that this change represents a polymorphism (Leach et al., 1993. Cell 75:1215-1225).

This change is unlikely to affect the recognition of other nucleotides in the same acceptor site, as mutational analysis of the polypyrimidine tract in which it resides suggests that these nucleotides are independently recognized by the spliceosome (Stephens & Schneider, 1992. J. Mol. Biol. 228:1124-1136; Roscigno et al., 1993. J. Biol Chem. 268:11222-11229). One hundred ninety six normal human sites were found having the same or lower information content as the hMSH2 acceptor containing this substitution. 51 of these contain cytosine at position -5. Either the true mutation lies elsewhere, in this or another gene (Leach et al., 1993. Cell 75:1215-1225; Bronner et al., 1994. Nature 362:258-261; Papadopoulos et al., 1994. Science 263:1625-1629), or the change indicates that this base is involved in a genetic control mechanism other than mRNA splicing (Amrein et al., 1994. Cell 76:735-746).

To summarize, inference of genetic mutations in splice junction recognition sites based on consensus sequences may be inaccurate, whereas information analysis of sequence variants can distinguish between polymorphic nucleotides and mutant sites. True mutations are expected to reside in positions where the sequence conservation in bits significantly exceeds the background variation and where the base frequency decreases significantly.

A similar approach may be applied to the analysis of other conserved transcriptional and translational signals or protein motifs in human sequences.

Example 2

Roscigno et al. (1993 J. Biol. Chem., 268:11222-11229) determined the effect of making changes in the polypyrimidine tract of adenovirus 2 intron of the major late promoter Leader 1 and Leader 2 splicing unit (GenBank accession J01917 coordinate 7100, Adenovirus type 2 DNA). They mutated this site and measured the splice product RNA divided by the wild-type product produced. These data and their standard deviations were measured from their graphs (See Roscigno FIGS. 3 and 4) in millimeters. The logarithm of these values were plotted against the predicted R.sub.i values (FIG. 21). One case of zero splicing was removed because the logarithm cannot be taken, and because small amounts of splicing may have occurred but were not reported. The correlation coefficient is 0.81. This case demonstrates that the R.sub.i analysis can predict the strength of a splice site within the experimental error.

Example 3

This example demonstrates the use of the present invention as a tool for identifying binding sites and manipulating the affinity of a binding site by specific changes within positions of the sequence.

In this example, Fis binding sites are analyzed. Fis is a bacterial protein which functions by binding to specific binding sequences on DNA and bending DNA in sitespecific recombination systems. The resulting information content model is used to locate previously unidentified sites adjacent to known ones. DNA mobility shift experiments were then performed to determine if the predicted sites are bound by Fis in vitro.

Searching sequences with the Fis individual information matrix model. The programs Scan, Xyplo and DNAplot were used to study Fis binding site on the fis promoter. At the transcription initiation site of the fis promoter, there are 6 strong Fis sites (Ball et al., 1992 J. Bact., 17: 8043-8056; Ninnemann et al., 1992 EMBO J. 11: 1075-1083). The Scan results show up to 13 additional sites in the immediate region of the promoter, but few elsewhere on the gene (FIG. 8). Presumably these correspond to the weaker sites noted by Ball et al. (1992 J. Bact., 174: 8043-8056).

In the bacteriophage P1 cin, bacteriophage P7 cin, and E. coli el4 pin enhancers, a potential overlapping site occurs 7 base pairs (.about.1/2 helical turn) to the left of the previously identified proximal site. (FIG. 12, right three graphs). Since this potential site is outside the region between the proximal and distal sites, we named it the "external" site. When a new site is on the right, it is 11 bases from the previously identified site, while a new site on the left is 7 bases from the previously identified site. We do not know if this correlation is coincidental. We also observed that a pattern corresponding to site III in gin (Koch et al., 1991 Nuc. Acids. Res., 19:5915-5922) appears in all other enhancers scanned except hin that in three cases a weaker potential site falls exactly between the distal site and that site with spacings of 10 to 13 base pairs (FIG. 12). Because the nomenclature for binding sites is already obscure, we decided not to name these sites.

Two Fis sites have been identified in the E. coli oriC locus at coordinates 202 (8.2 bits) and 283 (5.7 bits) (Filutowicz et al., 1992 J. Bact., 174:398-407). There is another strong potential Fis site exactly 11 bases from the 202 site at coordinate 213 (8.0 bits). Footprinting data in Filutowicz (1992, FIG. 5b, c, site "I") shows DNase I protection that covers both sites.

Total sequence conservation at Fis sites. The total number of Fis sites in the E. coli genome is not known, so the information needed to locate those sites (R.sub.frequency) cannot be calculated. However, the total sequence conservation at the binding sites is 8.5 bits, which suggests that there is one site every 2.sup.85 =362 bases or an average of 4 sites at each of the 3239 genes of the 4,673,000 bp genome. It also implies that about 1300 Fis molecules would be needed to fill the Fis sites. When we searched Ecoseq7, which contains 60% of the known E. coli sequences (Rudd, et al. 1993, ASM News, 59: 335-341), for Fis sites with more than 1 bit of sequence conservation we found 36,000 sites, so there should be 60,000 possible Fis sites in the entire genome. These estimates are comparable to the number of Fis molecules per cell, which ranges from close to zero in stationary cells to between 50,000 and 100,000 molecules per cell during the transition to exponential growth or an increase in nutrients.

Example 4

This example provides an illustration of designing binding sites with the method of the present invention.

We chose 32 bases of the hin sequence because according to the information-theory based search this region contains two overlapping Fis sites, one of which is the Fis site proximal to the recombination junction hixL (Bruist et al., 1987 Gens Dev., 1:762-772). We added 5 bases on each end--half a twist of DNA--to be sure we were not missing important components, although this region does not show up significantly in the sequence logo. Beyond these ends we added EcoRI and HindIII overhangs. We created three other sequences using the anticonsensus of the Fis sequence to destroy the proximal site, the newly identified "medial" site, or both sites (FIG. 13). The anticonsensus sequence is the sequence which should bind Fis the worst. It is predicted from the number of bases at each position (n(b, l) numbers matrix or the R.sub.i (b,l) weight matrix, FIG. 14) by noting which bases appear least frequently at each position of the site. In ambiguous cases we chose C or G when possible because these appear rarely in the logo (FIG. 15). We used the same rationale in designing the DNA from bacteriophage Pl cin.

These sequences and their complements were synthesized (Midland Oligos, Inc., Midland, Tex., USA) so that when annealed they provide sticky EcoRI and HindIII ends. Annealed oligos were ligated into plasmid pTS385 digested with EcoRI and HindIII and transformed into E. coli DH5.alpha. as previously described (Hengen & Iyer, 1992 Biotechniques, 13:57-62). Transformants were selected on LB media containing 50 .mu.g/ml kanamycin and 50 .mu.g/ml of ampicillin. When necessary, we transformed E. coli BL21/DE2 (Studier & Moffatt, 1986 J. Mol. biol., 189:113-130) and selected them on the same media containing 1 mM IPTG. We knew from previous experiments that the parental plasmid pTS385 is conditionally lethal to this strain because a strong T7 promoter is positioned between the EcoRI and HindIII sites. Induction of T7 RNA polymerase with IPTG thus provided a strong selection for recombinant plasmids containing the intended insert Fis DNA, eliminating all but a few wild-type pTS385 plasmids from the lot of transformants. The resultant plasmids were screened by restriction analysis and PCR amplified using primers flanking the inserted DNA pTS37fl 5' acatttcccgaaaagtgc 3' and pTS37rl 5' cggaacacgtagaaagcca 3'. When recombinants were identified, plasmid DNA was transformed into and maintained within E. coli DH5.alpha.. The sequence between the EcoRI and HindIII sites was then confirmed by dideoxy sequencing with an ABI model 373A automated sequencer (Hunkapiller et al., 1991 Science, 254:59-67).

For gel mobility shifts (Fried & Crothers, 1981 Nuc. Acids Res. 9(23):6515-6525), we used Fis protein cloned and purified from E. coli obtained as a gift from R. Johnson 1986 Cell, 46:531-539). Plasmid DNA from the 8 clones was purified by the method of Birnboim and Doly (1979 Nuc Acids Res. 7:1513-1523) or Hengen (1995 Biotechniques, 13:57-62), digested with EcoRI, end-filled with biotin-11-dUTP using the Klenow fragment of E. coli DNA polymerase I, and linearized with BglII, which cleaves 369 bp from the EcoRI site. The 369 bp DNA fragment was purified away from the larger plasmid fragment by electrophoresis through SeaPlaque GTG agarose (FMC, Rockland, Me., USA), sliced from the gel, and extracted using a freeze-and-spin method through Costar Spin-.RTM. centrifuge tubes containing 0.2 .mu.m pore size nylon filters. Purified DNA was extracted with an equal volume of isoamyl alcohol to remove residual ethidium bromide, digested with HindIII, heated to 65.degree. C. for 30 minutes to inactivate the HindIII enzyme, and cooled to room temperature.

Binding assays were accomplished by incubating DNA at approximately 1 nM with various concentrations of Fis protein ranging from 125 to 1000 nM at room temperature for 15 minutes in 25 mM Tris HCl (pH 7.6), 80 mM NaCl, 1 mM EDTA, 2 mM DTT 100 .mu.g/ml acetylated bovine serum albumin, and 100 .mu.g/ml calf thymus DNA. Gel shift analysis was done by separation of the different species on a 8.0% polyacrylamide gel in 1.times.TBE. The DNA was electro-transferred onto Tropilon-Plus.TM. nylon membrane (Tropix, Inc. Bedford, Mass., USA) with a Hoefer Semi-Phor Model TE70 semi-dry transfer unit for 30 minutes at 30 mA, and crosslinked by exposure to 254 nm UV light for 10 minutes on a UVP Model T5-15 transilluminator (Ultra-Violet Products, Inc., San Gabriel, Calif., USA). Biotinylated DNA was detected using a Southern-Light.box-solid. chemiluminescent kit using the CSPD.RTM. substrate (Tropix, Inc., Bedford, Mass., USA) and exposure to Kodak BioMax MR film.

Strong Fis sites separated by 11 and 7 base pairs were designed by selecting the most frequent base at each position in the Fis sequence logo (FIG. 6, FIG. 9). These were then merged with the same sequence shifted by 11 or 7 base pairs. 5 extra bases were added to the ends and the DNAs were made self complementary (FIG. 8). They were synthesized with biotin on the 5' end and gel purified (Oligos Etc. Wilsonville, Oreg., USA). To insure complete annealing, they were heated, and slowly cooled to room temperature.

To determine whether overlapping Fis sites can be simultaneously bound by Fis, we synthesized strong Fis sites which overlap by 11 or 7 base pairs and tested their properties by gel shift. Fis protein shifts both DNAs, but the DNA with two Fis sites separated by 11 bases was shifted once, while the DNA with two Fis sites separated by 7 bases is shifted twice. This demonstrates that Fis molecules separated by 11 bases are on the same face of the DNA and collide with each other, while those separated by 7 bases are on different faces and do not collide. These results are consistent with many observations of Fis sites naturally separated by 7 or 11 bases, with molecular modelling and with a detailed analysis of the sequence logo structure which reveals that the Fis sites are internally redundant at spacings of 7 and 11 bases. We therefore propose that the collision and non-collision properties of Fis are used in genetic control systems as part of molecular flip-flops. Such flip-flops may be useful for constructing molecular computers.

Example 5

This example demonstrates the relationship between information content and binding ability.

Fis sites at inversion regions. When we scanned our Fis R.sub.i (b,l) model across DNA inversion regions, we discovered that each known proximal site had an overlapping sequence with the same characteristics as a Fis site (FIG. 12). To test whether the new sites exist, we performed gel shift experiments on DNAs in which we presumably had knocked out neither, one, or both of the sites.

Under our experimental conditions hin does have a second site as predicted, since the knockout of the stronger proximal site still allowed the DNA to shift (FIG. 16). However, more Fis protein was required to shift an equivalent amount of DNA than for the wild-type proximal site, indicating that Fis binds weakly to the medial site. This is consistent with the weaker sequence conservation of the medial site (4.5 bits) compared to the proximal site (9.0 bits).

For the cin experiment, the stronger proximal site was confirmed, but the weaker external site showed a barely detectable shift (visible on the original X-ray film). To our surprise, when both cin sites were destroyed, we still detected a weak shift, which is stronger than that of the external site. In the process of destroying both the external and proximal sites, we inadvertently created a new Fis binding site shifted one base to the left of the original external site. The new site is only 1.3 bits, but it still gives a band shift.

Example 6

EXAMPLE OF SPLICE MUTATIONS

Six donor site mutations have been examined, all of which cause beta+ thalassemia, i.e. there is some normal splicing.

Three of these mutations are in exon 1 of the beta globin gene and give a mild thalassemia phenotype. The normal intervening sequence 1 ("IVS1") donor at position 246 has 4.96 bits. There are two cryptic sites in the normal sequence that are apparently not used in vivo, but which are more likely to be used either if the position 246 site is mutated to become weaker or if mutations occur that make them stronger.

The cryptic sites are at positions 208 (7.69 bits) and 230 (8.73 bits), i.e. in exon 1. Mutations at position

228 (t.fwdarw.a) increases the site at 230 to 10.86 bits

232 (g.fwdarw.a) increases the site at 230 to 9.14 bits

235 (g.fwdarw.t) increases the site at 230 to 9.96 bits

The difference in information content between the normal and mutant sequences appears to be rather small, as is the phenotypic affect.

Conversely, mutations in the donor site itself, even ones that are somewhat removed from the splice site result in preferential splicing at these cryptic sites. At position 251 (i.e. +5), G.fwdarw.C results in a reduction to 1.01 bits and G.fwdarw.T results in a reduction to 1.04 bits. Patients with these mutations have beta-plus thalassemia, but splicing at this site is severely reduced compared to normal. In contrast, T.fwdarw.C mutation at position 252 (i.e. +6) results in a reduction to 3.54 bits. This mutation is not a severe beta+ thalassemia, with splicing of the normal message occurring at 50-70% of the wild-type splice site.

It may be useful to use Scan to analyze for cryptic splice sites in the normal sequence close to the splice donors and acceptors that are normally used for all of the human genes in the database. Then, a correlation can be made to the disease database of splice mutations with that list to see whether those splice mutations are more severe than others where no such cryptic sites can be found.

Two mutations have been found in intron 1 which activate cryptic acceptors: g355a and t362g, upstream from the one normally used. The site created by g355a has 4.89 bits and has a beta+ thalassemia phenotype. In monkey kidney cells--not erythroid cells, the cryptic site is used 90% of the time, the normal site 10%. In erythroid cells the abnormal message is not detected, but processed mRNA levels are lower than normal. The site created by t362g has 5.08 bits and the normal site is not used in the heterologous expression system. This would be interpreted as a beta-0 thalassemia, except that the cell type in which splicing is analyzed appears to be important, so it may not be possible to draw the inference of beta-0 thalassemia. There appears to be a minimum threshold of information required for choice of the splice acceptor, but as long as the cryptic acceptor falls within the normal range it can and will be used.

An interesting cryptic acceptor site in intron 2 has been identified. The normal intron sequence contains a splice acceptor site at 1177 that is stronger than the one adjacent to exon 3 (position 1448). The site at 1177 has 14.779 bits and the one at 1448 has 13.33 bits. An A.fwdarw.G mutation at 1447 has been described which has a beta-0 (no mature globin mRNA) phenotype. This mutation reduces information content to 5.17 bits at the normal splice site (curiously, one is created at 1446 with 7.046 bits). Note that both of these are in the normal range. However, neither can compete with the cryptic site at 1177, so that essentially all of the spliced message is untranslatable and unstable. This site is so strong that mutations that create new donor sites between 1177 and 1448 create an untranslatable exon with the 1177 as 5' end (then, the IVS IVS2 donor splices to 1177 instead of 1448).

These two examples are paradoxical. In the first intron, the cryptic sites are weaker by Ri analysis than the normal acceptor but they are preferred. In the second intron, the cryptic acceptor is stronger than the "mutant" site in the normal acceptor and is preferred even though the "mutant" site has respectable information content. These results are reconciled in that the spliceosome processively reads the sequence until it finds an acceptable site (from 5'.fwdarw.13') and makes a lariat.

Example 7

Mutations at the +3 position of the donor splice site in different genes were analyzed. Specifically, the sequence alterations were G.fwdarw.T in Von Willebrand Factor mRNA, A.fwdarw.G in Ornithine Transcarbamylase mRNA ("OTC") (exon 7), and G.fwdarw.C in CD18 (beta integrin). The first causes a form of hemophilia, the second--congenital hyperammonemia, mental retardation and usually infantile death, and the third, recurrent often fatal infections due to deficient expression of leukocyte adhesion glycoproteins. The severity of these different diseases is or appears to be correlated with the splice site mutation present.

In the Von Willebrand Factor mutation, there is exon skipping because the splice site is not recognized in some instances (and because there are no cryptic site in the neighborhood--which is confirmed by the scan). The normal site has 10.07 bits while the mutation has 5.97 bits. Experimentally it appears that the affected homozygote (by RT-PCR) makes similar amounts of mutant and normal transcripts. Clotting, however, is markedly reduced in the homozygote due to low levels of factor present. This may be related more to the turnover and stability of the factor (which is found in plasma).

In the OTC patient, the substitution does not change the information content of the site very significantly. The previous normal site has 6.954 bits, the "mutant" has 6.554. The Northern and Western blots do not demonstrate a reduced expression and sythesis of OTC in this patient. Also reported is a T.fwdarw.C substitution at the invariant +2 site which does abolish expression experimentally and has -11.138 bits of information. The OTC (+3) change, represents another polymorphism.

The last mutation, in CD18, was found in a set of related individuals with moderate deficiency phenotype. This mutation does not completely abolish splicing: however the level was measured to be 3% of that seen in normal individuals. The normal splice site has 9.179 bits and the mutant site has 4.78 bits, which appears to be towards the low end of the distribution.

Example 8

Steroid 21-hydroxylase gene splice site substitutions in intron 2 (IVS21).

Mutations in the Cytochrome P450 (C21) (which encodes Steroid 21-hydroxylase), cause congenital adrenal hyperplasia ("CAH"), a recessive disorder. Patients with this disorder display a virilizing phenotype or a salt wasting phenotype. Virilization is more apparent in females, in males it can result in precocious puberty and hypersexualization. Most of the mutations characterized to date result form gene conversion of the B gene by the neighboring A gene, which is non-functional pseudogene. These two genes are very similar in sequence, there are numerous nucleotide substitutions in the A gene that when introduced into the B gene by gene conversion result in a non-functional P450(C21)B allele. Depending on the extent of the gene conversion event, the mutated sequences may affect the entire B gene or a subset of sequences in this gene.

Higashi et al. (Proc. Nat Acad Sci. USA, 85: 7486-90, 1988) described two patients with CAH that exhibited substitutions in the acceptor sequence of IVS2 of P450 (C21)B. Patient 10 was a virilized female with a C.fwdarw.G transversion at position -12 of the normal splice site (pos. 2333 of GenBank locus M12792). Patient 7 was a salt-wasting male with a C.fwdarw.A transversion at the same site. The substitutions in both of these individuals arose by gene conversion of the 5' or amino-terminal domain of the B gene by the A-pseudogene. The 3' terminal segment of the CAH gene was not involved in the gene conversion event. This led the investigators to suspect hat these amino-terminal nucleotide substitutions may have been responsible for inactivation of these CAH alleles. S1 nuclease protection studies show that the C.fwdarw.G substitution abolishes mRNA splicing at this acceptor and results in the exclusive use of 3 preexisting and new cryptic acceptor sites upstream of the normal site, and premature termination of translation.

Individual information analysis of these substitutions is consistent with the S1 nuclease protection experiments. The C.fwdarw.G substitution creates an adequate cryptic acceptor with 7.99 bits of information from a site with 0.70 bits. The normal acceptor decreases slightly from 12.1 to 10.5 bits, however it is within the range of functional sites. In order to explain the preference for using the cryptic acceptor over the normal site even though it has a lower R.sub.i value, it would appear that the cryptic site is detected by the spliceosome prior to seeing the normal site. This preference for a weaker, but adequate cryptic acceptor has been seen at similar mutations in several other genes that we have been examined and may be a consequence of processivity of the spliceosome in recognizing acceptor sites.

In contrast, the C.fwdarw.A allele does not create a new splice recognizer sequence at position -12 (there is a small decrease in the Ri at this site compared to the normal sequence to 0.41). It does not appreciably reduce the information content of the normal acceptor site either (form 12.1 to 10.0 bits), which is within the range of functional sites. This analysis indicates the C.fwdarw.A is a genetic polymorphism independent of the S1 nuclease digestion result. The prevalence of this substituter in patients with CAH is therefore unrelated to the diagnosis. We would predict that if a similar number of normal individuals without evidence of this disorder were examined, this substitution would also be detected frequently.

Example 9

Fis sites at the nrd promoter demonstrate prediction of sites for which footprinting data exist.

By using a degenerate consensus pattern, 5 Fis binding sites were found upstream of the transcriptional start site of the nrd operon of E. coli (Augustin, et al., 1994 J. Bact., 176:378-387). When we scanned for potential Fis binding sites, about 8 more sites were identified (FIG. 20).

These sites were easily confirmed to be true sites since Cu-phenantroline footprinting of this region had been carried out by (See Augustin et al., FIG. 3), corresponding well with our predictions even though none of these sites were used in the generation of the R.sub.i (b,l) weight matrix model used in this analysis.

The DNA sequence was from the GenBank; accession number K02672 (Carlson, et al., 1984, Prac. Natl. Acad. Sci. USA, 81:4294-4297).

Transcription begins at position 0 (GenBank coordinate 3395) and proceeds to the right. Potential Fis sites (R.sub.i .gtoreq.2 bits) relative to the start of transcription are at:

-349 (2.0 bits), 348 (6.6 bits), -327 (8.7 bits), -283 (13.8 bits), -272 (5.8 bits), 230 (8.9 bits), -221 (2.2 bits), 209 (3.2 bits), 202 (8.8 bits), -173 (2.6 bits), -158 (4.4 bits), -129 (6.4 bits) and -17 (4.8 bits).

Five Fis sites were identified by Augustin et al. to be in the ranges: -310 to -328 (probably site 327), -268 to -285 (site -272), -187 to -204 (probably site -202), -142 to -160 (probably site -158), and -122 to -139 (site -129) relative to the start of transcription. These are indicated by filled squares in FIG. 20.

Those sites which were located by the Scan program and visible on the footprinting data of Augustin et al. (Augustin, et al. FIG. 4, lanes 4 and 5) but not previously described, are indicated by filled circles.

The two DNA sites found by Augustin et al. are at -52 and -40 and indicated by open squares.

We extracted this sequence using a new feature of the Delila program which allows sequences to be renumbered by giving the command "default coordinate zero" followed by instructions of the form "get from 3395 -4000 to 3395 +4000;". Thus when this sequence was searched with the Scan program, the reported locations were relative to the transcriptional initiation point.

This plot also differs from those of FIG. 12 in that the individual information scores are drawn as lines from the bottom up, rather than from zero bits up or down. This is set by using a switch within the DNAplot parameter file.

It will thus be seen that the objects set forth above, among those made apparent from the preceding description, are effectively obtained. Since certain changes may be made above system and method without departing from the scope of the invention, it is intended that all matter contained in the above description and shown in the accompanying drawing shall be interpretive as illustrative and not limiting. It is also understood that the following claims are intended to cover all of the generic and specific features of the invention described herein, and all statements of the scope of the invention which, as a matter of a language might be said to fall between. ##SPC1##

Claims

We claim:

1. A method of determining sequence conservation of a binding site on a test nucleic acid sequence using a set of sequences, said set including nucleic acid base identity, base position and alignment information, comprising: (a) extracting base identity and alignment information to generate an aligned set of sequences corresponding to at least one pre-determined criterion; (b) generating in a computer an information weight matrix for said extracted information to provide an information model for said nucleic acid sequences, said matrix being calculated in accordance with the formula: R.sub.i (b,l)=2-(-log.sub.2 f(b,l)+e(n(l)), where .function.(b,l) is the frequency of each base b at position l in the aligned sequences and e(n(l)) is a sample size correction factor for the n sequences in .function.(b,l); (c) applying in a computer the individual information weight matrix R.sub.i to the test sequence, thereby assigning each position within the test sequence an information-based numerical value; and (d) determining the sequence conservation of the binding site on the test sequence by adding together the information-based numerical values.

2. The method of claim 1 further comprising the step of changing at least one base in said test sequence to be analyzed and repeating the steps (c)-(d).

3. The method of claim 2 wherein said model is a splice site donor sequence model.

4. The method of claim 2 wherein said model is a splice site acceptor sequence model.

5. The method of claim 1, further comprising outputting the evaluation of the matrix in at least one of the following: (i) standard deviation of the individual base information from the wild-type n sequences, and (ii) the one tailed probability.

6. The method of claim 1, further comprising: outputting said matrix to a graphic interface.

7. The method of claim 6 wherein said sequence is DNA and the amplitude of the numerical value signal displays a sinosoidal wave having maximum values between 1 and 2 bits high having a period of 10.6 positions.

8. A method of claim 1, further comprising the step of displaying nucleic acid sequence information as a Walker, said information including base identity, base position, DNA helix angle and information weight matrix for said sequence, comprising displaying the sequence in a set of letters, wherein the amplitude of each letter at each base corresponds to the individual information weight of the letter at the base position.

9. The method of claim 8 wherein the display includes an indicia at one base position, a portion of said indicia defining the base position, individual information, and Z score of the Walker.

10. The method of claim 9 wherein a comparison between the individual information content signal for the region and a mean value signal of information is used to generate a Z score.

11. The method of claim 8 wherein the indicia include a letter representation of the base at that position.

12. The method of claim 8 wherein the Walker amplitude is positive for energetically-favored bases.

13. The method of claim 8 wherein the base position number of at least one base is displayed digitally.

14. The method of claim 13 wherein the individual information content of at least one base position is displayed.

15. The method of claim 14 wherein said display is digital.

16. The method of any claim 8, 9, 11-15 wherein said Walker is colorized.

17. A computer program product comprising: a computer usable medium embodying computer readable program code means for analyzing a test nucleic acid sequence using aligned nucleic acid sequences, the program code means comprising: program code means for calculating the individual information weight matrix for said aligned sequences in accordance with the formula: R.sub.i (b,l)=2-(-log.sub.2 f(b,l)+e(n(l)), where .function.(b,l) is the frequency of each base at position l in the aligned sequences and e(n(l) is a sample size correction factor for the n sequences in .function.(b,l).

18. The computer program product of claim 17, wherein said medium is a tape medium.

19. The computer program product of claim 17, wherein said medium is a CD-Rom medium.

20. The computer program product of claim 17, wherein said medium is a random access memory.

21. The computer program product of claim 17, wherein said program code means is readable by a digital computer.

22. The computer program product of claim 17, wherein said program code means is in the PASCAL programming language.

23. A computer system having a central processing unit under the control of the program of claim 17.

24. A method for characterizing a binding site, utilizing a processor capable of generating an R.sub.i (b,l) information weight matrix signal and a display means capable of displaying screens associated with a plurality of functions comprising the steps of: (a) applying said information weight matrix signal to a region of a sequence corresponding in size to the binding site; (b) assigning each position within the region, an information-based numerical value signal based upon the information weight matrix signal; (c) displaying said numerical value signal on said display means; and (d) characterizing the binding site by accumulating the information-based numerical value signal for each position within the region, producing an individual information content signal.

25. The method of claim 24 wherein said applying step is performed in multiple regions corresponding in size to the binding site.

26. The method of claim 1 or 24, wherein the binding site is characterized as having a polymorphism.

27. The method of claim 1 or 24, wherein the binding site is characterized as having a deleterious mutation.

28. A method of identifying binding sites utilizing a processor having a data entry means comprising the steps of: (i) generating an information weight matrix signal based upon known binding site sequences; (ii) applying a region corresponding in size to said binding sequence of an unknown sequence signal to the information weight matrix signal; (iii) analyzing each position within the unknown sequence signal to determine information content signal at said position; (iv) adding the information content signals together thereby generating an individual information content signal; and (v) identifying a binding site upon generation of a positive individual information content signal.

29. The method of claim 28 wherein the individual information weight matrix signal is applied step wise to each position within a sequence thereby generating a series of individual information content signals.

30. The method of claim 29 further comprising a display means wherein said series of individual information content signals are displayed in graphical form.

31. An article of manufacture for use in determining sequence conservation of a binding site on a test nucleic acid sequence using an aligned set of sequences comprising a computer readable medium, said medium containing individual information weight matrix generation signals for forming an information model corresponding to at least one predetermined criterion, said individual information weight matrix being defined by the formula: R.sub.i (b,l)=2-(-log.sub.2 f(b,l)+e(n(l), where .function.(b,l) is the frequency of each base b at position l in the aligned sequences and e(n(l) is a sample size correction factor for the n sequences in .function.(b,l); said medium further comprising individual information weight matrix generating signals for generating an individual information weight matrix from the aligned sequences containing the binding site in accordance with said formula, thereby assigning each position within the aligned sequences an information-based numerical value; said medium further containing characterizing signals for determining the sequence conservation of the binding site in the test nucleic acid sequence based upon the individual information weight matrix.

32. The article of manufacture according to claim 31, wherein said model is a splice site donor sequence model.

33. The article of manufacture according to claim 31, wherein said model is a splice site acceptor sequence model.

34. The invention of claim 31 , wherein said model is a promoter sequence model.

35. The invention of claim 31, wherein said model is a regulatory sequence model.

36. The article of manufacture of claim 31, wherein the binding site is characterized as having a polymorphism.

37. The article of manufacture of claim 31, wherein the binding site is characterized as having a deleterious mutation.

38. A method of displaying nucleic acid sequence information as a Walker, said information including base identity, base position, DNA helix angle and information based numerical values for said sequence, comprising: displaying the sequence in a set of letters, wherein the amplitude of the letter at each base corresponds to the individual information weight of each base at the base position.

39. The method of claim 38 wherein the Walker includes an indicia at a base position.

40. The method of claim 39 wherein the Walker includes the alpha letter indicia corresponding to the base.

41. The method of claim 38 wherein the Walker amplitude is positive for energetically-favored bases.

42. The method of claim 38 wherein the base position number of at least one base is displayed digitally.

43. The method of claim 42 wherein the information content of at least one base position is displayed.

44. The method of claim 42 wherein said display is digital.

45. The method of claim 38-44 wherein said Walker is colorized.

46. A computer program product comprising a computer usable medium embodying computer readable program code for displaying nucleic acid sequence information as a Walker, said information including base identity, base position, DNA helix angle and information weight matrix for said sequence, the program code means comprising: program code means for displaying the sequence as a Walker, wherein the amplitude of the Walker at each base corresponds to the individual information weight at the base position.

47. The computer program product according to claim 46, wherein said medium is a tape medium.

48. The computer program product according to claim 46, wherein said medium is a CD-Rom medium.

49. The computer program product according to claim 46, wherein said medium is a random access memory.

50. The computer program product according to claim 46, wherein said medium is an optical disk.

51. The computer program product according to claim 46, wherein said program code means is readable by a digital computer.

52. The program product of claim 46 further including program code means for displaying an indicia at a base position.

53. The program product of claim 46 further including program code means for displaying the alpha letter indicia corresponding to the base.

54. The program product of claim 46 further including program code means for displaying the Walker amplitude as positive for energetically-favored bases.

55. The program product of claim 46 further including program code means for displaying the base position number of at least one base is digital.

56. The program product of claim 46 further including program code means for displaying the information content of at least one base position.

57. The program product of claim 46 further including program code means for displaying said Walker in color.

Non-Patent Literature (13)

Stephens et al., "Features of Spliceosome Evolution and Function Inferred from Analysis of the Information at Human Splice Sites", Journal of Medical Biology, 1992, vol. 228, pp. 1124-1136.
Schneider et al., "Sequence logos: a new way to display consensus sequences", Nucleic Acids Research, 1990, vol. 18, pp. 6097-6100.
Staden, R., "Computer methods to locate signals in nucleic acid sequences", Nucleic Acids Research, 1984, vol. 12, pp. 505-519.
Schnider, T. "Poster on Information Theory and Individual Information," Web page at http://www-lmmb.ncifcrf.gov/.about.toms/just.say.no/index.html, Feb. 27, 1996--updated Oct. 18, 1996.
Thomas D. Schneider, "Sequence Logos, Machine/Channel Capacity, Maxwell's Demon, and Molecular Computers: A Review of the Theory of Molecular Machines," Nanotechnology vol. 5, pp. 1-8, 1994.
Nathan D. Herman, et al., "High Information Conservation Implies that at Least Three Proteins Bind Independently to F Plasmid incD Repeats," Journal of Bacteriology, pp. 3558-3560, Jun. 1992.
Thomas D. Schneider, et al., "Sequence Logos: A New Way to Display Consensus Sequences," Nucleic Acids Research, vol. 18, No. 20, pp. 6097-6100, Sep. 12, 1990.
Otto G. Berg, et al., "Selection of DNA Binding Sites by Regulatory Proteins: The Binding Specificity of Cyclic AMP Receptor Protein to Recognition Sites," Journal of Molecular Biology, No. vol. 200, No. 4, pp. 709-723, 1988.
Otto G. Berg, "Selection of DNA Binding Sites by Regulatory Proteins: The LexA Protein and the Arginine Repressor Use Different Strategies for Functional Specificity," Nucleic Acids Research, vol. 16, No. 11, pp. 5089-5105, Apr. 26, 1988.
Thomas D. Schneider, "Information and Entropy of Patterns in Genetic Switchs," Maximum-Entropy and Bayesian Methods in Science and Engineering, vol. 2, pp. 147-154, 1988.
Otto G. Berg, et al., "Selection of DNA Binding Sites by Regulatory Proteins: Statistical-Mechanical Theory and Application to Operators and Promoters," Journal of Molecular Biology, vol. 193, No. 4, pp. 723-750, Feb. 20, 1987.
Rodger Staden, "Computer Methods to Locate Signals in Nucleic Acid Sequences," Nucleic Acids Research, vol. 12, No. 1, pp. 505-519, 1984.
Gary D. Stormo, et al., "Use of the `Perceptron` Algorithm to Distinguish Translational Initiation Sites in E. coli," Nucleic Acids Research, vol. 10., No. 9, pp. 2997-3011, 1982.

Abstract

Description