HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Full Text (PDF) All Versions of this Article: rna.374907v1 13/4/478    most recent Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Kim, N. Articles by Schlick, T. Search for Related Content PubMed PubMed Citation Articles by Kim, N. Articles by Schlick, T. Social Bookmarking                   What's this?

A computational proposal for designing structured RNA pools for in vitro selection of RNAs

¹ Department of Chemistry, New York University, New York, New York 10003, USA
² Courant Institute of Mathematical Sciences, New York University, New York, New York 10012, USA

	ABSTRACT

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION APPENDIX ACKNOWLEDGMENTS REFERENCES

 
Although in vitro selection technology is a versatile experimental tool for discovering novel synthetic RNA molecules, finding complex RNA molecules is difficult because most RNAs identified from random sequence pools are simple motifs, consistent with recent computational analysis of such sequence pools. Thus, enriching in vitro selection pools with complex structures could increase the probability of discovering novel RNAs. Here we develop an approach for engineering sequence pools that links RNA sequence space regions with corresponding structural distributions via a "mixing matrix" approach combined with a graph theory analysis. We define five classes of mixing matrices motivated by covariance mutations in RNA; these constructs define nucleotide transition rates and are applied to chosen starting sequences to yield specific nonrandom pools. We examine the coverage of sequence space as a function of the mixing matrix and starting sequence via clustering analysis. We show that, in contrast to random sequences, which are associated only with a local region of sequence space, our designed pools, including a structured pool for GTP aptamers, can target specific motifs. It follows that experimental synthesis of designed pools can benefit from using optimized starting sequences, mixing matrices, and pool fractions associated with each of our constructed pools as a guide. Automation of our approach could provide practical tools for pool design applications for in vitro selection of RNAs and related problems.

Keywords: in vitro selection; RNA pool design; mixing matrix; sequence-structure map; graph theory

	INTRODUCTION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION APPENDIX ACKNOWLEDGMENTS REFERENCES

 
In vitro selection is an experimental approach that allows the screening of large random-sequence libraries of nucleic acid molecules (10¹⁵) for a specific function, such as binding or catalysis (Ellington and Szostak 1990; Tuerk and Gold 1990; Wilson and Szostak 1999; Jäschke 2001; Storz 2002). In recent years, numerous target-binding nucleic acid molecules (aptamers) have been identified; targets include organic molecules, antibiotics, proteins, and whole viruses (Wilson and Szostak 1999; Hermann and Patel 2000). In addition, in vitro selection experiments have led to novel RNA enzymes (ribozymes) and have ramifications for biomolecular engineering, for example, the design of allosteric ribozymes and biosensors (Soukup and Breaker 1999a,b, 2000) and aptamers for functional genomics (Famulok and Verma 2002). Other emerging applications of engineered RNAs include RNA synthetic biology, where designed RNAs are used to control cellular functions (e.g., regulate gene expression) (Isaacs et al. 2006). These exciting advances offer new investigative and application tools for molecular biology, proteomics, molecular medicine, and diagnostics (Breaker 2004). As applications of in vitro selection technology expand, the demands for efficient selection of complex RNA motifs increase in importance.

Many RNAs identified from random pools have simple structural motifs (e.g., stem–loop, stem–bulge–stem–loop) (Lee et al. 2004). Indeed, our graph-based analysis of random pools (25–100 nucleotides [nt]) showed that different RNA secondary topologies are far from uniformly distributed, with low yields for multiply branched structures, although complex structures gradually become more frequent as RNA length increases (Gevertz et al. 2005). Interestingly, recent experimental findings suggest that enhancing the structural diversity of RNA pools increases the possibility of obtaining novel RNAs with high activity (Carothers et al. 2004, 2006). Specifically, GTP aptamers with high-binding affinities are found to be more complex structurally than low-binding-affinity aptamers. The principal reason for the lack of structural diversity in random pools is due to incomplete and inefficient sampling of the astronomical size of the sequence space; random sequence sampling is inefficient because structures are not uniformly distributed in sequence space (Gevertz et al. 2005). To overcome this problem, heuristic approaches have been used to enhance the structural diversity of RNA pools. For example, structured pools have been synthesized by maintaining a constant stem–loop (GTP aptamer selection) (Davis and Szostak 2002) and by introducing random segments in existing RNA structures (e.g., purine nucleotide synthase and domains of group I ribozymes) (Jaeger et al. 1999; Ohuchi et al. 2002, 2004; Lau et al. 2004; Yoshioka et al. 2004). Recent works have also investigated the effects of sequence length (Legiewicz et al. 2006) and nucleotide composition (Knight et al. 2005) on recovery of specific RNAs. In addition, different functional classes of single-stranded RNAs have been found to have similar nucleotide compositions, implying evolutionary convergence (Schultes et al. 1997). The success of heuristic approaches depends on the details of the introduced sequence biases and the RNA function targeted. It is thus a challenge to develop systematic pool design approaches based on deeper understanding of pool sequence and structural complexity for the discovery of novel and complex RNAs.

To enhance in vitro selection experiments, RNA pools must possess sufficient sequence and structural complexity to ensure that the target RNA property exists in the pool. Given that we know little about the distribution of active RNAs in sequence and structural space (Carothers et al. 2004), an important goal of pool design is to maximize sequence and structural diversity without synthesizing all possible sequences. Even if complete coverage of sequence space is possible, not all regions of the space are likely to be productive for finding novel RNAs. This was suggested by recent analysis showing that the properties of GTP aptamers are correlated with their sequence/structural information content (Carothers et al. 2004). Unlike sequence space, the complexity of RNA structure space is more difficult to characterize quantitatively. At the secondary structural level, structural distributions of RNA pools can be analyzed using graph theory (Gan et al. 2003; Kim et al. 2004). Such an analysis shows that random pools are not structurally diverse (Gevertz et al. 2005), suggesting that pool structural diversity depends on how the sequence space is sampled. Indeed, understanding the relationship between sequence and secondary/tertiary structure spaces is essential for the design of effective pools for in vitro selection of RNAs. Thus, developing methodologies for generating and analyzing sequence pools possessing diverse RNA sequences and structures could enhance in vitro selection technology. Ultimately, a deeper understanding of the distribution of active RNAs in sequence and structural space will emerge through productive interactions between theoretical analysis and experiment.

Here we develop a computational approach for improving pool sequence and structural diversity by sampling sequences representing diverse regions of sequence space. We show that effective sampling of sequence space regions can be performed using nucleotide base "mixing matrices" for nucleotide transition rates applied to chosen starting sequences. Mixing matrices applied to given sequences are essentially generators of sequence pools and can be used to guide the reactants during in vitro selection experiments. Since we show that different regions of the sequence space are associated with distinct structural distributions, designed pools with specified target secondary structures can be obtained by optimizing a set of mixing matrices and starting sequences to approximate the target structural distributions. Figure 1 illustrates the relations among pool sequence/structure analysis, mixing matrix and starting sequence, and pool synthesis.

View larger version (54K): [in this window] [in a new window]   FIGURE 1. Modeling the RNA pool generation process using mixing matrices and analysis of pool structural distributions using tree graphs. The mixing matrix applied to any starting sequence specifies the mutation rates for all nucleotide bases. The matrix elements of each row represent nucleotide base (A, C, G, U) composition in a vial or synthesis port. Mixing matrices and starting sequences can be optimized to yield target structured pools.

	MATERIALS AND METHODS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION APPENDIX ACKNOWLEDGMENTS REFERENCES

 
Defining mixing matrices for generation of nonrandom sequence pools
Understanding pool synthesis strategies provides clues for improving in vitro selection technology. The standard experimental protocol involves synthesizing DNA sequences and then transcribing them into RNA sequences by RNA polymerases. Current sequence synthesis strategies include chemical synthesis of short sequences (<150 nt); enzymatic assembly of short strands; synthesis of designed sequences with constant (double-stranded) and variable (single-stranded) regions; and synthesis of sequences around a designed sequence by random and biased mutations (Wilson and Szostak 1999). Random pools with short sequences are normally synthesized chemically, whereas longer sequence pools can be assembled enzymatically from shorter, chemically synthesized strands using techniques such as ligation or template-directed polymerization (Jäschke 2001; Stuhlmann and Jäschke 2002).

Our strategy is to define substrate pools with both random and biased sequence mutations around specific (starting) sequences to generate designed sequences with constant and variable regions. For this purpose, we introduce the mixing matrix M, whose elements specify mixing (or "contamination" or "doped") in the four phosphoramidite (A, C, G, and U [i.e., T in DNA synthesizer]) vials; applying mixing matrices to starting sequences leads to designed sequence pools. (Note that since successive nucleotide bases are selected independently from the vials, sequence synthesis methods correspond to the zeroth-order Markov process.) This representation of pool generation or synthesis using mixing matrices enables computational analysis of pool sequence and structural diversity. Our goal is to define, via computational analysis, an optimal set of starting sequences, mixing matrices, and associated weights for a given target structural distribution in the pool. Figure 1 illustrates pool design and synthesis via mixing matrix and analysis of sequence/structure space.

For pool synthesis using four vials or ports, the corresponding mixing matrix M is a 4x4 matrix that specifies the molar fractions of nucleotide components A, C, G, and U (T) in the four vials. Thus, the "ij" element of M (i.e., M _ij ) denotes the molar fraction of base j in the vial "for base i." It describes how we can dope that vial for i by introducing other bases ji into it as well. The design problem involves selecting those doping ratios and starting sequences. For example, M _AU is the fraction of U (T) nucleotides in the vial for A, M _AA is the fraction of A in the vial A, and M _UA is the fraction of A in the vial U (T). Thus, the elements of each row of the matrix sum to unity:

(1)

If the DNA synthesizer is to produce a fixed sequence, then a vial for base i has 100% base i and zero fraction of other bases (i.e., M _ii =1 and M _ij =0 for ij). If M _ii <1 and M _ij 0, contaminations are introduced, as specified by the off-diagonal elements of M. The expected number of mutations in a synthesized sequence is determined by

(2)

where N_j is the number of nucleotides of type j in the original sequence.

Ideally, we would like to determine the mixing matrices M for a target structural distribution or on the basis of specific biological-motivated contamination protocols. In practice, the inverse design problem—specify M and analyze the resulting structural distribution—is much easier to perform. We thus construct different mixing matrices motivated by biological covariance mutations and analyze their coverage of sequence space via a standard clustering method. Direct modeling of pools generated by a specific mixing matrix can be made by exploiting correlations between bases in folded RNAs. For example, the bases in paired and unpaired regions are correlated, allowing assignment of matrix elements for mutating bases in such regions.

Our biological motivation for choosing the mixing matrix classes is as follows. We consider mixing matrices with symmetric elements, M _AU=M _UA, M _CG=M _GC, M _GU=M _UG, to preserve base pairs. Such matrices cover the sequence subspace approximating covariance mutations (e.g., AU to UA, CG to GC, GC to UA). Covariance mutations have been used to analyze the secondary structure and sequence consensus of RNA sequence families. For example, this approach has been successfully applied to search for tRNA-related sequences and other small RNAs (Eddy and Durbin 1994). Alternatively, to disrupt stems and generate new structures, we can consider mixing matrices that do not preserve base pairs. Such matrices include asymmetric matrices without the property of covariance mutations. Noncovariance mutations, including random mutations, are commonly used to generate sequence pools for in vitro selection applications.

To sample the sequence space, we define five classes of mixing matrices motivated by biological considerations, based primarily on sequence transformations associated with covariance mutations. The mixing matrix classes are characterized by the following matrix elements: (A) varying diagonal elements M _ii with the condition M _AA=M _CC=M _GG=M _UU; (B) M _CC=M _GG=1; (C) M _AA=M _UU=1; (D) M _AC=M _UG=1; and (E) M _CA=M _GU=1. Within each class, several mixing matrices are constructed whose elements are distributed uniformly in steps of 0.25. A total of 22 mixing matrices representing the five classes are displayed in Figure 2. The matrix classes to which they belong are as follows: class A matrices 1–6, class B matrices 7–10, class C matrices 11–14, class D matrices 15–18, and class E matrices 19–22. Note that in vitro experiments effectively use random pools generated by a constant 4x4 mixing matrix, where all 16 elements are 0.25; this corresponds to our matrix 4.

View larger version (50K): [in this window] [in a new window]   FIGURE 2. Our five classes of 22 mixing matrices (MM) for generating diverse sequence pools. The matrix classes are developed based on alteration of diagonal elements (class A) and covariance mutations (classes B–E). For pool synthesis using four vials, the mixing matrix is a 4x4 matrix specifying the molar fractions of nucleotide components A, C, G, and U in the four vials. The columns represent the molar fraction of the four bases in vial for each base denoted in each row.

In RAG, RNA graphs are organized into n-vertex families, and members of a family are ordered using a topological index (i.e., Laplacian eigenvalues) (Fera et al. 2004; Gan et al. 2004). Structural complexity can be measured by the graph's vertex number (V) and the second smallest Laplacian eigenvalue (₂). For example, a linear chain has a smaller eigenvalue than a branched structure. The number of motifs in an n-vertex family increases with n, the number of vertices. For example, the 6-vertex tree family has six distinct trees and the 7-vertex family has 11 trees (see RAG Web resource at http://monod.biomath.nyu.edu/rna). For easy reference, each tree motif is labeled by vertex number and ordering within the family; for example, members of the 6-vertex family are labeled 6₁, 6₂, ..., 6₆. Since vertex number (V) is related to RNA sequence length L, these are constant length pools; in fact, we found empirically that L=20(V–1) (Gan et al. 2003). Our pool design will focus on tree structures because RNA folding algorithms for tree structures are efficient; computationally demanding pseudoknot folding algorithms are also available (Rivas and Eddy 1999; Ren et al. 2005). RNA graph theory also provides a complete set of pseudoknot and nonpseudoknot motifs for more general assessment of pool structural diversity.

Starting sequences for pool generation
The six starting sequences with distinct tree structures (Fig. 3) are 70S (chain F) (80 nt), tRNA (81 nt), P5abc domain of group I intron (56 nt), GTP-binding aptamer (69 nt), modified P5abc domain (51 nt), and modified GTP-binding aptamer (54 nt). As shown in Figure 3, distinct tree structures are represented as graphs by converting stems to edges and other structural elements (e.g., loop, bulge, etc.) to vertices according to tree graph rules developed previously (Gan et al. 2003). As shown in Figure 3, the Laplacian eigenvalue (₂) indicates the structural complexity of starting sequences. Generally, the starting structure allows exploration of the structural neighbors of that structure via mutations. For random mutation rates (constant matrix elements of 1/4), the generated pools have no memory of the starting sequence. We generate pools with all possible combinations of 22 mixing matrices and six starting sequences for pool structured designs.

View larger version (27K): [in this window] [in a new window]   FIGURE 3. Starting sequences and their secondary structures for pool synthesis using mixing matrices. Displayed are the secondary structures and corresponding tree graphs for four existing and two modified existing RNAs. Laplacian eigenvalue (2) of the tree graph indicates the structural complexity.

(1)

where the four-component vector X _i , i=1, 2, ..., n, represents the nucleotide base:

(2)

The matrix Y represents the sequence pool generated by M with starting sequence S. For example, if

(3)

and

(4)

then Y is given by

(5)

The probability of finding a new sequence S' in the pool can be calculated from Y. For example, P(S'=ACU)=Y ₁₁•Y ₂₂•Y ₄₃=0.1•0.1•0.7 and P(S'=GUC)=Y ₃₁•Y ₄₂•Y ₂₃=0.1•0.1•0.1. Similarly, we can calculate the frequency of sequences with a specified base-pairing scheme using this mathematical formulation. However, a rigorous mapping of sequence space (Y) to shape space (possible RNA graphs) requires an RNA folding algorithm, as described in our previous work (Gan et al. 2003).

The challenge in computational pool design is to find an optimal set of mixing matrix (M), starting sequence (S), and weight (or pool fraction) for generation of target-structured pools. In principle, the mixing matrix can be calculated using statistical thermodynamics from the distribution of shapes in the designed pool. Assuming that the designed pool consists of N noninteracting RNA molecules, the probability of finding topology t in the pool is

(3)

where E(S_i ) is the energy of sequence S_i , =1/kT, (E_i ) is the density of states, and S _i, =1/kT, (E _i) is an RNA topology operator defining tree or pseudoknot shapes enumerated by RNA graph theory. Equation (3) defines the relation between sequence pool {S_i } and target structural distribution P(t). Recently, we calculated P(t) distributions for 25–100 nt random pools using a folding algorithm and a program for converting secondary structures into tree graphs (Gevertz et al. 2005). Thus, the goal is to determine the sequence pool {S_i }, or mixing matrices generating that pool, to produce the target distribution P(t). In the Appendix, we describe a practical protocol for finding optimal mixing matrices approximating the target P(t) based on analyses of sequence space and pool structural distribution. Alternative pool design methods may also be developed based on Equation (3).

Pool sizes
For practical reasons, our computations used relatively small pools of 10,000 sequences. To show the effect of pool size, Figure 4 plots the frequency of several tree motifs (4₁, 4₂, 5₁, 5₂, 5₃, and 6₁) for pool sizes of 5,000–60,000 sequences using mixing matrix 4 (MM4) and the initial tRNA sequence. We see that the pool fractions for distinct tree motifs saturate beyond 5000 sequences, indicating that the error due to sample size is small. The rapid saturation of pool fraction stems from mapping secondary structures using simple graphs. If detailed motif features (size of loops, stems, etc.) are incorporated into the mapping, larger pool sizes will certainly be required.

View larger version (15K): [in this window] [in a new window]   FIGURE 4. Effect of pool size on pool fractions of selected tree motifs. The pools are generated using random mixing matrix 4 and starting tRNA sequence in Figure 3B.

	RESULTS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION APPENDIX ACKNOWLEDGMENTS REFERENCES

 
Coverage of sequence space regions generated by mixing matrix classes and starting sequences is distinct from random pools
We consider our set of starting sequences (Fig. 3) and 22 mixing matrices (Fig. 2) to explore the sequence/structure space of sequence pools and to optimize target pools.

To analyze the clustering patterns in sequence space, we cluster all sequences generated by the mixing matrices using a standard clustering technique (e.g., MDS), allowing visualization of sequence similarity/dissimilarity properties (Cox and Cox 1994). A similar procedure is commonly used for investigating the diversity of chemical compound libraries (Xie et al. 2000). Such analysis helps establish the relation between each mixing matrix and the generated sequence space. Given a pool of sequences, we define Hamming distances (number of dissimilar bases) between all pairs of sequences (see Materials and Methods), allowing data projections in 2D, 3D, and higher dimensions.

Figure 5, A and B, shows the 2D and 3D clustering of sequences generated by 22 mixing matrices using starting sequences for the modified P5abc (Fig. 3E) and 70S (Fig. 3A) RNAs, respectively. In Figure 5A, we see that the sequences generated by the five mixing matrix classes and the P5abc starting sequence cover distinct regions of the sequence space, especially the boundary and central regions. The boundaries are spanned by matrix classes B–E, and the central region by matrix class A. Intriguingly, the random MM4 produces sequences that are localized in sequence space, showing that the standard approach does not provide an efficient sampling of diverse regions of sequence space in agreement with observations. More adequate sampling of sequence space is provided by the nonrandom mixing matrices. The 70S starting sequence yields similar global 2D and 3D clustering patterns: the five matrix classes yield clusters in distinct sequence regions (Fig. 5B). Although the 22 mixing matrices of Figure 2 provide a comprehensive coverage of the sequence space, some regions remain sparsely populated, indicating that the matrix classes must be expanded for more complete coverage. Still, the chosen matrix set is diverse enough for initial assessment of our structured pool design concept.

View larger version (49K): [in this window] [in a new window]   FIGURE 5. Two- and three-dimensional clustering plots using the MDS transformation for sequences generated from 22 mixing matrices (labeled 1–22) (A) starting with a modified P5abc domain (Fig. 3E) and (B) with 70S (Chain F) (Fig. 3A). The distance between any sequence pair is the Hamming distance, a measure of the number of dissimilar nucleotide bases. Axes represent two or three largest components of the projection. Each color represents a sequence pool generated by one of the 22 mixing matrices; the "x" mark on the left represents result for an invariant sequence transformation corresponding to diagonal matrix M ii =1. The mixing matrices are grouped into five classes (A–E) according to their matrix properties (Table 1).

View larger version (21K): [in this window] [in a new window]   FIGURE 6. Pool fraction distributions for six tree motifs in pools generated from 22 mixing matrices (labeled 1–22), starting with a modified P5abc domain (Fig. 3E), which has a 42 tree motif. The results for the random pool 4 (marked with arrow) are displayed as dot-filled histograms.

Second, we note a pattern in the correlation between 4₂ and 5₂ trees and between 3₁ and 4₁ trees. With the Figure 3E starting sequence, Figure 6 and Table 2 show that sequence pools from matrices 7–13 have a large proportion of 4₂ and 5₂ tree structures compared with the random pool 4. Similarly pools from matrices 14–22 possess >30% 3₁ and 4₁ tree structures. This pattern emerges because the 4₂/5₂ and 3₁/4₁ tree-motif pairs are related by an internal loop or bulge, which can be created by a few mismatched base pairs.

Third, the structural distributions generated by a tRNA sequence (5₃ tree motif; Fig. 3B) differ from those for the modified P5abc domain (4₂ tree motif; Fig. 3E) in one important respect (Table 2). The most likely motifs are the simpler 5₁ and 5₂ trees rather than the starting 5₃ motif. MM1, for example, generates only 5% 5₃ motif, but 26% 5₂ tree. In contrast, MM7, which preserves C and G bases, generates 23% 5₃ trees, while other combinations of mixing matrices and starting structures yield almost no 5₃ trees (Table 2). The mean mutation rate for MM7 from the starting tRNA sequence (Fig. 3B) is 0.1 (8.5 positions among 81 nt). Thus, the 5₃ motif populations produced by matrices 1 and 7–11 are much higher than in random pools (1.3%). The difficulty of generating significant populations with the tRNA-like 5₃ tree motif likely stems from the lower thermodynamic stability of 5₃ compared to 5₁ and 5₂ trees. Our analysis shows that matrices 7–9 generate sequences that are favorable for stabilizing the 5₃ motif because these matrices preserve energetically favorable CG base pairs.

To increase the population of complex folds like the tRNA-like 5₃ tree motif, we consider refining the mixing matrices 7–9. Since class B matrices 7–9 produce a higher frequency of 5₃ tree, we search for matrices in the neighborhood of this class by exhaustively varying the elements in each row with M _ij =0.2, yielding 56 possible cases. Assuming that each row is independent, the total number of mixing matrices around the class B matrix region is 56² or 3136, since two rows (second and third) are identical and the other two rows (first and last) have a total of 56 cases each. We filter the 3136 trial mixing matrices, yielding better than a 23% tRNA-like 5₃ tree structure. Remarkably, 12 of the 3136 mixing matrices for tRNA-like topology fulfill our requirement forming 5₃ motifs. We use these "MMT" matrices (Fig. 7) to generate graph-structural distributions with tRNA shapes, as shown in Table 3. For example, MMT6 generates a 51% tRNA-like 5₃ tree motif with 15 mutations out of 81 bases.

View larger version (32K): [in this window] [in a new window]   FIGURE 7. Twelve refined or variants of class B mixing matrices for enhancing pools with the tRNA-like (53) structure (MMT); they generate pools with at least 23% 53 tree motif.

Sequence/structure correlations exist in designed pools
The above survey of tree structural distributions provides an analysis of RNA shapes in designed pools. We now analyze sequence/structure distributions at the nucleotide base level generated by the 22 mixing matrices starting with a 51-nt P5abc domain (Fig. 3E). In Figure 8, we use sequence Hamming and tree edit distances to quantify sequence and structure distances, respectively, as defined in Materials and Methods. Recall that the Hamming distance is the number of differing letters between two equal-length RNA sequences aligned end to end (Hamming 1987). The tree edit distance between two (full) tree secondary structures measures the minimum sum of the cost (insertion, deletion, and replacement of nodes) along an edit path for converting one tree into another (Hofacker 2003).

View larger version (25K): [in this window] [in a new window]   FIGURE 8. Contour plots of sequence/structure relationships using Hamming distance versus tree edit distance for pools generated by 22 mixing matrices, starting from a modified P5abc domain (Fig. 3E). Note that the X and the Y axes are always 0–100 and 0–60, respectively, and that each intensity bar indicates the frequency of joint distance distributions (the frequency outside the box is 0). There are 10,000 sequences in each pool.

As the mixing matrices are altered, the distribution of the tree edit distances also changes. Generally, tree edit distance increases with mutation rate. For example, Figure 8 shows that tree edit distances become larger as diagonal elements of matrix classes A, B, and C decrease. Changing C and G bases (matrix class C) has a larger effect on the starting structure than changing A and U bases (matrix class B), as evident from the pool distances from the origin. This is due to lower free energies associated with GC base pairs compared to AU base pairs. Thus, Figure 8 indicates that the distribution of sequence/structure distances from the initial sequence is controlled by the elements of the mixing matrices. Although the patterns of sequence/structure distributions are not sensitive to the starting sequences (data not shown), the densities within the localized regions are markedly changed. Figure 8 shows that the pools generated by most mixing matrices (except for 1, 8, and 9) and starting sequences of a modified P5abc domain produce a single cluster. We find that contour plots with string edit distance or base-pair distance (data not shown) show somewhat less information about pool structural properties than those with tree edit distance. Other secondary structure measures may also be used to capture structural differences among folds in the same vertex or tree class. For example, it is informative to know the distribution of stem, loop, and bulge sizes (Fontana et al. 1993).

Parameter optimization can lead to design of structured RNA pools
The preceding analysis of sequence space and assessment of structural distributions generated by nonrandom mixing matrices allow design of target structured pools. Here we use the pool design algorithm (Appendix) to develop several structured pools by selecting an optimal combination of starting sequences, mixing matrices, and associated weights {(S_i , M _i , _i )}. The best combination for a target pool is dictated by the frequency data (Fig. 6; Tables 2, 3).

To illustrate, Table 4 shows four examples of designed pools that are rich in specific tree structures (e.g., 4₁, 5₁, 5₂); also displayed are their pool characteristics (mixing matrix weights and tree motif frequencies). Specifically, our target pools are: Pool T_A with 4₁ and 4₂ structures; Pool T_B with 5₁, 5₂, and 5₃ structures; Pool T_C with 4₂, 5₂, and 5₃ structures; and Pool T_D with 4₁, 4₂, and 5₃ structures. Pools T_A and T_B are 4- and 5-vertex pools, respectively, and Pools T_C and T_D are pools with mixed n-vertex structures. Each designed pool represents an optimal combination of starting sequences, mixing matrices, and associated weights derived using Step 5 of our design algorithm (see Appendix). Briefly, we initially choose pool fractions T ₁, T ₂, ..., T_n for target motifs and the number of mixing matrices to approximate the target pool. We then use Equation (6) in the Appendix to calculate the weight ₁, which depends on T ₁, starting sequence S ₁, and mixing matrix M ₁. Next, we minimize the error between the target and estimated target pool fractions, Equation (8) in the Appendix, over all pools generated by starting sequence/mixing matrix pairs {(S_i , M _i )}. This procedure yields optimized starting sequences, mixing matrices, and weights; the mean mutation rate is calculated based on these sequences, mixing matrices, and their weights.

Pool T_B, targeting 20% each of the 5₁, 5₂, and 5₃ structures, is generated from the Figure 3A sequence for 5₁ with matrix 13 at a weight of 0.18 and the Figure 3B sequence for 5₃ with matrix T12 at a weight of 0.82. For this pool optimization, we expanded our mixing matrix/starting sequence repertoire to include those in Table 3 (the 12 mixing matrices for generating pools with a high frequency of 5₃ motifs, which are extremely rare in random pools). Thus, this optimization was performed over the set of 144 (22x6+12x1) mixing matrix/starting sequence pairs.

Resulting Pool T_B contains 20% each 5₁, 5₂, and 5₃ tree motifs, compared with 23%, 16%, and 0% for the random pool (MM4), matching the target exactly. We found that using the 12 MMT matrices dramatically increases the population of the 5₃ motif (at a cost of decrease of 3₁, 4₁, and 4₂ motifs). The 6₁ structure (9%) is the next highest species in Pool T_B.

Target Pools T_C and T_D are mixed pools with both 4- and 5-vertex tree structures, designed from our 144 mixing matrix/starting sequence pairs. The targets for Pool T_C are 4₂, 5₂, and 5₃ tree motifs, and those for Pool T_D are 4₁, 4₂, and 5₃ tree motifs (20% for each). Pool T_C is generated by the Figure 3E sequence (MM9, 0.60) and the Figure 3B sequence (MMT2, 0.40), and Pool T_D is produced by the Figure 3B sequence (MMT6, 0.329) and the Figure 3F sequence (MM13, 0.608). The results are as expected: Pool T_C has frequencies for 4₂, 5₂, and 5₃ motifs of between 19% and 20%; Pool T_D has frequencies for 4₁, 4₂, and 5₃ trees of 17%, 20%, and 20%, respectively, all within 3% of the target.

Our designed pools above, involving three of the 22 mixing matrices and two of the 12 MMT matrices, only touched the surface of possibilities. Still, in practice, it might be preferable to approximate a target pool using a small number of mixing matrices. Once our algorithm is automated (Appendix), exploration of pool design can be routinely performed.

A designed pool improves the selection of GTP aptamers
We now apply our pool design approach for enhancing GTP-binding aptamers. Szostak's group recently found that the GTP aptamer's binding affinity is correlated with the informational complexity (Carothers et al. 2004, 2006). Informational complexity is correlated with structural complexity (e.g., number of stems, vertex number of tree graph). As the information content and binding affinity decrease (Carothers et al. 2004, see their Fig. 1, panels A and B), the aptamers have simple structures such as 2₁ or 3₁ tree motifs. Specifically, a high-affinity GTP aptamer with high informational complexity (Carothers et al. 2004, see their Fig.1, panel C) has the 5₂ tree structure (Fig. 3D). Interestingly, no GTP aptamer with a 4₂ tree structure (Fig. 3F) has been reported, although it is structurally similar to the 5₂ tree. Because the frequency of the 4₂ motif is only 12% in the random Pool T_F (Table 4), we propose designing a GTP aptamer pool by enriching the pool with 5₂ and 4₂ motifs. Our target pool fractions (T_i ) are 20% for 4₂ and 26% for 5₂. Our optimization yields Pool T_E (Table 4) as a combination of two subpools: the Figure 3D sequence (MM13, 0.625) and the Figure 3F sequence (MM10, 0.375). The frequencies of 4₂ and 5₂ trees in the designed pool are 21% and 26%, respectively, nearly as desired and very different for the 12% and 16% distributions of these motifs in the random Pool T_F. The sequence/structure contour plots in Figure 9 show differences between the designed and random pools; the designed pool has a relatively high mean mutation rate of 0.349.

View larger version (27K): [in this window] [in a new window]   FIGURE 9. Comparison of designed GTP (upper) and random (lower) pools using contour plots of Hamming distance versus tree edit distance. The GTP pool is generated by 62.5% MM13 starting with the 52 motif (Fig. 3D) and 37.5% MM10 starting with the 42 motif (Fig. 3F). The random pool is generated using the starting sequence in Figure 3D.

	DISCUSSION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION APPENDIX ACKNOWLEDGMENTS REFERENCES

Thus, our work suggests that designing pools for enhancing in vitro selection can follow several research avenues. Maximizing sequence and structural diversity broadly can increase the probability of finding a given RNA property using nonrandom mixing matrix/starting sequence pairs. An advantage of this approach is that designed pools can be directly implemented in pool synthesis. Alternatively, we can target a specific structural distribution by determining optimal mixing matrices and starting sequences without explicit sequence/structure mapping. Our targeted pool design can be applied to known structures, novel motifs (Kim et al. 2004), complete sets like n-vertex pools, or perhaps submotifs of RNAs (Zorn et al. 2004). Of course, a more comprehensive set of mixing matrices covering wider regions of the sequence space should be sought systematically. For example, matrices conserving noncanonical base pairs (AC, CA, GA, AG, etc.) can complement our current set, which conserves canonical base pairs; there are 12 such classes from a total of 16 possible base pairs.

Another design theme is enrichment of pools with structures resembling a target-active molecule. We illustrated this approach using GTP-binding aptamers. The conventional approach – designing pools in the sequence neighborhood of a target molecule (Lau et al. 2004; Ohuchi et al. 2004; Yoshioka et al. 2004), however, does not ensure that the designed pools will cover the structural neighbors of the target molecule, unless sequence mutations are made to localized sequence segments, as is commonly done in many experiments. In contrast, our optimized pool design approach (Appendix), allows enrichment of pools with specific RNA topologies or structures (e.g., tRNA-like 5₃ tree). In addition, novel tree topologies in the neighborhood of the target molecule, as suggested by structural enumeration (Kim et al. 2004), could be similarly engineered.

Clearly, further developments of sequence/structure analysis techniques are needed to improve the pool design and overcome specific limitations. Understanding the sequence/structure relationship is one of the most fundamental biological problems not only for RNA but also for proteins. In our analysis, we are limited by the usage of numerical secondary structure folding algorithms, which are still imperfect and inefficient for predicting pseudoknot structures. However, our risk has been reduced here by "folding" of small RNAs (<100 nt) only and focusing on statistical properties (e.g., frequencies of topologies). A general strategy for improving structure prediction is to consider many suboptimal structures using, for example, the Boltzmann sampling method (Ding and Lawrence 2003).

Ultimately, RNA tertiary and higher-order folding is essential to understand RNA function. Perhaps progress on this problem will be realized in the near future. For now, we offer our sequence pools possessing diverse RNA secondary structures as an approach to enhance in vitro selection technology.

Our pool design algorithm can be fully automated given target RNA shapes (and possibly starting sequences). We are developing a publicly available Web server to allow experimentation of pool design and analysis of RNA pool properties (e.g., base composition, size distribution of stems, bulges, etc.), and to define optimal mixing matrices for pool synthesis. Experimental synthesis of designed pools (specific structural motifs and their frequency) can be performed by using optimized starting sequences, mixing matrices, and associated weights. When available, location of this server will be noted on our group Web site (http://monod.biomath.nyu.edu). We hope that this tool will help stimulate the productive interaction between theoretical and experimental efforts.

	APPENDIX

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION APPENDIX ACKNOWLEDGMENTS REFERENCES

 
An algorithm for designing structured RNA pools
Our pool design algorithm is based on analyses of sequence and structure spaces to allow design of specific structures, including novel RNA-like motifs identified using graph theory analysis (Kim et al. 2004). The algorithm below assumes that we have available reference data such as shown in Tables 2 and 3 that relate mixing matrices and starting sequences to motif distributions in resulting pools. The sequence space regions are mapped via various mixing matrices using a standard clustering method; the structural distribution is computed by converting secondary structures into tree graphs. By knowing the structural distributions of various sequence space regions, we then can optimize the choice of starting sequences and mixing matrices to approximate the target structured pool for future work.

Our pool design algorithm involves the following steps:

1. Specify a target distribution of topologies/shapes.
   
2. Define candidates for starting sequences and mixing matrices that aim to cover the sequence space. The mixing matrices have been constructed, for example, based on covariance mutations. The mixing matrices and starting sequences may remain the same for different structured pool designs. We "visualize" the diversity of a set of RNA sequences using a standard sequence similarity/dissimilarity clustering based on Hamming distance (number of dissimilar bases) between any pair of aligned sequences. In this study, we used mainly six starting sequences and constructed 22 mixing matrices to cover the sequence space (see Results).
   
3. Compute shape frequency distributions corresponding to all starting sequence/mixing matrix pairs, as discussed below and detailed in our previous study (Gevertz et al. 2005). This step analyzes pool structural diversity.
   
4. Choose the number of mixing matrices to approximate the designed pool.
   
5. Find an optimal combination of starting sequences (S_i ) and mixing matrices (M _i ) and associated weights ( _i ) to approximate the target RNA shape distribution. The mathematical procedures for this step are detailed below.
   

The designed pool is composed of k smaller subpools defined by the set {(S_i , M _i , _i )}, i=1, 2, ..., k. The above pool design algorithm can be fully automated given target RNA shapes (and possibly starting sequences). We are planning to make publicly available a Web server to allow experimentation of pool design and analysis of RNA pool properties, and to obtain mixing matrices for pool synthesis. Experimental synthesis of designed pools can be performed by using trial S_i , M _i , and _i .

In Step 3, the pool structural distribution is calculated by mapping RNA secondary structures into graph space. This is done by predicting secondary structures of all sequences using the Vienna RNAfold package and then converting them into tree graphs, as described elsewhere (Gevertz et al. 2005). It is known that 73% of known base pairs are predicted by free-energy minimization algorithms such as RNAfold for sequences with <700 nt (Mathews and Turner 2006). For greater accuracy, the Boltzmann sampling method can be used to generate a set of 1000 suboptimal structures (Ding and Lawrence 2003), although at a higher computational cost (1000 times pool size). Specifically, base-pairing information in the .ct file generated by the RNAfold program is used to convert a secondary fold into a tree graph. The topologies of the folds are determined using Laplacian eigenvalues of tree graphs as implemented in our RNA Matrix Program (Gan et al. 2004) (server available at http://monod.biomath.nyu.edu/rna). Specifying tree topologies using eigenvalues is inexact because different topologies can have the same spectrum; the assignment error rate is a few percent for small tree topologies (<10 vertices). This step is similar to the RNAshapes program, which uses bracket notations for representing secondary structures (Giegerich et al. 2004; Steffen et al. 2006). Unless stated otherwise, each sequence pool has 10,000 sequences, which is adequate for assessing structural distributions using simple tree graphs. Structure prediction and conversion to tree graphs for 10,000 80-nt sequences require 1 h on an SGI 300 MHz MIPS R12000 IP27 processor.

In Step 5, we approximate a target structural distribution by optimizing a set of starting sequence/mixing matrix pairs based on pool structural frequency data. Generally, we consider a designed pool composed of k subpools, each generated with a mixing matrix/starting sequence pair and associated with a weight _i : p(S ₁, M ₁, ₁), p(S ₂, M ₂, ₂), ..., p(S_k , M _k , _k ), where ₁+₂+...+ _k =1 and p(S_i , M _i , _i ) denotes synthesizing the _i fraction of the pool sequences using starting sequence S_i and mixing matrix M _i . Optimization of the three pool parameters S_i , M _i , and _i can be formulated as follows: If the nx1 matrix T is the target distribution with T_i fractions of structures 1, 2, ..., n and F_l (S_i , M _i ) is the pool fraction of structure l generated by starting sequence S_i and mixing matrix M _i in Tables 2 and 3, the pool parameters (S_i , M _i , _i ) can be optimized by the following equation:

(4)

where =(_{1, 2,} ..., _k ) subject to the conditions ₁+₂+...+ _k =1 and _i 0. Since experimental implementation of pool synthesis is simpler with fewer mixing matrices, we consider the solution of for k=2 below; the optimization procedure can be generalized. Formula (4) with only two mixing matrices M ₁ and M ₂ reduces to

(5)

The solution for the only weight is

(6)

The estimated pool fractions for the other shapes or topologies 2, 3, ..., n are derived from the known ₁, F ₁(S ₁, M ₁), and F ₁(S ₂, M ₂) as follows:

(7)

We then optimize (S ₁, M ₁), and (S ₂, M ₂) by minimizing the error

(6)

The above procedure will allow us to obtain the optimized parameters ₁, (S ₁, M ₁), and (S ₂, M ₂) for a target distribution T. The convergence of the procedure depends on the number of mixing matrices and starting sequences, or coverage of the sequence/structure space.

	ACKNOWLEDGMENTS

TOP ABSTRACT