Systems and methods for high throughput aptamer structure classification and assessment of aptamer structure functionality

Abstract

Systems and methods are disclosed for high throughput classification and functionality assessment of aptamer structures. A non-limiting exemplary method includes: receiving, by a processor, a plurality of aptamer sequences; generating, by the processor, a structure for each aptamer sequence of the plurality of aptamer sequences to determine a plurality of structures; identifying, for each structure, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures; classifying, by applying the set of parameters of each structure into a machine learning model, the plurality of structures into one or more aptamer clusters; and determining a target binding characteristic for each aptamer of the one or more aptamer clusters.

Description

SYSTEMS AND METHODS FOR HIGH THROUGHPUT APTAMER STRUCTURE CLASSIFICATION AND ASSESSMENT OF APTAMER STRUCTURE FUNCTIONALITYRELATED APPLICATION

[0001] This application claims priority to U. S. Provisional Patent Application Serial No. 63 / 717,735, filed on November 7, 2024, which is incorporated herein by reference in its entirety for all purposes.STATEMENT REGARDING FEDERALLY SPONSORED R& D

[0002] This invention was made with government support under grant / contract number 510776 awarded by The Simons Foundation; grant / contract number 2027277 awarded by the National Science Foundation; grant / contract number 2318817 awarded by the National Science Foundation; grant / contract number 2404470 awarded by the National Science Foundation; and grant / contract number R61MH135106 awarded by the National Institutes of Health - National Institute of Mental Health (NIH-NIMH). The government has certain rights in the invention.INCORPORATION OF THE SEQUENCE LISTING

[0003] This application contains a Sequence Listing, which is incorporated herein by reference in its entirety. The accompanying Sequence Listing xml file, named “UCLA. P0221WO_SL” was created on November 6, 2025 and is 24,581 bytes.BACKGROUND

[0004] The embodiments disclosed herein are generally directed towards high throughput systems, processes and methods that determine, classify and assess the functionality of aptamer structures.

[0005] Aptamers are oligonucleotide structures formed of oligonucleotide sequences, such as single stranded DNA, that are typically selected and used for binding to targets (e.g., ions, organic molecules, peptides, proteins, and cells) with high affinity. The aptamers’ small size, reproducible chemical synthesis, biocompatibility, low immunogenicity, and structural stability relative to the targets make aptamers particularly advantageous in biosensing, therapeutics, and other applications. The identification of effective oligonucleotide sequences for aptamer generation is thus crucial to the formation of effective aptamers. Aptamer candidate sequencesare typically obtained from a SELEX screening process and sequenced using next-generation sequencing (NGS), after which they are organized hierarchically based on the number of hits (i.e., the number of times each sequence occurs in the NGS output). Aptamer candidate sequences typically form secondary structures (e.g., 2 dimensional structures) and tertiary structures (e.g., 3 dimensional structures) based on bonds linking one nodes (e.g., nucleotides) of one or more regions of the aptamer candidate sequence to other nodes, often of other regions. When an aptamer binds its target, conformational rearrangements facilitate energetically favorable interactions with the target and surrounding solution ions. This results in global changes in aptamer secondary and tertiary structure intramolecular rearrangements. Modeling and studying the secondary and tertiary structures formed by aptamer candidate sequences is therefore crucial to the selection of effective aptamer candidate sequences and production of effective aptamers.

[0006] Although conventional NGS and aptamer selection techniques can yield large numbers (e.g., thousands) of candidates, very few candidates (e.g., a few dozen, typically having the highest hits) are advanced for further analysis of target binding properties and suitability for application use. Thus, the range of options for aptamer candidate sequences is low and experimentally intractable. Therefore, numerous potentially useful, even exemplary candidates are left unexplored. In particular, the current range of options for aptamer candidate sequences are insufficient for producing aptamers suitable for targeting small molecules. Selecting aptamers for targeting small molecules is notoriously difficult due to the limited number of functional groups, affording less opportunity for inciting specificity via aptamer-target hydrogen bonds, electrostatic interactions, and hydrophobic interactions.

[0007] Conventional techniques for modeling and studying the secondary and tertiary structures formed by aptamer candidate sequences are not as effective, especially for the study of aptamers to target small molecules. For example, such techniques are restricted by the limited selection of aptamer candidate sequences, leaving out potentially useful aptamer candidate sequences for small molecule binding. Furthermore, conventional techniques often present inconsistent and conflicting results of predicted secondary or tertiary structure, and often require large amounts of computation and causing computational resource burden. The limited range and aforementioned defects of aptamer candidate structures also comprises subsequent analysis of the aptamer candidate structures (e.g., for binding characteristics). Furthermore, there is an inadequate integration of the analysis of the aptamer candidates with earlier sequencing and structure prediction steps, making aptamer research time consuming and computationally inefficient, particularly for high throughput analysis of aptamer candidates.

[0008] There is thus a desire and need for more comprehensive systems and methods for determining and modeling structures for aptamers that expand the possibilities of aptamer candidates, especially for targeting small molecules. Furthermore, there is a desire and need for such systems and methods to be high throughput, integrated, and computationally efficient in delivering insights (e.g., regarding binding characteristics) for a large number of aptamers.Various embodiments of the present disclosure address one or more of the above described shortcomings.SUMMARY

[0009] Recognized herein is a need for high-throughput, integrated, and computationally efficient systems and methods for the classification and functionality assessment of aptamer structures, particularly in the context of targeting small molecules. Addressing shortcomings of the conventional approaches, the present disclosure introduces a machine learning-driven approach that processes large volumes of aptamer sequences to generate their structural models, assess inter- structural similarities, classify aptamers into functional clusters, and determine target binding characteristics, thereby expanding the pool of viable aptamer candidates and improving the precision and scalability of aptamer discovery workflows.

[0010] The present disclosure relates to relates to a high-throughput computational system and method for analyzing, classifying, and selecting aptamer sequences based on their predicted secondary structures. The presently described system efficiently generates structural models of aptamer candidates and represents these structures using mathematical descriptors that capture both topological and energetic features. These descriptors can then be processed using machine learning techniques to identify patterns, cluster structurally similar aptamers, and prioritize candidates with high potential for target binding. By enabling the systematic downselection of aptamers from large sequence pools, the present disclosure facilitates scalable and interpretable discovery of functional aptamers, including those targeting small molecules.

[0011] Accordingly, provided herein is a computer-implemented method for high throughput classification and functionality assessment of aptamer structures, the method comprising: (a) receiving, by a processor, a plurality of aptamer sequences; (b) generating, by the processor, a structure for each aptamer sequence of the plurality of aptamer sequences to generate a plurality of structures; (c) identifying, for each structure of the plurality of structures, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures; and (d) classifying, by applying the set of parameters of the each structure into a trained machine learning model, the plurality of structures into one or more aptamerclusters. In some embodiments, the computer-implemented method can further comprise (e) determining a target binding characteristic for an aptamer of the one or more aptamer clusters.

[0012] In some embodiments, the computer-implemented method can further comprise generating a visualization of the one or more aptamer clusters. In some embodiments, the structure, for the each aptamer sequence, can be a secondary structure or a tertiary structure based on the each aptamer sequence. In some embodiments, the generating the structure for each aptamer sequence can comprise: optimizing a selection of a set of bonds between pairs of nodes of the each aptamer sequence to achieve a minimum free energy. In some embodiments, the generating the structure for each aptamer sequences can comprise: (i) identifying a plurality of nodes in the each aptamer sequence; (ii) determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the each aptamer sequence; and (iii) selecting, among the plurality of admissible bonds, a set of bonds between the each pair of nodes to form the structure by optimizing to achieve the minimum free energy for the structure. In some embodiments, the selecting the set of bonds can comprise: (i) determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to a face class comprising one or more template graphs of the structure; and (ii) filtering, from the selecting the set of bonds, the admissible bond that fails to conform to the one or more template graphs via a two-level stacking region template. In some embodiments, the face class can comprise a hairpin loop, a stacking region, a bulge loop, an interior loop, or a multi-branch loop.

[0013] Also provided herein is a system for high throughput classification and functionality assessment of aptamer structures, the system comprising: (a) a processor; and (b) memory storing instructions that, when executed by the processor, cause the processor to: (i) receive a plurality of aptamer sequences; (ii) generate a structure for each aptamer sequence of the plurality of aptamer sequences to generate a plurality of structures; (iii) identify, for each structure, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures; (iv) classify, by applying the set of parameters of the each structure into a trained machine learning model, the plurality of structures into one or more aptamer clusters; and (v) determine a target binding characteristic for an aptamer of the one or more aptamer clusters.

[0014] In some embodiments, the memory storing instructions, when executed, can further cause the processor to: generate a visualization of the one or more aptamer clusters. In some embodiments, the structure, for the each aptamer sequence, can be a secondary structure or a tertiary structure based on the each aptamer sequence. In some embodiments, the instructions, when executed, can cause the processor to generate the structure for the each aptamer sequence by: optimizing a selection of a set of bonds between pairs of nodes of the each aptamer sequenceto achieve a minimum free energy. In some embodiments, the instructions, when executed, can cause the processor to generate the structure for the each aptamer sequence by: (i) identifying a plurality of nodes in the each aptamer sequence; (ii) determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the each aptamer sequence; and (iii) selecting, among the plurality of admissible bonds, a set of bonds between nodes to form the structure by optimizing to achieve the minimum free energy for the structure.

[0015] In some embodiments, the instructions, when executed, can cause the processor to select the set of bonds by: (i) determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to a face class comprising one or more template graphs of the structure; and (ii) filtering, from the selection of the set of bonds, the admissible bond that fails to conform to the one or more template graphs via a two-level stacking region template. In some embodiments, the wherein the face class can comprise a hairpin loop, a stacking region, a bulge loop, an interior loop, or a multi-branch loop.

[0016] Further provided herein is a non-transitory computer-readable medium (CRM) having stored thereon computer-readable instructions executable to cause performance of operations comprising: (a) receiving, by a processor, a plurality of aptamer sequences; (b) generating, by the processor, a structure for each aptamer sequence of the plurality of aptamer sequences to determine a plurality of structures; (c) identifying, for each structure, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures; (d) classifying, by applying the set of parameters of the each structure into a trained machine learning model, the plurality of structures into one or more aptamer clusters; and (e) determining a target binding characteristic for each aptamer cluster of the one or more aptamer clusters.

[0017] In some embodiments, the computer-readable instructions can be executable to cause performance of operations further comprising: generating a visualization of the one or more aptamer clusters. In some embodiments, the structure, for the each aptamer sequence, can be a secondary structure or a tertiary structure based on the each aptamer sequence. In some embodiments, the generating the structure for each aptamer sequence comprises: optimizing a selection of a set of bonds between pairs of nodes of the each aptamer sequence to achieve a minimum free energy for the structure. In some embodiments, the generating the structure for each aptamer sequences can comprise: (i) identifying a plurality of nodes in the aptamer sequence; (ii) determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the aptamer sequence; and (iii) selecting, among the plurality of admissible bonds, a set of bonds between nodes to form the structure by optimizing to achieve the minimumfree energy for the structure. In some embodiments, the selecting the set of bonds can comprise: (i) determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to a face class comprising one or more template graphs of the structure; and (ii) filtering, from the selection of the set of bonds, the admissible bond that fails to conform to the one or more template graphs via a two-level stacking region template.

[0018] Additional aspects and advantages of the present disclosure will become readily apparent to those skilled in this art from the following detailed description, wherein only illustrative embodiments of the present disclosure are shown and described. As will be realized, the present disclosure is capable of other and different embodiments, and its several details are capable of modifications in various obvious respects, all without departing from the disclosure. Accordingly, the drawings and description are to be regarded as illustrative in nature, and not as restrictive.

[0019] All publications, patents, and patent applications mentioned in this specification are herein incorporated by reference to the same extent as if each individual publication, patent, or patent application was specifically and individually indicated to be incorporated by reference. To the extent publications and patents or patent applications incorporated by reference contradict the disclosure contained in the specification, the specification is intended to supersede and / or take precedence over any such contradictory material.BRIEF DESCRIPTION OF FIGURES

[0020] The features of the present disclosure are set forth with particularity in the appended claims. A better understanding of the features and advantages of the present disclosure will be obtained by reference to the following detailed description that sets forth illustrative embodiments, in which the principles of the disclosure are utilized, and the accompanying drawings of which:

[0021] FIG. 1 is a flowchart of a non-limiting exemplary method for determining, classifying, and assessing the functionality of aptamer structures. The flowchart describes high-throughput determination of aptamer secondary structures, followed by machine learning-based analysis to prioritize and downselect from thousands of candidate aptamer sequences, which are originally obtained via SELEX, into a smaller, more tractable subset of sequences. These selected candidates are then subjected to experimental validation to assess target affinity and selectivity.

[0022] FIG. 2A is a flowchart illustrating a non-limiting exemplary computer-implemented method for determining, classifying, and assessing the functionality of structures of aptamers according to non-limiting embodiments of the present disclosure.

[0023] FIG. 2B is schematics illustrating non-limiting exemplary template graph consisting of four nodes. Such graph is referred to as a “stack.”

[0024] FIG. 2C is schematics illustrating non-limiting exemplary face classes for calculating minimal free energy configurations, wherein each face class comprises one or more template graphs. In all of the face classes, the rest of the aptamer continues out from each stack.

[0025] FIG. 2D is schematics of visual representation of descriptors used to represent the secondary structure of single-stranded DNA sequences. A non-limiting exemplary secondary structure is depicted as a graph (upper left corner). The structural matrix and Motzkin path descriptors can capture the topological features of the secondary structure graph, while the Bag of Faces (BoF) descriptors convey information about the energetic configuration of the secondary structure.

[0026] FIG. 3 is a block diagram illustrating a computer system 300 upon which embodiments of the present disclosure can be implemented

[0027] FIG. 4 is serotonin aptamer showing secondary structure at Step 1 described in Example 1 (left), the removal of edges (i, j) with \j - i\ < 4 (middle), and Step 2 (right). This is a matrix of all candidate pairs for the secondary structure. The strategy was to eliminate pairings until there was at most one interaction per nucleotide, with the final configuration corresponding to the lowest free energy. In Step 2, all candidate nodes that did not participate in at least one stack connection were filtered out, which greatly reduced the search space prior to Step 3.

[0028] FIG. 5 is a graph illustratingnorm of the adjacency matrix (twice the number of edges) before (solid) and after (dashed) filtering for 1000 random sequences for lengths from 10 to 150, and the ratio of the two (top left). Error bars show the standard deviation.

[0029] FIG. 6 is a graph illustrating running time differences between SeqFold and the presently improved SeqFold on randomly generated single- stranded DNA sequences of various lengths. Ten sequences were computed per string length with a variance of about 0.5%, not shown. Algorithms were performed on AMD Ryzen 95900HX clocked to 3.80 GHz (Base clock 3.30 GHz, boosted to 3.80), NVIDIA RTX 3070 Laptop GPU, 16 GB system ram.

[0030] FIG. 7 is schematics illustrating comparison between SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold, as implemented within the UNAfold software, for an exemplary sequence. The energies were computed by inputtingeach folded structure into the mfold software. The sequence is 5’-GGGACGACGGGGCACATTGTGCT ATTCAGTTGTTCCGCAGGAGAGTCGTCCCGCCTAGCTATTCAGTTGTTCCGCAGGAG AGTCGTCCC-3’ (SEQ ID NO; 6). Structure visualization generated using ViennaRNA, specifically forgi.

[0031] FIG. 8 is schematics illustrating comparison folding of two exemplary sequences between SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold. The figure illustrates that both SeqFold and the presently improved SeqFold consistently fail to predict the junctions represented in the secondary structures produced by the presently described machine learning algorithm and mfold. Two sequences from Table 1 were considered: SI (top row) with a three-way junction and S2 (bottom row) with a four- way junction. Note that Table 1 also mentions the mfold energies associated with each illustrated secondary structure. Structure visualization generated using ViennaRNA, specifically forgi.

[0032] FIG. 9 is schematics illustrating an example reflecting either a difference in coaxial stacking energies or a difference in the multibranch energy function. Comparison folding of an exemplary sequence between SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold as implemented within the UNAfold software. The sequence is an aptamer for Sgc-3b. The sequence being folded is 5’-TTTACTTATTCAATT CCCGTGGGAAGGCTATAGAGGGGCCAGTCTATGAATAAGTTT-3’ (SEQ ID NO: 7). Structure visualization generated using ViennaRNA, specifically forgi.

[0033] FIG. 10 is schematics illustrating an example with non-canonical base pair in one of the stacks. Comparison between SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold as implemented within the UNAfold software. This sequence has been found to be an aptamer for Theophylline. The sequence being folded is 5’-GACGACGATTGTGGTCTATTCATAGGCGTCCGCTGAGTCGTC-3’ (SEQ ID NO: 8). Structure visualization generated using ViennaRNA, specifically forgi. The first three methods do not find the non-canonical base pair.

[0034] FIG. 11 is schematics illustrating a visualization of the chording procedure introduced in Example 5. The sequence structure and Motzkin path found in FIG. 2D. Fewer line gives the associated parenthetical sequence. Red lines connect interacting nucleotides and have corresponding parentheses in the parenthetical sequence.

[0035] FIG. 12 is schematics illustrating two-dimensional representation using t-SNE of Motzkin path descriptors for the dataset of 4450 unique single-stranded DNA sequences.Motzkin path descriptors characterize the topology of the sequences’ secondary structure graph. Each plotted point represents a unique single-stranded DNA sequence. The brighter the color, the higher the count of the sequence. Illustrated are also the secondary structures of sequences x, x2, x3, and x4with counts in the top 0.1 percentile. Next to each of these is plotted an exemplary secondary structure of a single-count sequence in the dataset with the same topology.

[0036] FIG. 13 is schematics illustrating a comparison of the folding algorithms using the high-count norepinephrine aptamer determined by the SELEX process. An initial stem length of four nucleotide pairs was forced in the folding. Notably, the presently described machine learning algorithm gives an intermediate structure between SeqFold and mfold / UNAfold. The sequence is 5 ’ -ACGACGGGGCAC ATTGTGCTGTTCATCTGTTCCGCAGGAGAGTCGT-3 ’ (SEQ ID NO: 9). Structure visualization generated using ViennaRNA, specifically forgi.

[0037] FIG. 14 is bar graphs illustrating the Bag of Faces (BoF) energetic configurations on the x-axes, and the number of sequences with each energetic configuration on the y-axes, for sequences inare the sets of sequences with secondary structure topologically equivalent to one of the four high-count sequences we analyze: q with count 77352, x2with count 67049, x3with count 28432 and x4with count 12126. The heights of the bars are shown in log scale. The numbers on top of each bar are the numbers of sequences with the respective energetic configuration. The bars with stripes are associated with the energetic configurations of high-count sequences.

[0038] FIG. 15 is a t-SNE two-dimensional representation of topic distribution descriptors for our dataset of 4450 unique single- stranded DNA sequences, and a bar plot illustrating the computed clusters and the number of sequences per cluster. In the t-SNE, topic distribution descriptors represent information related to the energetic configuration of the secondary structures of the sequences. Each plotted point represents a unique single- stranded DNA sequence. Each color is associated with a different cluster. The alphabetic characters and the dotted circles highlight the clusters containing at least one sequence with an associated count in the top 0.1 percentile. For each such cluster, we provide an exemplary low-count secondary structure. In the bar graph, bars with stripes are associated with sequences with counts in the top 0.1 percentile. The numbers above each bar are the counts of sequences included in the respective cluster. The heights of the bars are shown in log scale.

[0039] FIG. 16 is schematics illustrating single- stranded DNA secondary structure with inner loop from single base-pair mismatch. Highlighted in the boxes are the terminal mismatches of the internal loop. The notation WX / YZ indicates that W and X are consecutive nucleotides on the DNA strand, and Y and Z are the corresponding consecutive nucleotides on thecomplementary portion of the DNA strand. That is, W pairs or mismatches with Y, and X pairs or mismatches with Z.

[0040] FIG. 17 is schematics illustrating the Cocaine aptamer with sequence 5’-ACAGCTG GGTGAAGTAACTTCCTAAAAGGAACAGAGGG-3’ (SEQ ID NO: 10). Secondary structure as determined by SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold as implemented within the UNAfold software. Structure visualization generated using ViennaRNA, specifically forgi. This example illustrates the different ways each method addresses co-axial stacking.

[0041] FIG. 18 is schematics illustrating comparison folding of an exemplary sequence between SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold as implemented within the UNAfold software. This aptamer binds to a membrane protein of the glioma cell line SHG44. The sequence being folded is 5’-CACAGGTTCCAGGTAATACCTAAGGGTATGCTCTCGCCTATTATATGGAGCAC-3’ (SEQ ID NO: 11). Structure visualization generated using ViennaRNA, specifically forgi.

[0042] FIG. 19 is schematics illustrating comparison folding of an exemplary sequence between SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold as implemented within the UNAfold software. The sequence is an aptamer for Botulinum neurotoxin type A. The sequence being folded is 5’-TTTTATTTTATTTTATTTT AAAAGGCGAATTCAGGGGACGTAGCAATGACTGCC-3’ (SEQ ID NO: 12). Structure visualization generated using ViennaRNA, specifically forgi.DETAILED DESCRIPTION

[0043] The present disclosure describe various exemplary embodiments of systems, software and methods for determining, classifying, and assessing the functionality of structures of aptamers. The disclosure, however, is not limited to these exemplary embodiments and applications or to the manner in which the exemplary embodiments and applications operate or are described herein.

[0044] Provided herein is an improved approach to high-throughput DNA secondary structure determination based on an open-source code (https: / / zenodo.org / records / 7986470). The code has been modified to reduce computational requirements and generate secondary structures more consistent with those predicted by mfold. The present disclosure addresses limitations inherent to the original algorithm through an alternative method based on subgraph matching.

[0045] Provided herein is an open-source Python implementation capable of computing high-throughput minimum free energy (MFE) secondary structures for thousands of DNA sequenceswithin minutes, which is described in detail infra. This implementation can be compared with existing tools to illustrate various computational challenges associated with aptamer secondary structure prediction. The implementation can serve as the initial computational step in a workflow illustrated in FIG. 1 in which aptamer candidates can be subsequently classified according to one or more similarity metrics based on their predicted secondary structures.Evaluating the similarity or dissimilarity among these structures can be essential for applying machine learning techniques to detect structural patterns (e.g., clustering) and infer potential binding properties. Based on these similarity metrics, several machine learning methods can be applied, including dimensionality reduction for data visualization, spectral clustering for group identification, and topic modeling for uncovering common structural features. For topic modeling, the present disclosure adapts a classical “bag of words” framework into a “bag of faces” representation, where structural regions of DNA sequences can be treated analogously to words in a document. The present disclosure applies this computational framework to large pools of single- stranded DNA sequences, including those derived from selection campaigns targeting specific molecules, such as, but not limited to, the neurotransmitter and hormone norepinephrine.

[0046] The present disclosure provides an effective and efficient open-source alternative to existing nucleic acid secondary structure prediction tools. The presently disclosed methods can generate and subsequently analyze a large number of secondary structures using modem machine learning techniques to identify clusters of sequences that exhibit similarity to high-frequency sequences identified through next- generation sequencing (NGS). Through this approach, aptamer motifs associated with target recognition can be revealed, thereby advancing the identification and selection of aptamers.

[0047] Unless otherwise defined, scientific and technical terms used in connection with the present teachings described herein shall have the meanings that are commonly understood by those of ordinary skill in the art. Further, unless otherwise required by context, singular terms shall include pluralities and plural terms shall include the singular.

[0048] As discussed, identifying effective aptamer candidate sequences and accurately modeling structures (e.g., secondary and tertiary structures) from the aptamer candidate sequences is crucial for effective aptamer production. However, current selection processes of aptamer candidate sequences yield a particularly limited range of possibilities of aptamer candidates for further analysis, and this limited range is particularly problematic for the formation of aptamers for targeting smaller molecules. Furthermore, conventional modeling techniques for modeling the secondary and tertiary structures of aptamers are also restricted bythe limited range of aptamer candidate possibilities, and are also inconsistent, computationally inefficient, and lack integration with other processes for aptamer selection, such as the sequencing of aptamer candidate sequences and the analysis of aptamer properties.

[0049] As discussed, desire and need for more comprehensive systems and methods for determining and modeling structures for aptamers that expand the possibilities of aptamer candidates, especially for targeting small molecules. Furthermore, there is a desire and need for such systems and methods to be high throughput, integrated, and computationally efficient in delivering insights (e.g., regarding binding characteristics) for a large number of aptamers.

[0050] The present disclosure describes a novel and nonobvious solution to the aforementioned shortcomings by describing a system and method for accurately, reliably, and efficiently determining, classifying, and assessing the functionality of structures of aptamers. In various embodiments, the systems and methods overcome the limitations of conventional techniques by more comprehensively selecting aptamer sequences based on structuring candidate sequences via subgraph matching to optimize for minimal free energy. Furthermore, the systems and methods provide a computationally efficient and integrated pipeline for the study of aptamers by relying on a machine learning model to cluster large numbers of aptamer candidates and deliver visualization of these clusters along with their properties (e.g., binding characteristics).

[0051] FIG. 1 is an illustration of an example methodology for determining, classifying, and assessing the functionality of structures of aptamers, according to example embodiments of the present disclosure. The example methodology can include, but is not limited to, sequencing 104 a large number of oligonucleotides from a library (e.g., oligonucleotide library 102) for aptamer candidate sequences. In some embodiments, the sequencing 104 can involve binding the oligonucleotides, partitioning (e.g., via beads), and amplification. In some embodiments, the sequencing 104 can involve next generation sequencing (NGS) methods, systems, and / or platforms. As previously discussed, selection of aptamer candidate sequences thus formed 108 for further analysis can be based on the more comprehensive techniques for determining the structure, as described herein. Thus, the exemplary methodology can further include, at block 110, determining (e.g., the modeling) the secondary (e.g., 2 dimensional) structure of each of a large plurality of aptamer sequences, yielding a large plurality of structures (e.g., about 8000 structures). At block 112, the example methodology can include relying on a machine learning model (e.g., a classification model) to down select from the large plurality of candidate aptamer sequences to a tractable and more relevant number of aptamer clusters. The clustering can therefore improve computational efficiency as subsequent analysis can be based on clusters as opposed to individual aptamer sequences. Furthermore, the one or more clusters can be analyzedfor further target affinity and selectivity determination (referred to herein as “binding characteristics) (block 114). In some embodiments, the exemplary methodology can include generating a visualization indicating the binding characteristics and / or the one or more clusters.

[0052] Further details of the present disclosure are provided in the Examples 1-12.I. Overview of Aptamer Identification

[0053] Aptamers can be identified from large combinatorial libraries, theoretically consisting of 4N sequences (> 1 billion), where N is the number of nucleotides in the sequence random variable region. The identification process, called SELEX, generates thousands of aptamer candidate sequences per selection (FIG. 1).

[0054] Aptamer candidates can be sequenced by next- generation sequencing (NGS) and organized hierarchically by numbers of ’hits,’ i.e., the number of times each sequence occurs in the NGS output. While SELEX in combination with NGS is considered high-throughput, yielding thousands of candidates, from a practical standpoint, only a few dozen of the highest copy number sequences can be typically advanced for experimental determination of target binding properties and suitability for application use. As a result, the aptamer candidate sequence space often becomes experimentally intractable, leaving numerous potentially valuable - even highly promising - candidates unexplored.

[0055] In solution, aptamers form secondary and tertiary structures determined by their primary sequences. When an aptamer binds its target, conformational rearrangements facilitate energetically favorable interactions with the target and surrounding solution ions. This results in global changes in aptamer secondary and tertiary structure via intramolecular rearrangements.

[0056] Certain artificial intelligence (Al) tools have been developed for clustering aptamer sequence. Aptacluster uses statistical methods to group aptamer sequences by similarity. In Aptacluster, similarity is defined by the number of nucleotide edits needed to convert one sequence to another. AptaTrace detects motifs associated with minimum face energy (MFE) binding properties (vide infra). Bashir et al., Nat. Commun. 12:2366 (2021) used particle display to partition a library of aptamers by affinity and then train machine learning models on these data to predict affinity. Sun et al., TrAC Trends Anal. Chem. 157:116767 (2022) reviewed computational tools for aptamer clustering. And, Kato et al., BMC Bioinformatics 21:263 (2020) developed a string-based method for cluster analysis. Nonetheless, existing methods do not focus on aptamer secondary structure, information that is critical to the target binding potential.

[0057] Modem machine learning methods enable sorting, comparing, and learning from massive amounts of data, with examples ranging from the analysis of text documents to high-dimensional video sequences (Lai et al., IMA J. Appl. Math. 81(3):409-431 (2016); Merkurjev et al., IEEE International Conference on Image Processing (ICIP) 689-698 (2014)) to large heterogeneous knowledge graphs (Ostaszewski et al., Mol. Syst. Biology 17(10):el0387 (2021)).

[0058] In view of the large number of candidate sequences generated during aptamer selection, there is a significant need for a high-throughput pipeline that incorporates secondary structure analysis. To support such a pipeline, secondary structures must be computed efficiently (FIG. 1). Once determined, machine learning techniques can be applied to organize and classify aptamer candidates based on their secondary structures, and potentially tertiary structures, including those relevant to target binding. Of particular interest are approaches capable of identifying low-frequency aptamer sequences that exhibit high target affinity, selectivity, or other desirable characteristics, which might otherwise be excluded during conventional SELEX workflows.A. Aptamer Chemistry and the SELEX Process

[0059] Aptamers can be identified by an in vitro directed evolution selection method called SELEX (Ellington et al., Nature 346:818-822 (1990); Ellington et al., Nature 355:850-852 (1992)). A large combinatorial library comprising billions of oligonucleotide sequences can be synthesized, wherein each sequence includes a randomized region flanked by two constant regions. The constant regions can be utilized for polymerase chain reaction (PCR) amplification and sequence identification via next- generation sequencing (NGS). The randomized region can serve as the binding site for potential target interaction. The length of the randomized region correlates with the theoretical diversity of secondary structures that can be formed by the oligonucleotides. While longer randomized regions can increase structural diversity and the likelihood of identifying aptamer candidates with target affinity, practical limitations associated with the synthesis, handling, and sampling of longer sequences may necessitate a balance in sequence length selection.

[0060] Large molecular targets, such as proteins, can be immobilized onto a stationary phase, and the oligonucleotide library can be passed through the stationary phase under controlled conditions. Sequences exhibiting low affinity for the target typically pass through, while sequences with higher affinity can be retained by the target-functionalized stationary phase. These higher-affinity sequences can then be eluted through competitive binding with the target in solution and subsequently collected. The collected sequences can be amplified using polymerase chain reaction (PCR), forming an enriched pool that serves as the selection library for subsequent rounds. Selection stringency can be progressively increased with each round to enhance theenrichment of sequences exhibiting improved target affinity and selectivity, including the use of counter- selection steps to remove sequences with undesired binding to non-target species. Upon reaching the desired stringency criteria, next-generation sequencing (NGS) can be performed to identify the primary nucleotide sequences of the candidate aptamers within the final enriched pool (Kohlberger et al., Biotechnol. Appl. Biochem. 69(5): 1771-1792 (2022)). This sequencing step can leverage high-throughput platforms to enable rapid and comprehensive analysis (Komarova et al., Int. J. Mol. Sci. 21(22):8774 (2020)).

[0061] Selection of aptamers for small molecule targets, such as, but not limited to, neurotransmitters, hormones, metabolites, and ions, can be particularly challenging due to the limited number of functional groups available on such targets. These small molecules provide fewer opportunities for noncovalent interactions with aptamers, such as hydrogen bonding, electrostatic interactions, 7t-stacking, and hydrophobic interactions. To address this challenge, a solution-phase SELEX methodology can be employed. In this approach, the selection library can be immobilized to a stationary phase, rather than the target. A short sequence complementary to a portion of one of the constant regions can be affixed to the stationary phase, and the library can be designed such that both constant regions are complementary to one another. The target molecule in solution can then be introduced to the immobilized library. Sequences exhibiting conformational changes upon target binding can be released from the stationary phase and collected. These sequences can undergo further rounds of selection under increasingly stringent conditions. Solution-phase SELEX can preferentially identify stem-loop aptamers, which are characterized by complementary 3’ and 5’ termini that form a stem structure upon target interaction, as demonstrated by the selection strategies described herein infra.

[0062] The in vitro SELEX process, while effective, can be labor-intensive and does not guarantee successful identification of high-affinity aptamers. Computational, or in silico, methods have been developed to support and enhance SELEX workflows both before and after empirical selection steps. However, the vast combinatorial space inherent to oligonucleotide libraries used in SELEX remains beyond the reach of comprehensive empirical or computational exploration. Furthermore, structural data on aptamers that bind to small molecule targets is comparatively limited relative to those targeting larger biomolecules, such as proteins, thereby constraining the utility of existing data-driven modeling approaches. Although various computational techniques have been investigated, they have not yet been fully utilized for small molecule targets in areas such as ab initio prediction of oligonucleotide folding, modeling of aptamer-target interactions, or de novo prediction of functional aptamers directly from primary sequence data.B. Existing Methods for Aptamer Structure

[0063] Several computational approaches can be employed to predict the secondary structures of single- stranded oligonucleotides. In some embodiments, the computational approach can involve determining conformations that minimize the overall free energy of the molecule (Zuker and Stiegler, Nucleic Acids Res. 9(1): 133-148 (1981)). In some embodiments, the computational approach can utilize partition functions to compute base-pairing probabilities across all possible structures (McCaskill, Biopolymers 29(6-7): 1105- 1119 (1990)). Both of these approaches can be typically implemented using dynamic programming algorithms originally introduced in Waterman, Adv. Math. Suppl. Stud. 1:167-212 (1978). Additional, less commonly used methods can include maximum matching algorithms (Nussinov et al., SIAM J. Appl. Math. 35(l):68-82 (1978)), which aim to identify structures with the greatest number of base pairs, and heuristicbased approaches that apply simplified rules based on assumed folding behaviors (Martinez, Nucleic Acids Res. 12(lPartl):323-334 (1984)).

[0064] Multiple software tools are available for predicting the secondary structure of aptamers, including mfold (Zuker, Nucleic Acids Res. 31(13):3406-3415 (2003)), NUPACK (Zadeh et al., J. Comput. Chem. 32(1): 170- 173 (2011)), RNA Composer (Sarzynska et al., Proteins: Struct. Funct. Bioinform. 91(12):1790-1799 (2023); Popenda et al., Nucleic Acids Res.40(14):ell2 (2012)), web 3DNA (Lu and Olson, Nucleic Acids Res. 31(17):5108-5121 (2003); Lu and Olson, Nat. Protoc. 3(7): 1213- 1227 (2008); Li et al., Nucleic Acids Res. 47(W1): W26-W34 (2019)), and ModeRNA (Rother et al., Nucleic Acids Res. 39(10):4007-4022 (2011), among others.

[0065] Mfold is generally regarded as the benchmark tool in this domain. It applies free energy minimization algorithms to identify the lowest-energy folded conformation of a given nucleotide sequence, with user-defined inputs such as temperature and ion concentrations.Validation of predicted structures typically requires empirical confirmation through experimental methods. The mfold software is implemented in C / C++ and is available under a free license for academic and nonprofit use, while commercial use may require a paid license.

[0066] SeqFold is an open- source alternative to mfold that also utilizes the dynamic programming approach originally developed by Zuker and Stiegler, Nucleic Acids Res. 9(1): 133-148 (1981). In contrast to mfold, SeqFold is implemented in Python, which is widely used for the development and deployment of machine learning (ML) algorithms. This implementation facilitates integration into computational workflows that process large sets of single-strandedDNA sequences, such as aptamer candidates, by enabling the generation of corresponding minimum free energy (MFE) secondary structures and subsequent application of machine learning-based analytical methods.

[0067] The current implementation of SeqFold includes certain limitations, such as inconsistencies between the dynamic programming algorithm originally proposed by Zuker and Stiegler and the manner in which the algorithm is implemented in the codebase. These inconsistencies can result in the prediction of minimum free energy (MFE) secondary structures that deviate significantly from those generated by mfold, which is treated as ground truth for comparative purposes.

[0068] The present disclosure demonstrates that machine learning methods can serve as effective tools for analyzing large collections of DNA secondary structures. Machine learning techniques can be applied to a variety of tasks, including visualization of the structural configuration space, identification of clusters of structures with shared properties, and prediction of target binding activity based on structural features.C. Multigraphs and Graph Matching for DNA

[0069] The present disclosure incorporates concepts from Yang et al., IEEE Trans. Netw, Sci. Eng. 10(4):1846-1862 (2023); Moorman et al., IEEE Trans. Netw, Sci. Eng. 8(2):1367-1384 (2021); and Tu et al., IEEE International Conference on Big Data (Big Data) 2575-2582 (2021) in subgraph matching for multiplex networks, also referred to as labeled directed multigraphs.

[0070] A labeled multigraph can be defined as:G = (V, E, L, C),where V represents the set of nodes, and E represents the set of edges, which are undirected in this context. The labeling function L assigns a label to each node, corresponding to nucleotide bases such as adenine (A), thymine (T), guanine (G), or cytosine (C). The channel assignment function C assigns a channel to each edge, where channels can represent different types of molecular interactions, such as backbone connectivity or nucleotide base-pairing interactions.

[0071] In this context, perfect matching can be understood as a subgraph isomorphism from a template graph Gtto a background graph G (sometimes referred to as the world graph).Subgraph isomorphism is known to be NP-complete. Subgraph matching algorithms date back to Ullmann, J. ACM 23(l):31-42 (1976), which utilizes a tree search approach, maintaining a search state and traversing the tree of possible search states, with backtracking employed when the end of a branch is reached. Due to the size of the search tree, computational complexity can be reduced by refining the search space at each step to avoid unnecessary branches. Other treesearch methods include VF2 (Cordelia et al., IEEE Trans. Pattern Anal. Mach. Intell. 26(10): 1367-1372 (2004)), and its variants VF2 Plus (Carletti et al., VF2 Plus: An Improved Version of VF2 for Biological Graphs, GbRPR, Springer 9069:168-177 (2015)), VF3 (Carletti et al., Introducing VF3: A New Algorithm for Subgraph Isomorphism, GbPRP 10310:128-139 (2017)), VF2++ (Jiittner and Madarasi, Discrete Appl. Math. 242:69-81 (2018)), and for certain graph structures, RVRI-DS (Bonnici et al., BMC Bioinformatics 14(7): S 13 (2013)). Techniques such as constraint propagation and filtering have been shown to significantly reduce computation time. The subgraph matching problem can be generally recognized as having combinatorially complex solution spaces; in Yang et al., IEEE Trans. Netw, Sci. Eng. 10(4):1846-1862 (2023), instances of subgraph matching problems have been shown to exhibit 10100or more isomorphisms.

[0072] In the present disclosure, the aptamer secondary structure folding problem can be viewed as a bipartite incomplete graph matching problem in which nodes in the DNA graph are matched to conjugate nodes within the same structure. A secondary structure can be represented by creating a second copy of the DNA primary sequence and searching for a bipartite matching between the DNA sequence and its copy. A constraint can be applied such that any base in one copy can only pair with a single conjugate base in the second copy. As such, this constitutes an incomplete or inexact matching problem, with the optimal matching defined by the minimum free energy (MFE). Existing algorithms for DNA secondary structure prediction can be constructed similarly to tree search and backtracking method in Ullmann, J. ACM 23(1):31-42 (1976). The present disclosure introduces additional techniques from constraint propagation to improve the speed of DNA secondary structure algorithms. The presently described approach can then be used to analyze thousands of secondary structures simultaneously using modern machine learning methods.II. Computer Implemented Methods of the Disclosure

[0073] FIG. 2A is a block diagram illustrating an example computer-implemented method for determining, classifying, and assessing the functionality of structures of aptamers according to non-limiting embodiments of the present disclosure. One or more blocks or processes described in the blocks of FIG. 2 can be performed by one or more computing devices (e.g., such as but not limited to computing system 300 as will be described herein). For example, the one or more blocks or processes can be performed by the processor 304 based on instructions provided by any one of, or a combination of memory components 306 / 308 / 310 and user input (e.g., provided via the input device 314), as are discussed herein. Refer to Examples 1-12 for further informationregarding the computer implemented methods provided herein, in accordance with various embodiments.

[0074] In various embodiments, at block 202, a processor (e.g., processor 304) can receive a plurality of aptamer sequences. For example, the aptamer sequences can be generated using an oligonucleotide library. In some embodiments, the aptamer sequences can be generated via NGS of a biological sample. The aptamer sequences can include a large plurality of aptamer sequences (e.g., about 8000) for further analysis.

[0075] At block 204, the computing device can generate a secondary structure for each aptamer sequence of the plurality of aptamer sequences to determine a plurality of secondary structures. In some embodiments, generating the secondary structure can include optimizing a selection of a set of bonds between pairs of nodes of the aptamer sequence to achieve a minimum free energy. In some embodiments, generating the secondary structure can include identifying a plurality of nodes (e.g., nucleotides) in the aptamer sequence, and then determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the aptamer sequence. In some aspects, an admissible bond can be based on chemical rules pertaining to allowable set of nucleotides that can bind to another nucleotide (e.g., an adenine (A) can bind with thymine (T), a cytosine (C) can bind with guanine (G), etc.). Furthermore, a set of bonds between nodes, among the plurality of admissible bonds, can be selected to form the secondary structure by optimizing to achieve a minimum free energy for the secondary structure. In some embodiments, selecting the set of bonds can include: determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to one or more template graphs of the secondary structure; and filtering, from the selection of the set of bonds, the admissible bond that fails to conform to the one or more template graphs.

[0076] As shown in FIG. 2C, the one or more face classes comprising one or more template graphs can include but are not limited to a hairpin 220, an inner loop 230, a bulge 240, or a multibranch 250. In some embodiments, each face class can be distinguishable based on the positioning of a stacking region 262. In at least one embodiment, the stacking region 262 may refer to a grouping of four nodes (e.g., nucleotides) as shown in FIG. 2B, from which linear aptamer sequences can extend from. Refer to Examples 1-12 for further information regarding the face classes comprising one or more template graphs.

[0077] At block 206, the computing device can identify, for each secondary structure, a respective set of parameters assessing a similarity criteria between the secondary structure and another secondary structure of the plurality of secondary structures. For example, the set of parameters can include or comprise a set of vector- valued descriptors of the secondary or tertiarystructures. The vector- valued descriptors represent the structures with a set of numerical values that summarize information about the structures’ properties. In at least one embodiment, the set of parameters can refer to distances between vector- valued representations of each structure to define the similarity criteria between each aptamer structure (e.g., the larger the Euclidean distance between two vector valued descriptors of respective aptamer structures, the less similar the aptamer structures. It is also contemplated that other metrics in addition to or alternative to Euclidean distance can be used to quantify the similarity between the set of parameters (e.g., vector valued descriptors), such as but not limited to Manhattan distance and cosine similarity.

[0078] In some embodiments, for example, as shown in FIG. 2D, the set of parameters (e.g., vector valued descriptors) for a given structure 260 can include but are not limited to the following: an adjacency matrix of a graph representing the structure (referred to herein as the structural matrix) 272; a Motzkin path vector 274, or combinatorial objects to study lattice paths; and bag-of-faces (BoF) descriptors 276 used for natural language processing.

[0079] The structural matrix descriptor 272 can comprise of the adjacency matrix of the structure graph, where, in at least one embodiment, the structural matrix descriptors can provide information on the topology of the structures.

[0080] In at least one embodiment, Motzkin path vector 274 can be a path of length / Fin H X H starting at (0, 0) and ending at (0, / F) with each step being one of the three following types (+1, +1), (+1, 0), (+1, —1). The path can be recorded as a balanced parenthetical sequence with rests, i.e., a sequence vG {(,., )}F which is balanced in the parentheses. One can go from such a sequence to a path by the map (■-> (+1, +1),. ■-> (+1, 0), ) ■-> (+1, -1) and the balanced condition ensures that one never cross below the x axis. As shown in FIG. 2D, each path step steps one to the right. The Motzkin path vectors can further represent the topology of the structure based on a path indicating steps across a sequence of nucleotides. By comparing the Motzkin path vectors, one can identify structural similarities between different folded sequences. For example, the cumulative sum of the Motzkin path vector can be used to describe the secondary structure of the sequences in numerical experiments. In some embodiments, various dimensions can be reduced where the vectors are identical across all data points. This additional step can eliminate redundant information from the Motzkin path descriptors, resulting in a more compact representation.

[0081] In some embodiments, the bag of faces (BoF) descriptors 276 can be based on the folded structure of the aptamer candidates being characterized by one or more faces and their associated energies. The BoF descriptor is thus based on a count of the number of faces and their respective energy configurations in a given folded sequence of a structure.

[0082] In some embodiments, the set of parameters for the similarity criteria can indicate a similarity to structures that are known to have ‘good’ binding properties (e.g., those structures are associated with a high count). The similarity criteria can thus be structural similarity and energetical similarity. Structural similarity can pertain to the topology of folded sequences of the candidate aptamer, while energetical similarity can refer to the energetical configurations of the structures, such as the energy associated with the faces composing the secondary structure.

[0083] Using the Motzkin path 272 and / or the structural matrix 274 descriptors, sequences of the structure can be identified with other structures topologically equivalent to that of a given sequence of interest. Using the bag of faces (BoF) descriptors 276, such sequences can be identified with other structures with similar energetical configurations.

[0084] At block 208, the computing device can classify, by applying the set of parameters of each secondary structure into a machine learning model, the plurality of structures into one or more aptamer clusters. The machine learning model can be an unsupervised machine learning model configured to cluster a plurality of graphs (e.g., structures) and / or sets of parameters respectively associated with the graphs into one or more clusters. Example unsupervised machine learning models can include but are not limited to: spectral clustering, k-means clustering, hierarchical clustering, mixture of Gaussians, and the like. For example, in at least one embodiment, the computing device can perform the classification by building latent topic models using the BoF descriptors 276 of the plurality of structures. The computing device can then identify topic distributions based on the set of parameters of each structure to describe and cluster the plurality of structures. In some embodiments, the topic modeling can rely on Nonnegative Matrix and / or other linear approaches. For example, occurrences of face-energy configurations across the structures (e.g., as measured by the BoF descriptors 276) as well as the respective sequence can be inputted into NMF. The clusters of aptamer candidates can thus be identified based on being described by similar topic mixture distributions. The distance between any two secondary structures can be measured by comparing how dissimilar their topic mixture distributions are.

[0085] In some embodiments, the datasets (e.g., sets of parameters) applied to the machine learning model can be improved or can undergo dimensionality reduction for more efficient computation. In particular, the dimensionality reduction can ensure the maintenance of the relevant intrinsic geometry of the aptamer structures. In at least one embodiment, the dimensionality reduction can be performed via t-distributed stochastic neighbor embedding (t-SNE) and / or via principal component analysis (PC A).

[0086] In some embodiments, more than one machine learning model can be used. The machine learning model can be unsupervised or semi-supervised (e.g., a clustering algorithm, a support vector machine, etc.).

[0087] In an alternative embodiment, the machine learning model can be trained based on a training data set comprising a plurality of reference sets of parameters with known or labeled clustering results. In some embodiments, machine learning models can include logistic regression techniques, linear discriminant analysis, linear regression analysis, artificial neural networks, machine learning classifier algorithms, or classification / regression trees. In various other embodiments, machine learning systems can employ Naive Bayes predictive modeling analysis of several varieties, learning vector quantization artificial neural network algorithms, or implementation of boosting algorithms such as Adaboost or stochastic gradient boosting systems for iteratively updating weighting to train a machine learning classifier to determine a relationship between an influencing attribute, such as received environmental data, and a system or environmental characteristic and / or a degree to which such an influencing attribute affects the outcome of such a system or environmental characteristic.

[0088] At block 210, the computing device can determine a target binding characteristic for each aptamer cluster. In various embodiments, the target binding characteristic can refer to an ability of an aptamer structure or cluster to bind to a target antigen, protein, and / or molecule of interest. For example, the target binding characteristic can include but is not limited to a measurement of target affinity or a target selectivity. In some embodiments, the target binding characteristic can be determined in silico, based on stored information for a target antigen, protein, and / or molecule of interest. Also or alternatively, the target binding characteristic can be determined via assaying techniques performed manually.

[0089] At block 212, the computing device can generate a visualization of the one or more aptamer clusters. For example, the computing device can display the structure or structural characteristics of each aptamer clusters. In some embodiments, variations within a structure of an aptamer cluster as it relates to individual aptamer structures of an aptamer cluster can be noted in the visualization. In some embodiments, functional regions of the aptamer cluster can be displayed in the visualization. In some embodiments, the target binding characteristic of the aptamer cluster or of one or more functional regions of the aptamer cluster can be displayed in the visualization.A. Subgraph Matching and Free Energy Minimization for DNA Folding

[0090] Provided herein is a full algorithm for DNA secondary structure prediction. General ideas from Yang et al., IEEE Trans. Netw, Sci. Eng. 10(4): 1846-1862 (2023); Moorman et al., IEEE Trans. Netw, Sci. Eng. 8(2):1367-1384 (2021); and Tu et al., IEEE International Conference on Big Data (Big Data) 2575-2582 (2021) on subgraph searches in multiplex networks can be followed. Each DNA strand can be treated as a linear string of nodes, with node labels corresponding to the primary sequence of nucleobases - adenine (A), thymine (T), guanine (G), and cytosine (C). The graph-matching problem can be solved using hierarchical filters applied to candidate base-paired nodes, followed by an exhaustive tree search over the remaining solution space. The tree search can be accelerated using principles such as structural equivalence and node cover (Yang et al., IEEE Trans. Netw, Sci. Eng. 10(4): 1846-1862 (2023)), or by ordering the search according to decreasing energy (Tu et al., IEEE International Conference on Big Data (Big Data) 2575-2582 (2021)).

[0091] The DNA secondary structure problem can be characterized as an annotated incomplete subgraph matching problem, in which a graph is matched to itself through internal interactions between nodes (nucleotides). The objective can be to identify the optimal number and arrangement of such interactions that minimize the overall free energy of the structure.

[0092] Non-limiting exemplary methods on an elimination scheme is described in detail in Example 1.

[0093] In some embodiments, Step 1 can involve representing a nucleic acid strand as a linear graph, with edges defined between adjacent nodes corresponding to elements of the primary sequence. A mapping can then be performed between the strand and a duplicate of itself according to predetermined rules. For purposes of graph matching, sequence elements (e.g., nucleotides such as adenine, thymine, cytosine, and guanine) can be used as labels on the graph nodes, assigned in the order corresponding to the primary sequence.

[0094] A list of candidate nodes in the target graph to which each node can potentially be matched can be generated. The matching rules can be defined by canonical base pairing principles. For example, in some embodiments, for each adenine (A) in the primary sequence, all thymine (T)-labeled nodes in the duplicate can be identified as potential matches; similarly, guanine (G) nodes can be matched with cytosine (C)-labeled nodes. This process can be extended for other bases as well, establishing a complete set of candidate matched nodes. These candidate pairings can then be refined through systematic elimination steps until an optimal structure can be determined. Non-limiting exemplary output from this step is illustrated in FIG.4 (left panel).

[0095] In some embodiments, Step 2 can function as a topological filtering stage, analogous to the “topology” filter described in Moorman et al., IEEE Trans. Netw. Sci. Eng. 8(2): 1367-1384 (2021). During this step, constraints can be applied to eliminate candidate node matches that do not preserve the first-order connectivity present in the reference or template graph, as required by the subgraph matching problem.

[0096] For structures such as nucleic acid aptamers, the filtering process can be conducted based on predefined structural configurations, such as stacks. In certain embodiments, a stack can be defined by two adjacent matched node pairs, with connectivity verified through both primary sequence adjacency and complementary matching. Candidate node pairs that are not part of a valid stack configuration can be excluded. For instance, in some implementations, matched edges can be classified as “admissible” if they satisfied these topological constraints.

[0097] Illustratively, a stack can include two adjacent nodes in the primary structure connected by backbone edges, along with complementary base-pair edges linking these nodes to corresponding nodes in the duplicated strand. In some embodiments, as illustrated in FIG. 2B, primary edges can be shown in black and base-pair edges in red. The remaining secondary edges after this step can be designated as admissible, and their distribution can be visualized, for example, in the right panel of FIG. 4. In some cases, this filtering process can remove over half of the initial candidate edges, as depicted in FIG. 5.

[0098] In some embodiments, Step 3 can involve refining the set of matched candidate node pairs to identify a configuration corresponding to the minimum free energy, with the constraint that at most one candidate pairing can be selected per nucleotide. This optimization can require knowledge of structural motifs (faces) and their associated energy values. Five distinct face types can be considered: hairpin loops, stacking regions, bulge loops, interior loops, and multibranch loops, which can serve as subgraph groupings to guide the selection of matched pairs. Each face type can be assigned a scalar energy, and the total energy of a predicted structure can be defined as the sum of the energies for all faces present.

[0099] The goal of this step can be to determine the configuration of admissible base-pair interactions that minimizes the total energy across the sequence. In some implementations, this can be accomplished through an exhaustive or guided search over all substrings of the sequence. For each candidate base pair (i, j), the algorithm can evaluate the energy of all admissible faces that could include the pair, beginning with shorter substrings and extending to larger regions of the sequence. The pairing with the minimal energy contribution can be identified and designated as the “last pair,” ensuring that no additional faces exist beyond the bounds of the selectedsubstring. The set of nucleotide pairs involved in the optimal energy configuration can then be selected as the final structural prediction.

[0100] While energy optimization for simpler face types, such as hairpins, interior loops, and bulges, can be performed using relatively straightforward heuristics, optimization involving multibranch structures can require exploring all valid subsets of admissible edges, creating a significantly larger and more complex search space. In some embodiments, this can be addressed using recursive dynamic programming techniques, which can be further detailed in Examples 10-12.

[0101] A recursive strategy, originally proposed by Zuker and Stiegler, Nucleic Acids Res. 9(1): 133— 148 (1981), can be employed to compute the minimum free energy configuration. This approach can also be used in software implementations such as mfold and SeqFold. Structural information can be stored in data structures (e.g., caches), such as Svand Sw, which can map index pairs (i, j) to their corresponding energy values and associated base-pairing configurations. In some embodiments, these caches can be implemented as lists of lists, with each element containing an energy value and a structural configuration.

[0102] The cache Sv(i, j) can represent the minimum energy configuration of the substring [i, j] under the assumption that i and j form a base-pair interaction. Conversely, sw(i, j) can represent the optimal configuration of the same substring without enforcing any interaction between i and j. If the optimal configuration includes a base-pair interaction between i and j, then the values in the two caches can be equal, i.e., Sv(i, j) = Sw(i, j). The recursive computation for Swcan be defined as follows:Formula (1)

[0103] To compute Sv(i, j), the algorithm can evaluate all valid face configurations that end at the (i, j) position and can select the one with the lowest total energy. This recursive process can facilitate efficient identification of the minimum free energy structure.1. Improved SeqFold Algorithm

[0104] Provided herein is an improved version of the SeqFold code, which can implement multiple enhancements for increased computational efficiency and structural accuracy.

[0105] Cache efficiency can be achieved by avoiding unnecessary recursive operations through a defined order of cache computation. One approach involves calculating Sv(i, j) for which li - j\ = k, beginning with k = 4 and proceeding incrementally through k = 5, 6,..., n. Alternatively, for each j = 5, 6,..., n, S w(i, j) can be computed for i = j - 4, j - 5,..., 1. Both approaches rely on the principle that, to compute Sw i, j), only values of Sw i', j') for all i < i' < j'< j are required. By computing cache values in one of these orders, recursive calls can be minimized, resulting in a reduction in computational cost on the order of n.

[0106] Graph-matching-based filtering can be used in the presently improved SeqFold to enhance prediction accuracy. In contrast to the original SeqFold implementation, a subgraph matching algorithm can be applied to identify all embeddings of the stack graph. This yields a restricted set of nucleotides and allowable edges, based on the assumption that isolated interactions do not form. In earlier algorithms, such isolated interactions were penalized by assigning high energy values (e.g., 1600 kcal / mol). This matching process can be executed with a complexity of O(n2), which offers computational advantages over previous methods. By constraining calculations to these filtered edges, the number of operations involved in loop processing, multibranch energy calculations, and the evaluation of S\ can be reduced. If a given pair (i, j) is not included in the valid edge set, it can be assumed that Sv(i, j) > Sw(i, j), and the calculation can be omitted.

[0107] Energy corrections have also been incorporated into the presently improved SeqFold. Specific inaccuracies in energy assignments present in the original SeqFold code were identified and corrected. Further implementation details regarding these corrections are provided in Example 9.

[0108] In some embodiments, the presently improved SeqFold and the presently described machine learning algorithm can be primarily based on the functions utilized in SeqFold for computing energy values associated with various structural features, such as hairpins, bulges, and internal loops. However, improvements can be introduced to these energy functions.Specifically, the computation of energies associated with internal loops generated by single basepair mismatches and with junctions, which can be multi-branch structures having no unpaired nucleotide between branches, can be modified.

[0109] The energy associated with internal loops resulting from single base-pair mismatches can be determined based on the energy contributions of the terminal mismatches of the loop. Terminal mismatches, which can occur at the ends of double-stranded loop regions, can comprise four nucleotides: two base-pairing nucleotides and two unpaired nucleotides internal to the loop. A left and right terminal mismatch can be associated with each loop, depending on the 5’ to 3’ orientation. In some embodiments, visual representations such as diagrams can be used to highlight these mismatch regions. The free energy associated with such loops can be computed as the sum of the energies of the left and right terminal mismatches. These energy values can be derived from enthalpy and entropy data using a specified formula. Lookup tables (e.g., Tables 2-5 herein) can be employed to assign energies to terminal mismatches, with different tables used depending on whether the mismatch includes sequence-terminal nucleotides.

[0110] In some embodiments, conventional implementations such as SeqFold can compute the energy of single base-pair mismatch loops inaccurately due to an incorrect determination of the left terminal mismatch. For example, such implementations can identify the terminal mismatch using the base pairs at the ends of stacking regions rather than the actual unpaired loop nucleotides. This issue can be addressed by reassessing and correcting the mismatch identification methodology.

[0111] With respect to junction energy computation, conventional approaches such as SeqFold can associate a junction with an energy of the form < Em= 4.6 + Est, where Est represents the total energy of individual branches or stacks. In contrast, the presently described approach can associate a lower junction energy, such as < Em= 0.6 + Est, thereby treating these structures as more thermodynamically stable. This adjustment can be based on heuristic tuning rather than experimental validation and can be designed to emulate behaviors observed in alternative folding software, such as mfold.

[0112] In some embodiments, it can be observed that conventional tools, including SeqFold and the improved SeqFold, can fail to predict the presence of junctions in cases where other tools, such as mfold and the presently described machine learning algorithm, succeed.Comparative analyses and examples can be provided to illustrate these discrepancies. It can be further speculated that such failures are not solely due to energy computation methods but also to limitations in the underlying dynamic programming implementation. For example, certain algorithms might not consider junctions as valid structural features (or “faces”), which can limit their ability to accurately predict complex secondary structures.2. Energy Optimization and Algorithmic Complexity

[0113] In some embodiments, each step of the energy minimization algorithm can involve making decisions about the type and size of structural elements, referred to as faces, that can form between nucleotide pairs. It can be assumed that a face can terminate at a specific nucleotide pair within a sequence substring. If no internal interactions are present, the face can be categorized as a hairpin loop. If a single internal interaction exists, it can be classified as an inner loop, bulge, or stack. Structures with multiple internal interactions can be identified as multibranch loops.

[0114] To determine the optimal configuration, the lowest energy arrangement can be found by evaluating all admissible face types, including hairpin loops, inner loop / bulge / stackconfigurations, and multibranch structures. The hairpin can require minimal computation due to its single configuration. Inner loops, bulges, and stacks, which can be characterized by two internal base pairs, can be evaluated through a search over admissible interactions, often resulting in polynomial time complexity. However, in some implementations, the number of required evaluations can be significantly reduced by using prior graph-matching steps to filter infeasible interactions.

[0115] Optimization for multibranch loops can be more complex, as it can involve exploring many combinations of internal interactions. In some implementations, assumptions about decomposing the sequence into smaller subsections can be used to reduce the computational burden. Recursive algorithms, similar to those used in prior methods, can enable efficient determination of multibranch configurations by reusing previously computed substructures.

[0116] In some embodiments, a parameter can be introduced to restrict the number of branches in multibranch loops, allowing further control over computational cost. The reduction in computational load can also benefit from graph-based pruning steps conducted earlier in the algorithm.

[0117] Overall, the time complexity of face optimization can depend on the sequence length and the allowed complexity of multibranch structures. In practice, runtime performance can vary between implementations depending on algorithmic choices and preprocessing steps. Empirical evaluation can reveal performance differences between methods, as illustrated by experimental comparisons.

[0118] Non-limiting exemplary implementation of this section is described in detail in Example 2.

[0119] In some embodiments, each step of an energy minimization algorithm can involve a decision regarding the type and size of a structural feature, or “face,” to be formed. A face can be assumed to terminate at a nucleotide pair (i, j) within a given substring [i, j] of a nucleic acid sequence. If no internal interactions are present, the face can be classified as a hairpin loop. If a single internal interaction is present, the face can be considered an inner loop, bulge loop, or stack. When two or more internal interactions are present, the face can be classified as a multibranch loop. Determining an optimal face configuration can require identifying the most energetically favorable arrangement of internal base pairs. The present disclosure describes systematic processes for computing these configurations and their associated energy values.

[0120] In some embodiments, under the assumption that nucleotides i and j interact, the lowest energy face structure can be determined by evaluating three types of face configurations: a hairpin loop, possible inner loop / bulge / stack combinations, and potential multibranch structures.For hairpins, the configuration space can be limited to a single structure, which can be computed with sublinear complexity.

[0121] The inner loop, bulge, and stack structures can be treated as a single class, as each can involve two interior base pairs. A search over all admissible edges - identified during a graph matching phase - within the substring (i, j) can be conducted to identify the configuration with the lowest energy. Although the worst-case computational complexity for this evaluation can be on the order of O(\j - il2), the number of evaluations can be significantly reduced in some embodiments due to earlier pruning of edge candidates.

[0122] In other embodiments, multibranch optimization can involve evaluating various combinations of admissible base pairs within the substring [i, j] that can give rise to multibranch configurations. This can result in a combinatorially large search space. However, existing approaches, such as those inspired by Zuker and Stiegler, can reduce the search complexity by assuming the existence of an intermediate nucleotide position k within [i, j] such that an optimal multibranch configuration can be formed by combining optimal substructures from [i+1, k] and [k, j-1], Algorithms incorporating this approach can reduce multibranch search complexity to linear time with respect to substring length.

[0123] In some embodiments, a control parameter m can be introduced to limit the number of branches in a multibranch structure, which can constrain the search to a subset of edge combinations involving at most m branches. This bounded complexity can be further reduced through prior graph-based filtering.

[0124] The combined complexity of computing an optimal face for each nucleotide pair (i, j) can, in some embodiments, be on the order of O(\j - il2) for traditional algorithms based on the Zuker-Stiegler approach and O(\j - il2m) for algorithms incorporating the above-described improvements.

[0125] In some implementations, this energy minimization can be executed across all admissible nucleotide pairings. In the absence of a precise estimate for the number of such pairings, a conservative upper bound based on all possible (i, j) pairs can be applied. This results in overall complexity estimates on the order o) for traditional implementations and O(n2m+2) for implementations using the improved method. In one non-limiting embodiment, the parameter m can be set to 4.

[0126] Experimental comparisons can demonstrate differences in runtime performance between implementations. For example, empirical scaling behavior can show that conventional implementations scale with a computational complexity of approximately < (n3'696), whileimproved implementations described herein can scale more efficiently, such as on the order of 풪(n3.477), as described in details in Example 2.B. Mathematical Representation of Secondary Structures

[0127] In some embodiments, multiple vector-valued descriptors can be utilized to mathematically represent and analyze the secondary structure of single-stranded nucleic acid sequences, such as DNA or RNA. These descriptors can encode both the topological and energetic characteristics of folded sequences in a numerical form suitable for computational processing, machine learning, and clustering analyses. In certain implementations, a first descriptor can be derived from the adjacency matrix of a graph representation of the secondary structure and can be referred to as a structural matrix. A second descriptor can utilize Motzkin paths, which are combinatorial constructs traditionally applied in the analysis of lattice paths in mathematics. A third descriptor type can be based on concepts from natural language processing, such as the “Bag of Words” methodology, and can be adapted to quantify features of nucleic acid structures in a manner analogous to how word frequency is used in textual analysis.

[0128] Collectively, these descriptors can provide complementary insights: the structural matrix can describe connectivity, the Motzkin path can encode topology, and the Bag of Faces can characterize energetic and compositional aspects of the structure. In some embodiments, one or more of these descriptors can be used individually or in combination to improve the accuracy, interpretability, and computational efficiency of secondary structure prediction and analysis.1. Structural Matrix Descriptor

[0129] In some embodiments, the structural matrix descriptor can be defined as the adjacency matrix corresponding to the secondary structure graph of a single-stranded nucleic acid sequence, such as DNA or RNA. In this representation, each nucleotide can be modeled as a node, and each base-pair interaction between nucleotides can be modeled as an edge connecting the respective nodes. To isolate the structural information relevant to folding topology, matrix entries corresponding to primary covalent bonds along the nucleotide backbone can be set to zero, such that only the secondary (hydrogen-bond) connections remain represented. This filtered adjacency matrix can therefore encode the global topology of the folded structure while excluding trivial backbone connectivity.

[0130] In certain embodiments, the structural matrix can be a symmetric binary or weighted matrix, where each nonzero element represents the presence or strength of a pairing interaction between two nucleotides. The dimensionality of the matrix can correspond to the total length of the nucleotide sequence, with indices reflecting nucleotide positions along the 5 ’-3’ direction. Because this descriptor captures only the connectivity pattern between paired regions, it can serve as a purely topological representation of the secondary structure.

[0131] In some embodiments, the structural matrix descriptor can be used to identify structural motifs or topological similarities among different folded sequences. For example, sequences that differ in nucleotide composition but share similar folding topologies can exhibit comparable adjacency structures. This allows for structure-based comparison independent of sequence identity, which can be particularly valuable in identifying functionally analogous aptamers or motifs across diverse sequence libraries.

[0132] In some embodiments, the adjacency matrix representation can be processed computationally to derive higher-level structural metrics, such as degrees of connectivity, clustering coefficients, or topological invariants. These metrics can further be incorporated into downstream analyses, such as structural similarity searches, clustering, or machine learningbased classification. Thus, the structural matrix descriptor can provide a compact yet information-rich representation that captures the essence of the folding topology without requiring direct reference to base composition or energetic values.2. Motzkin Paths Descriptor

[0133] In some embodiments, a second descriptor can be based on Motzkin paths, which are combinatorial mathematical objects used to describe non-crossing paths in a lattice grid. In this context, the Motzkin path can encode the secondary structure of a nucleic acid sequence by assigning directional values to each nucleotide position that reflect whether a nucleotide is paired, unpaired, or forms part of a base-pair bridge. The path can be represented as a sequence or vector of integers (e.g., 1, 0, -1) or as cumulative height values corresponding to the number of open pairings at each position. This representation can preserve information about the structural topology of the sequence (such as loops, stems, and unpaired regions) while avoiding redundancy and maintaining compatibility with geometric and statistical analysis methods.

[0134] In some embodiments, a nucleic acid secondary structure descriptor can be based on a Motzkin path formalism to capture the topology of the folded sequence. A Motzkin path of length N can be defined on the lattice N x N, starting at (0, 0) and ending at (0, N), with each step being one of (+1, +1), (+1, 0), or (+1, -1). Such a path can be recorded as a balancedparenthetical sequence s E {(,.,)}w, where the mapping (■-> (+1, +1),. ■-> (+1, 0), ) ■-> (+1, - 1) is used to translate between the parenthetical form and the lattice path form. Thebalanced-parentheses condition can ensure that the path never crosses below the x-axis. Since each step moves one unit to the right, the sequence can alternatively be represented as a vector v E {—1, 0, 1}W, where all partial sums are non-negative, or equivalently as a cumulative sum vector w, whereVj. A sequence scan be admissible if= sN= 0 and lsi+i— st I 1 for all i - conditions that mirror the non-crossing interactions in a nucleic acid secondary structure without pseudoknots.

[0135] In some embodiments, this Motzkin-path descriptor can be used to represent the secondary structure of a single-stranded DNA sequence by laying out the nucleotide backbone along an A-gon (or circle), placing open parenthesisat the nucleotide where a pairing begins, a dot at positions without pairing, and a close parenthesiswhere the pairing ends as one moves around the circle (e.g., clockwise from a fixed starting point). Each balanced parenthetical sequence then corresponds to a non-crossing chord diagram of the A-gon, where each chord corresponds to a nucleotide pair. Because DNA sequences in many cases avoid pseudoknots (i.e., crossing base-pair interactions), the mapping is bijective between the folding topology and the Motzkin path vector representation.

[0136] In some embodiments, the resulting descriptor can be further processed by excluding dimensions (positions) that are constant across all data points, for example, if an initial constant stem of defined length appears in all sequences, or if the last coordinate is always zero. This processing can lead to a more compact representation (e.g., reducing dimensionality) while preserving the meaningful variability in folding topology among the sequence set.

[0137] Further details of the Motzkin path descriptor are described in Example 5. The Motzkin-path descriptor can provide a compact, topology-oriented vector representation of a nucleic acid secondary structure that can be readily used for structural similarity comparisons, clustering, or as input features in machine learning analyses.3. Bag of Faces Descriptor

[0138] In some embodiments, a third descriptor, referred to as a Bag of Faces (BoF) descriptor, can be inspired by the “Bag of Words” approach commonly used in natural language processing. Instead of counting the frequency of words in a text, the BoF approach can quantify the frequency of structural features, or “faces,” and their associated energetic configurations in a folded nucleic acid structure. Each face class, such as hairpins, bulges, stacks, interior loops, and multibranch loops, can be assigned an energy value or energy range, and the number ofoccurrences of each class can be counted to form a feature vector. This descriptor can thus capture both the structural composition and energetic profile of the secondary structure, allowing for the classification, clustering, and comparison of large sets of candidate sequences.

[0139] The BoF approach can be conceptually derived from the Bag-of-Words (BoW) model commonly used in natural language processing (NLP). In NLP, BoW representations describe a body of text as a frequency count of word occurrences, discarding syntactic order in favor of statistical features. Analogously, in some embodiments, the BoF descriptor can record the frequency of occurrence of distinct structural face types (e.g., hairpin loop, bulge, internal loop, stack, or multibranch loop) paired with their corresponding energy levels as observed in a folded nucleic acid structure.

[0140] In some embodiments, the BoF descriptor can be implemented by first identifying all the individual faces present in the secondary structure of a sequence. Each face can then be paired with a scalar energy value determined from thermodynamic models (e.g., nearest-neighbor rules). These face / energy pairs can then be used to construct a vector where each element corresponds to a specific face type and energy range combination. For instance, one entry may represent the count of [stack, -1.5 kcal / mol], while another may correspond to [hairpin loop, 2.2 kcal / mol]. In this way, each folded structure can be encoded as a numerical vector capturing both its structural and energetic characteristics.

[0141] In some embodiments, the dimensionality of the BoF vector space can correspond to the number of distinct face / energy configurations observed across a given dataset of folded sequences. For example, a dataset containing many diverse sequences may result in a BoF representation spanning hundreds of unique face / energy combinations. Such high-dimensional vectors can be used for downstream machine learning tasks, including clustering of structurally similar sequences, anomaly detection, or as features in predictive models.

[0142] Further details of the BoF descriptor are described in Example 5. The BoF descriptor can provide a scalable, interpretable, and expressive feature representation of secondary structure configurations, facilitating structure-informed analysis of nucleic acid sequences in computational biology, synthetic biology, and biosensing applications.C. Machine Learning Algorithm of the Disclosure

[0143] Provided herein is an improved computational method for determining minimal free energy structures of nucleic acid sequences using an established energy function framework. Thepresent methods incorporate subgraph matching to identify all admissible interactions, thereby reducing the candidate search space and allowing greater flexibility in optimization approaches.

[0144] This approach offers advantages over existing algorithms by enabling more precise control over structure formation, improving prediction accuracy, and reducing unnecessary computations.

[0145] To determine the lowest energy multibranch configuration containing a given nucleotide pair (i, j), a more comprehensive search process can be conducted compared to previous methods. All possible combinations of up to m edges within the substring [i, j] can be evaluated. This enables finer control over the types of secondary structures produced by the algorithm and ensures that the entire configuration space can be examined, albeit with increased computational cost. Additional details of this procedure can be found in the Examples 10-12.

[0146] Cache requirements can be reduced by performing multibranch optimization, which can eliminate the need for maintaining a secondary cache. As a result, only the structure and energy information for pairs (i, j) that correspond to admissible interactions need to be stored, thereby minimizing overall cache size. The present method can also account for different temperature conditions, as described in Example 9.

[0147] Certain features can be incorporated into future implementations of the present method. For example, the code can assume a salt concentration of 1 M NaCl; adjustments to reference energies can be made for other salt concentrations. Coaxial stacking stabilization energies can be further refined, as current implementations can use simplified or averaged values, whereas experimentally determined, sequence-dependent values have been reported in prior literature. This difference can explain minor discrepancies observed in cases involving complex multibranch structures. In addition, only admissible edges, as defined in the subgraph matching step, are currently included. Non-canonical pairings, which are occasionally predicted by other algorithms, are not presently incorporated but could be added in future versions.

[0148] Provided herein is a method by which machine learning tools can be employed to analyze a large collection of DNA strands derived from the SELEX process. In some embodiments, a specific target-binding problem (e.g., norepinephrine) can be selected and the raw data obtained from next-generation sequencing after SELEX can be processed. Unlike conventional statistical methods, the approach can directly address the global geometry of DNA secondary structures in an unsupervised manner, with the goal of identifying aptamer candidates that may not have high sequence counts yet are structurally similar to high-count aptamers.Analogous pipelines have emerged in related areas of bioinformatics, such as Snekmer algorithm,which develops scalable fingerprinting of protein sequences based on amino acid recoding and hence links sequences with distant similarity.

[0149] Provided herein is newly developed algorithm for high-throughput processing of DNA sequence secondary structures. This algorithm can be directly paired with the SELEX directed evolution method to categorize large numbers of aptamer candidates using modem machine learning methods. The presently described workflow enables comparing, contrasting, and clustering aptamer candidates based on secondary structure and energetics. Namely, it allows sequences with similar structure and energetics to the few sequences with the highest amplification counts from the SELEX process to be identified.

[0150] In some embodiments, the presently described machine learning algorithm can be expanded to include those items mentioned in Section II, A of the present disclosure. The pipeline presented can be ready to be paired with experiments to determine aptamer candidate target-binding affinity and selectivity. In some embodiments, the machine learning methodologies presented herein can be expanded, as additional studies are carried out, to test additional clustering methods.

[0151] In some embodiments, aptamer 3D structures can be incorporated into the machine learning process. This can be more challenging due to the high dependency on the accuracy of the secondary structure information provided. Programs such as AMBER24 use the mechanical force field simulation of the target molecule, usually followed by energy minimization software to do molecular docking analyses. In some embodiments, these types of simulations can be helpful in identifying potential aptamer-target binding sites, which can be experimentally investigated via biophysical techniques, such as, but not limited to, crystallography, nuclear magnetic resonance spectroscopy, and / or cryo electron microscopy (cryo-EM).1. Dataset and Data Pre-Processing

[0152] In some embodiments, a solution-phase systematic evolution of ligands by exponential enrichment (SELEX) process can be performed in accordance with established protocols to select aptamers for a target molecule, such as, but not limited to, a small molecule neurotransmitter. In certain implementations, multiple oligonucleotide libraries can be used, with each library comprising a randomized region of defined nucleotide length (e.g., 48 or 58 nucleotides). The oligonucleotides and associated primers can be standard desalted nucleic acid sequences, and each library can contain identical flanking regions surrounding the randomized segment. In some embodiments, complementary sequences can be placed at the 5’ and 3’ ends to promote stemformation, which can preferentially select for aptamers exhibiting stem closure upon target binding - a feature that can be advantageous in biosensor design.

[0153] In some embodiments, iterative selection rounds can be conducted, followed by semi-quantitative polymerase chain reaction (PCR) to determine when to adjust the stringency of the selection process. Selection rounds can be performed in a buffered solution (e.g., phosphate-buffered saline with divalent cations such as MgCh) under controlled pH and temperature conditions. The PCR parameters can include initial denaturation, cyclic denaturation, annealing, and extension steps, followed by a final extension. Selection stringency can be increased based on comparative band intensities between prewash and wash steps. In certain embodiments, structurally similar molecules (e.g., analogs or competitors) can be used for negative selection to eliminate non-specific binders.

[0154] Once the selection process reaches a point where band intensities in prewash and wash steps show no increase, aptamer samples can be subjected to high-throughput sequencing, such as next- generation sequencing (NGS). In some embodiments, sequencing can be performed after a predefined number of SELEX cycles (e.g., 9th and 13th for a 48-mer library, and 12th and 16th for a 58-mer library). Extended sequencing primers can be used to improve sequencing efficiency.

[0155] In some embodiments, raw sequencing data can be preprocessed to remove primer regions that do not contribute to the structural or functional characteristics of the aptamers.Sequence filtering can involve extracting randomized regions flanked by specific complementary stem-forming segments. While expected lengths correspond to the designed library sizes, sequences of variable length can be retained if they persist through selection, suggesting functional relevance. Sequences can be ranked by read counts, although low-count sequences can also be of interest. In such cases, computational tools can be applied to mine structurally relevant, low-frequency sequences based on similarities to high-frequency sequences using high-throughput structural analysis and clustering.

[0156] In some embodiments, the data cleaning process can result in one or more files containing nucleic acid sequences and associated frequency metrics. Duplicate sequences detected across cycles can be deduplicated by retaining only the highest-count variant.Sequences below a length threshold or containing invalid characters can be excluded. Additional length-based outlier filtering can be applied using percentile cutoffs. The final dataset can consist of thousands of sequences of defined lengths, all of which can be subjected to secondary structure prediction using one or more computational folding algorithms. In some embodiments, base pairing between terminal regions can be enforced to simulate stem formation. Structuralpredictions for thousands of sequences can be performed rapidly using standard computational hardware.2. Training Phase

[0157] A machine learning module as described herein is configured to undergo at least one training phase wherein the machine learning software module is trained to carry out one or more tasks including data extraction, data analysis, and output generation.

[0158] In some embodiments of the presently described machine learning algorithm, the algorithm comprises a training module that trains the machine learning module. The training module is configured to provide training data to the machine learning module, said training data comprising, for example, but not limited to, unique single-stranded DNA sequences obtained from selection (e.g., via SELEX) designed to identify aptamers for the neurotransmitter and hormone norepinephrine. In additional embodiments, the training data is comprised of singlestranded DNA sequences. In some embodiments of a machine learning module described herein, a machine learning module utilizes automatic statistical analysis of data in order to determine which features to extract and / or analyze from single-stranded DNA sequences. In some of these embodiments, the machine learning module determines which features to extract and / or analyze from a single- stranded DNA sequence on the training that the machine learning module receives.

[0159] In some embodiments, a machine learning module can be trained using a data set and a target in a manner that might be described as supervised learning. In these embodiments, the data set is conventionally divided into a training set, a test set, and, in some cases, a validation set. A target is specified that contains the correct classification of each input value in the data set. For example, single- stranded DNA sequences can be repeatedly presented to the machine learning module, and for each sample presented during training, the output generated by the machine learning module can be compared with the desired target. The difference between the target and the set of input samples can be calculated, and the machine learning module can be modified to cause the output to more closely approximate the desired target value. In some embodiments, a back-propagation algorithm can be utilized to cause the output to more closely approximate the desired target value. After a large number of training iterations, the machine learning module output can closely match the desired target for each sample in the input training set. Subsequently, when new input data, not used during training, is presented to the machine learning software module, it can generate an output classification value indicating which of the categories the new sample is most likely to fall into. The machine learning module is said to be able to “generalize” from its training to new, previously unseen input samples. This feature of amachine learning module allows it to be used to classify almost any input data which has a mathematically formulatable relationship to the category to which it should be assigned.

[0160] In some embodiments of the machine learning module described herein, the machine learning module utilizes an individual learning model. An individual learning model is based on the machine learning module having trained on data from a single individual and thus, a machine learning module that utilizes an individual learning model can be configured to be used on a single individual on whose data it trained.

[0161] In some embodiments of the machine training module described herein, the machine training module utilizes a global training model. A global training model is based on the machine training module having trained on data from multiple sources and thus, a machine training module that utilizes a global training model is configured to be used on multiple sources.

[0162] In some embodiments, the use of training models can change as the availability of single- stranded DNA sequences changes. As additional data becomes available, the training model can change to a global or individual model. In some embodiments, a mixture of training models can be used to train the machine training module.

[0163] In some embodiments, unsupervised learning can be used to train a machine training module to use input data such as, for example, single- stranded DNA sequences. Unsupervised learning, in some embodiments, includes feature extraction which can be performed by the machine learning module on the input data. Extracted features can be used for visualization, for classification, for subsequent supervised training, and more generally for representing the input for subsequent storage or analysis. In some cases, each training case can consist of a plurality of aptamer sequences.

[0164] Machine learning modules that are commonly used for unsupervised training can include k-means clustering, mixtures of multinomial distributions, affinity propagation, discrete factor analysis, hidden Markov models, Boltzmann machines, restricted Boltzmann machines, autoencoders, convolutional autoencoders, recurrent neural network autoencoders, and long short-term memory autoencoders. While there are many unsupervised learning models, they all have in common that, for training, they require a training set consisting of biological sequences, without associated labels.

[0165] A machine learning module can include a training phase and a prediction phase. The training phase is typically provided with data in order to train the machine learning algorithm. Data that is inputted into the machine learning module can be used, in some embodiments, to construct a hypothesis function to determine clusters of aptamer sequences whose predicted secondary structure descriptors are structurally and energetically similar to high-count (presumedhigh-affinity) sequences. In some embodiments, a machine learning module is configured to determine if the outcome of the hypothesis function was achieved and based on that analysis make a determination with respect to the data upon which the hypothesis function was constructed. That is, the outcome tends to either reinforce the hypothesis function with respect to the data upon which the hypothesis functions was constructed or contradict the hypothesis function with respect to the data upon which the hypothesis function was constructed. In these embodiments, depending on how close the outcome tends to be to an outcome determined by the hypothesis function, the machine learning algorithm will either adopts, adjusts, or abandon the hypothesis function with respect to the data upon which the hypothesis function was constructed. As such, the machine learning algorithm described herein dynamically learns through the training phase what characteristics of an input (e.g., data) is most predictive in determining whether the features of a given aptamer sequence's secondary structure and energetic configuration are indicative of high affinity or high occurrence within the selection pool.

[0166] For example, a machine learning module is provided with data on which to train so that it, for example, is able to determine the most salient features of nucleotide sequences and their folded secondary structures. The machine learning modules described herein train as to how to analyze structural and energetic configurations of aptamers, rather than analyzing these configurations using pre-defined instructions. As such, the machine learning modules described herein dynamically learn through training what characteristics of an input signal are most predictive in determining whether the features of a folded structure correspond to a likely or functionally significant aptamer conformation.(a) Training Data Construction

[0167] In some embodiments, the training dataset can be derived from high-throughput aptamer selection experiments (e.g., SELEX) in which a library of single-stranded DNA sequences is interrogated for binding to a target molecule. After sequencing (e.g., NGS) of selected pools, each sequence can be associated with a “count” or read-frequency metric reflecting its prevalence in the enriched pool. Sequences can also be filtered to remove invalid nucleotides, duplicates, extremely short fragments, or outlier lengths. In this way, a curated set of unique DNA candidate sequences (e.g., thousands of sequences) can be assembled, each annotated with its observed count metric.

[0168] In some embodiments, each candidate sequence can then be computationally folded using folding algorithms to generate secondary structure predictions. These predictions can include face-energy computations, subgraph matching representations, Motzkin path descriptors,structural-matrix descriptors, or Bag-of-Faces (“BoF”) vectors that encode structural topology and energetic features of each folded sequence. Thus a feature matrix can be generated: rows correspond to individual aptamer candidate sequences, and columns correspond to descriptor features (e.g., structural vector entries, energy bins, face-counts, Motzkin path coordinates).

[0169] In some embodiments, the training labels or target values can be derived from the read-counts (e.g., high-count vs low-count), or other metrics of sequence enrichment or binding affinity (if available). The machine learning module can thus be trained to learn a hypothesis function that maps descriptor vectors (structural / energetic features) to a predictive output (e.g., likelihood of high count / enrichment or favorable binding). During model training, techniques such as dimensionality reduction, clustering, topic modelling, or supervised classification / regression can be applied to capture patterns that distinguish high-enrichment sequences from low-enrichment sequences.

[0170] In some embodiments, data preprocessing steps can include normalization of counts (e.g., converting to percentile ranks), removal of redundant descriptor dimensions (e.g., eliminating descriptor coordinates with constant value across all sequences), and splitting of the dataset into training and validation subsets. The resulting trained model can then be used in the prediction phase to score new or untested candidate sequences based on their structural / energetic descriptors, thereby prioritizing sequences for experimental follow-up that were not initially high-count.3. Prediction Phase

[0171] Following training, the machine learning algorithm is used to determine, for example, the most probable secondary structure configuration of a given nucleotide sequence, based on patterns on which the system was trained using the prediction phase. With appropriate training data, the system can identify energetically favorable and structurally relevant folding patterns. For example, appropriate data derived from known aptamer sequences and their experimentally validated secondary structures can be submitted for analysis to a system using the described trained machine learning algorithm. In these embodiments, a machine learning algorithm can detect structural motifs, such as hairpins, bulges, and multibranch loops, that contribute to the overall folding configuration. In some embodiments, the machine learning algorithm further ranks or scores predicted structures based on learned energetic or functional relevance.

[0172] The prediction phase uses the constructed and optimized hypothesis function from the training phase to predict the probability of a nucleotide sequence adopting a particular secondarystructure and / or any features or metrics computed from the sequence’s inferred folding configuration, such as structural stability, motif presence, or energy minimization outcomes.

[0173] In some embodiments, in the prediction phase, the machine learning module can be used to analyze data derived from nucleotide sequences, such as aptamers or genomic regions, independent of any system or device described herein. In these instances, the new data can provide structural insights or predictive metrics that are required for determining the most probable RNA or DNA secondary structures, stability, or binding potential.

[0174] In some embodiments, a probability threshold can be used in conjunction with a final probability score generated by the machine learning algorithm of the present disclosure to determine whether a given nucleic acid sequence matches a learned structural motif or stability profile. In some embodiments, the probability threshold can be used to tune the sensitivity of the trained machine learning model. For example, the probability threshold can be 1%, 2%, 5%, 10%, 15%, 20%, 25%, 30%, 35%, 40%, 45%, 50%, 55%, 60%, 65%, 70%, 75%, 80%, 85%, 90%, 95%, 98%, or 99%. In some embodiments, the probability threshold can be adjusted if the structural prediction accuracy, sensitivity, or specificity falls below a predefined adjustment threshold. In some embodiments, the adjustment threshold can be used to determine the parameters of the training phase. For example, if the predictive accuracy of the probability threshold falls below the adjustment threshold, the machine learning algorithm of the present disclosure can be configured to extend the training phase and / or require additional labeled structural data. In some embodiments, additional experimentally validated secondary structure measurements and / or sequence data can be included in the training dataset. In some embodiments, such additional data can be used to refine and improve the performance of the machine learning model of the present disclosure.4. Model Architecture

[0175] In some embodiments, the machine learning model architecture can include, but not limited to, a multi-stage pipeline comprising feature extraction, dimensionality reduction, and supervised prediction. Single- stranded DNA sequences can first be transformed into highdimensional structural and energetic descriptors, such as structural matrices, Motzkin paths, and Bag-of-Faces vectors, that capture secondary structure topology and thermodynamics. These descriptors can then be embedded into a lower-dimensional latent space using unsupervised learning techniques (e.g., topic modeling or clustering) to reveal structural patterns. The resulting representations can be input into a supervised learning model (e.g., logistic regression or a neural network) trained to predict whether a given sequence shares features with high-performing aptamers. The architecture can optionally include interpretability and feedback mechanisms to guide experimental validation and iterative refinement.(a) Feature Extraction Module

[0176] In some embodiments, the presently described machine learning model architecture can begin with a feature-extraction module designed to systematically convert raw singlestranded DNA sequence data into structured, vectorized formats that capture biologically meaningful information. This module can implement multiple transformation pipelines in parallel or sequentially to derive complementary perspectives on each input sequence’s structural properties. For example, a structural matrix descriptor can encode the topological layout of base pair interactions as an adjacency matrix, omitting primary structural bonds to emphasize higher- order folding. Simultaneously, a Motzkin-path vector can provide a condensed representation of the secondary structure’s nested or non-crossing base pairings using a balanced parenthetical formalism or its equivalent numerical encoding.

[0177] In parallel, a Bag-of-Faces (BoF) vector can be generated to characterize the energetic landscape of the folded structure by counting the frequency of face-energy configurations such as stacks, bulges, hairpin loops, and multi-branch junctions with associated thermodynamic values. These three types of descriptors can be concatenated or processed independently to form a composite feature set that reflects both the geometric and energetic properties of the secondary structure.

[0178] The resulting high-dimensional representation can then be normalized, filtered for redundant or invariant dimensions, and optionally passed through a dimensionality reduction module (e.g., via PC A or autoencoding) to improve computational efficiency and enhance downstream model performance. This feature-extraction process serves as a critical first stage in the presently described machine learning pipeline.(b) Dimensionality Reduction and Embedding Module

[0179] In some embodiments, following feature extraction, the machine learning pipeline can apply a dimensionality-reduction and embedding module to convert high-dimensional descriptor vectors into a more compact, informative latent space. This step can be essential for reducing noise, improving computational efficiency, and facilitating model generalization, especially when working with large-scale datasets of single-stranded DNA aptamer sequences. The dimensionality -reduction module can use unsupervised learning algorithms, such as NonnegativeMatrix Factorization (NMF), Principal Component Analysis (PCA), or spectral clustering methods, to uncover intrinsic structures within the data without the need for labeled outputs.

[0180] These unsupervised methods can reveal latent variables, sometimes referred to as “topics” or “clusters”, that correspond to frequently occurring structural or energetic patterns in the dataset. For example, sequences exhibiting similar stacking arrangements, hairpin configurations, or multibranch loop structures with consistent energy profiles can be grouped into a common latent cluster. Such groupings can offer biological interpretability, potentially aligning with known structural families or functional motifs.

[0181] In some embodiments, these low-dimensional embeddings can also act as compressed feature representations that are passed directly into downstream prediction engines, such as classifiers or regression models, to enable structure-based predictions (e.g., binding affinity, specificity, or stability). The learned embeddings can thus serve dual roles: facilitating interpretation of the molecular diversity in the sequence dataset and acting as a critical interface between raw structural descriptors and predictive modeling components within the presently described machine learning algorithm framework.(c) Supervised or Semi-Supervised Prediction Module

[0182] In some embodiments, following the embedding of high-dimensional structural descriptors into a latent space, a supervised or semi-supervised prediction module can be employed to assess the likelihood that a given sequence exhibits desired functional characteristics, such as target affinity, specificity, or prevalence in SELEX-derived pools. This prediction module can accept as input either the original descriptor vectors or the lowerdimensional embeddings produced by the prior module and apply a hypothesis function to generate a probabilistic or categorical output.

[0183] The hypothesis function used in this stage can be any number of conventional machine learning models, such as logistic regression, support- vector machines (SVMs), decision trees, ensemble models (e.g., random forests or gradient boosting), or shallow neural networks. These models can be trained using labeled examples, such as sequences experimentally confirmed to exhibit high binding affinity, or semi-supervised using partially labeled datasets, where only a subset of sequences are annotated.

[0184] In some embodiments, the prediction output can be a continuous score representing the probability that the sequence belongs to a class of high-affinity aptamers or is structurally similar to high-frequency sequences identified during the SELEX enrichment process. Alternatively, it can yield binary or multi-class predictions. This predicted score can then be used to rank or filterlarge pools of candidate sequences, effectively narrowing down the search space for experimental validation.(d) Interpretability Layer

[0185] In some embodiments, an interpretability layer can be incorporated into the presently described machine learning pipeline, which can significantly enhance user confidence and insight in the model outcomes. After the downstream prediction module assigns a score or probability to each candidate sequence, the interpretability layer can trace that prediction back to specific structural motifs or graph sub-structures identified in the feature-extraction and embedding stages.

[0186] In some embodiments, this layer can utilize subgraph-matching algorithms to locate instances of stems, loops, hairpins, bulge regions, or multi-branch junctions in the folded structure that correspond to latent “topics” or embedding clusters identified by unsupervised learning. By doing so, the system can reveal which particular motifs (e.g., a conserved stacking region or an unusual bulge configuration) drove a sequence’s high-ranking status.

[0187] In some embodiments, the interpretability layer can output not only the matched motif locations (node indices, nucleotide positions) but also provide comparative visualizations or heat-maps indicating how the candidate’s structure aligns with known high-performing aptamer structures. This enables users to examine, for each high-scoring sequence, the specific base-pair interactions or structural features that the model leveraged, thereby supporting informed decisions about which sequences to synthesize, experimentally validate, or discard. It can also facilitate hypothesis generation; for instance, users might identify a novel stem-loop pairing motif that recurs across multiple high-scoring sequences, suggesting a functional role for that motif in target binding.

[0188] Additionally, in some embodiments, the interpretability layer can support active feedback and model refinement. If experimental binding results are subsequently obtained, the layer can highlight discrepancies between predicted structural motifs and actual binding outcomes, enabling users to refine the descriptor definitions, adjust embedding parameters, or incorporate new structural templates into the subgraph-matching database. Thus, the interpretability layer can serve as a bridge between algorithmic predictions and empirical validation, enabling a closed-loop workflow in which structural insights gleaned by the presently described machine learning algorithm can guide experimental aptamer development and inform subsequent rounds of computational modelling.(e) Feedback Loop

[0189] In some embodiments, the architecture of the presently described machine learning can be enhanced by incorporating a feedback loop that links computational predictions with empirical validation. This feedback loop can operate by experimentally testing the top-ranked candidate sequences identified by the prediction module, such as, but not limited to, through high-throughput aptamer binding assays, surface plasmon resonance, fluorescence assays, or other biochemical evaluation techniques. The experimental outcomes (e.g., binding affinities, specificity profiles, false positives) can then be used as labeled data to refine the machine learning pipeline.

[0190] In some embodiments, sequences confirmed to have high binding affinity can be labeled as true positives, while those failing experimental validation can be labeled as false positives. This labeled dataset can be incorporated into supervised re-training of the model, allowing the hypothesis function to update its weights or decision boundaries in a data-driven manner. In some embodiments, feedback from experimental results can also be used to revise the feature extraction or embedding processes, for instance, by identifying new structural motifs associated with true binding or by re-weighting feature vectors to emphasize previously underrecognized face / energy combinations.

[0191] In some embodiments, this feedback loop enables the presently described machine learning system to remain adaptable to changes in experimental design or target diversity. As new small molecule targets or biologically relevant conditions are introduced, the system can incorporate those variations into its training data and continuously improve prediction accuracy. In some embodiments, active learning strategies can be employed, wherein the model selects uncertain or borderline candidate sequences for experimental testing to maximize learning efficiency.

[0192] In some embodiments, the inclusion of a feedback loop can transform the presently described machine learning module into a closed-loop discovery platform, one that not only predicts aptamer candidates but also evolves and adapts based on empirical outcomes.5. Visualization

[0193] In some embodiments, data visualization can be used to assist in exploring the configuration space of folded nucleic acid sequences by transforming large amounts of secondary-structure descriptor data into compact visual formats. This transformation can facilitate the identification of structural patterns, relationships, and outlier sequences that might otherwise remain hidden when analyzing only a small subset of data at a time.

[0194] In some embodiments, these visualizations can present two- or three-dimensional representations of vector-valued descriptors (e.g., Motzkin-path vectors, structural matrix descriptors, Bag-of-Faces vectors) derived from the sequence dataset. Such lower-dimensional embeddings can be obtained through dimensionality-reduction techniques that preserve key aspects of the high-dimensional data, such as pairwise similarities or intrinsic manifold geometry.

[0195] In some embodiments, non-linear dimensionality-reduction algorithms (such as, but not limited to, t-distributed stochastic neighbor embedding, or t-SNE) can be applied. These methods can map high-dimensional descriptors to a lower-dimensional space by optimizing the divergence between similarity distributions in the original and embedded spaces. Parameters such as perplexity can be adjusted to balance focus between local and global data structures, smaller values can emphasize nearby neighbors, while larger values can highlight broader structure.

[0196] In some embodiments, the selection of dimensionality-reduction parameters (such as perplexity) can critically influence the resulting visual layout. It can be noted that algorithms such as t-SNE can scale quadratically with dataset size, which can limit applicability when analyzing very large datasets. On moderate-sized datasets (e.g., a few thousand sequences), computations can complete in seconds on standard hardware.

[0197] In some embodiments, colors or other visual cues can be applied to each data point based on additional sequence metadata (e.g., read counts from NGS). Clusters of sequences with similar descriptor vectors can then become visually apparent, high-count sequences can cluster in distinct regions, while low-count sequences may appear in other regions or overlap with the high-count clusters. Based on this, in some embodiments, one can identify sequences with similar structural or energetic configurations to top-performing candidates, even when their counts are low, thereby guiding subsequent analysis or experimental validation.

[0198] More detailed information regarding visualization search are described in Example 6.6. Similarity Search

[0199] In some embodiments, a structural and energetic similarity search can be performed on a dataset of single-stranded DNA sequences in order to identify sequences whose secondary structures share both the topological and energetic configurations of high-count aptamer candidates, thereby enabling the discovery of potentially high-affinity sequences beyond those with the highest sequencing counts.

[0200] In some embodiments, when analyzing large datasets of candidate sequences from processes such as SELEX, the objective can be to identify those sequences whose folded structures resemble those hypothesized to exhibit favorable binding properties (for example,those associated with high sequencing counts). In some embodiments, two distinct similarity domains can be considered: structural similarity, which relates to the topology of the folded sequence, and energetic similarity, which relates to the energetic configuration of the secondary structure, such as the energy contributions of each structural “face”.

[0201] In some embodiments, Motzkin-path descriptors or structural-matrix descriptors can be used to identify sequences that are topologically equivalent to a given reference sequence. Concurrently, Bag of Faces (BoF) descriptors can be used to identify sequences that share similar energetic configurations with the reference structure.

[0202] In some embodiments, D =denote the set of all single-stranded DNA sequences under investigation, M ■■= {mj}=1denote the associated Motzkin-path descriptor vectors, and F ■■=denote the associated BoF descriptor vectors, where mi and ftare the Motzkin path and BoF representations of the sequences xt, respectively. Further, x denote a reference sequence (for example, a high-count sequence) with descriptors m and f. For nonnegative thresholds emand e^, define:= {xLE D such that ||mj — m||2<:={xtE D such that \\ft— f|| < e^}. In these embodiments, T^oand 8^0can represent the subsets of sequences that are structurally and energetically equivalent (respectively) to the reference sequence x.

[0203] In some embodiments, the similarity search can be applied to a dataset containing thousands of sequences. One goal of this search can be to identify sequences that share both topological and energetic features with sequences whose counts fall in the top 0.1 percentile of the dataset (z.e., whose count exceeds 99.9% of all sequences). To ensure structural diversity, only a selected subset of top-count sequences (for instance, the first, second, fourth and fifth highest-count sequences) may be studied as references; sequences sharing a topology identical to the third highest-count sequence may be excluded if they duplicate the topological configuration of the highest-count reference.

[0204] In some embodiments, for each reference sequence xt, the cardinality of the structural-equivalence set \T(i 0) | can be determined. In a non-limiting example, one reference topology could correspond to 103 matching sequences, another to 59, yet another to 34, and another to 36. Thereafter, a further analysis of energetic profiles within each structural set can reveal how many distinct energetic variants exist (e.g., nine distinct energy profiles for one topology, three for another, seven for a third, and eleven for a fourth). Within those profiles, striped bars in a bar-plot can indicate profiles matching high-count reference sequences; for instance, 80 sequences may share both topology and energy with the first reference, 56 with the second, 26 with the third, and 20 with the fourth.

[0205] In some embodiments, visualization of the results can reveal that only a subset (for example, fewer than 11%) of the dataset sequences fall within clusters containing high-count references, thereby suggesting that further research efforts can be focused on this promising subset while deprioritizing the remaining majority of sequences.

[0206] More detailed information regarding similarity search are described in Example 7.7. Topic Modeling and Spectral Clustering

[0207] In some embodiments, sequences can be clustered based on the energetic configurations of their secondary structures, where latent topic models are built usingBag-of-Faces (BoF) descriptors and the resulting topic distributions are used to describe and cluster the sequences.

[0208] In some embodiments, topic modeling can be used as an unsupervised machine learning technique to discover hidden structural or energetic motifs within a corpus of sequences, by analyzing occurrences of structural elements (analogous to “words” in text) across all sequences. In some embodiments, unsupervised methods such as Latent Dirichlet Allocation (LDA) or Nonnegative Matrix Factorization (NMF) can automatically associate each sequence (analogous to a “document”) with a distribution over latent topics (analogous to “topics” in NLP).

[0209] In some embodiments, to perform topic modeling, NMF can be used: a linear algebraic method advantageous for computational efficiency and scalability to large datasets. In some embodiments, the model can take as input a matrix X E Blvxt(where w is the number of distinct face-energy configurations and t is the number of sequences). In some embodiments, the matrix can be decomposed into two nonnegative matrices W E Blv,kand H E IKk,t, such that: X ~ W • H. The row of W can represent basis vectors (topics) and the columns of H can encode the contribution of each topic for each sequence. The hyperparameter k (number of topics) can be defined a priori by the user and the algorithm can be applied to compute W and H via optimization.

[0210] In some embodiments, after topic modeling, clusters of sequences with similar topic -mixture distributions can be identified using clustering techniques (e.g., spectral clustering). A similarity graph can be constructed from the topic -mixture vectors, where each sequence is a node and edges encode similarity of topic distributions. Spectral clustering can leverage the eigenvectors of the graph Laplacian to partition the sequences into distinct clusters. The number of clusters can be defined a priori in line with the number of topics or other design criteria, facilitating grouping of sequences that share structural / energetic motif composition.

[0211] In some embodiments, visualization (e.g., via t-SNE) of the resulting clusters can demonstrate that high-dimensional descriptor similarity is well represented in low-dimensional space, enabling identification of structural families and guiding selection of candidate sequences for further experimental validation.

[0212] More detailed information regarding topic modeling and spectral clustering are described in Example 8.III. Computer Implemented System of the Disclosure

[0213] In various embodiments, the high throughput systems and methods for determining, classifying, and assessing the functionality of structures of aptamers can be implemented via computer software or hardware. Refer to the Examples 1-12 for further information regarding the system, devices and methods provided herein, in accordance with various embodiments.

[0214] FIG. 3 is a block diagram illustrating a computer system 300 upon which embodiments of the present teachings may be implemented. In various embodiments of the present teachings, computer system 300 can include a bus 302 or other communication mechanism for communicating information and a processor 304 coupled with bus 302 for processing information. In various embodiments, computer system 300 can also include a memory, which can be a random-access memory (RAM) 306 or other dynamic storage device, coupled to bus 302 for determining instructions to be executed by processor 304. Memory can also be used for storing temporary variables or other intermediate information during execution of instructions to be executed by processor 304. In various embodiments, computer system 300 can further include a read only memory (ROM) 308 or other static storage device coupled to bus 302 for storing static information and instructions for processor 304. A storage device 310, such as a magnetic disk or optical disk, can be provided and coupled to bus 302 for storing information and instructions.

[0215] In various embodiments, computer system 300 can be coupled via bus 302 to a display 312, such as a cathode ray tube (CRT) or liquid crystal display (LCD), for displaying information to a computer user. In some embodiments, the display 312 may be enable user input (e.g., via a touchscreen). For example, the display 312 may enable a user (e.g., a research scientist) to enter a target antigen, protein, and / or molecule of interest for which information concerning a binding characteristic of aptamers is desired by touching the display 312. An input device 314, including alphanumeric and other keys, can also enable the same and can be coupled to bus 302 for communication of information and command selections to processor 304. Another type of user input device is a cursor control 316, such as a mouse, a trackball or cursor direction keys forcommunicating direction information and command selections to processor 304 and for controlling cursor movement on display 312. This input device 314 typically has two degrees of freedom in two axes, a first axis (z.e., x) and a second axis (z.e., y), that allows the device to specify positions in a plane. However, it should be understood that input devices 314 allowing for 3-dimensional (x, y and z) cursor movement are also contemplated herein.

[0216] Consistent with certain implementations of the present teachings, results can be provided by computer system 300 in response to processor 304 executing one or more sequences of one or more instructions contained in memory 306. Such instructions can be read into memory 306 from another computer-readable medium or computer-readable storage medium, such as storage device 310. Execution of the sequences of instructions contained in memory 306 can cause processor 304 to perform the processes described herein. Alternatively, hard-wired circuitry can be used in place of or in combination with software instructions to implement the present teachings. Thus, implementations of the present teachings are not limited to any specific combination of hardware circuitry and software.

[0217] The term “computer-readable medium” (e.g., data store, data storage, etc.) or “computer-readable storage medium” as used herein refers to any media that participates in providing instructions to processor 304 for execution. Such a medium can take many forms, including but not limited to, non-volatile media, volatile media, and transmission media.Examples of non-volatile media can include, but are not limited to, dynamic memory, such as memory 306. Examples of transmission media can include, but are not limited to, coaxial cables, copper wire, and fiber optics, including the wires that comprise bus 302.

[0218] Common forms of computer-readable media include, for example, a floppy disk, a flexible disk, hard disk, magnetic tape, or any other magnetic medium, a CD-ROM, any other optical medium, punch cards, paper tape, any other physical medium with patterns of holes, a RAM, PROM, and EPROM, a FLASH-EPROM, another memory chip or cartridge, or any other tangible medium from which a computer can read.

[0219] In addition to computer-readable medium, instructions or data can be provided as signals on transmission media included in a communications apparatus or system to provide sequences of one or more instructions to processor 304 of computer system 300 for execution. For example, a communication apparatus may include a transceiver having signals indicative of instructions and data. The instructions and data are configured to cause one or more processors to implement the functions outlined in the disclosure herein. Representative examples of data communications transmission connections can include, but are not limited to, telephone modemconnections, wide area networks (WAN), local area networks (LAN), infrared data connections, NFC connections, etc.

[0220] It should be appreciated that the methodologies described herein, flow charts, diagrams and accompanying disclosure can be implemented using computer system 300 as a standalone device or on a distributed network or shared computer processing resources such as a cloud computing network.

[0221] The methodologies described herein may be implemented by various means depending upon the application. For example, these methodologies may be implemented in hardware, firmware, software, or any combination thereof. For a hardware implementation, the processing unit may be implemented within one or more application specific integrated circuits (ASICs), digital signal processors (DSPs), digital signal processing devices (DSPDs), programmable logic devices (PLDs), field programmable gate arrays (FPGAs), processors, controllers, microcontrollers, microprocessors, electronic devices, other electronic units designed to perform the functions described herein, or a combination thereof.

[0222] In various embodiments, the methods of the present teachings may be implemented as firmware and / or a software program and applications written in conventional programming languages such as C, C++, Python, etc. If implemented as firmware and / or software, the embodiments described herein can be implemented on a non-transitory computer-readable medium in which a program is stored for causing a computer to perform the methods described above. It should be understood that the various engines described herein can be provided on a computer system, such as computer system 300, whereby processor 304 would execute the analyses and determinations provided by these engines, subject to instructions provided by any one of, or a combination of, memory components 306 / 308 / 310 and user input.

[0223] In describing the various embodiments, the specification may have presented a method and / or process as a particular sequence of steps. However, to the extent that the method or process does not rely on the particular order of steps set forth herein, the method or process should not be limited to the particular sequence of steps described. As one of ordinary skill in the art would appreciate, other sequences of steps may be possible. Therefore, the particular order of the steps set forth in the specification should not be construed as limitations on the claims. In addition, the claims directed to the method and / or process should not be limited to the performance of their steps in the order written, and one skilled in the art can readily appreciate that the sequences may be varied and still remain within the spirit and scope of the various embodiments. Similarly, any of the various system embodiments may have been presented as a group of particular components. However, these systems should not be limited to the particularset of components, nor their specific configuration, communication and physical orientation with respect to each other. One skilled in the art should readily appreciate that these components can have various configurations and physical orientations (e.g., wholly separate components, units and subunits of groups of components, different communication regimes between components).

[0224] Although specific embodiments and applications of the disclosure have been described in this specification, these embodiments and applications are exemplary only, and many variations are possible.EXAMPLES

[0225] These examples are provided for illustrative purposes only and not to limit the scope of the claims provided herein.EXAMPLE 1. Subgraph Matching and Free Energy Minimization for DNA Folding

[0226] This Example describes a non-limiting method on an elimination scheme built as follows:

[0227] Step 1: The DNA strand was represented as a linear graph with edges defined between adjacent nodes corresponding to the primary sequence. A mapping was performed between the strand and a duplicate of itself according to specified rules. For the purposes of graph matching, nucleotide bases (adenine (A), thymine (T), cytosine (C), and guanine (G)) were used as labels on the graph nodes. These labels were assigned in the order corresponding to the primary DNA sequence.

[0228] A list of candidate nodes in the target graph to which every other node could be matched was created. The rules for matching were defined by canonical nucleotide conjugate base pairs. Specifically, for each adenine (A) in the primary sequence, the matching candidates were identified as all nodes labeled with thymine (T). Similarly, for each guanine (G), the candidate list included all nodes labeled with cytosine (C). The list was completed by identifying all candidate nodes for thymine (T) and cytosine (C). At the conclusion of this step, a set of possible matched candidate nodes was established. These matched nodes were then systematically eliminated, as described below, until an optimal structure was obtained. An example of the output from Step 1 is shown in FIG. 4 (left panel).

[0229] Step 2: This step served as an analog of the “topology” filter described in Moorman et al., IEEE Trans. Netw, Sci. Eng. 8(2): 1367- 1384 (2021). A constraint was applied to eliminate matched candidate nodes that did not connect to nodes matching the first-order connections in thetemplate graph, as defined in the subgraph matching problem. For DNA aptamer structures, the search was conducted based on stack configurations, as illustrated in FIG. 2B. All matched candidate nodes that were not part of a stack - defined by two adjacent matched nodes - were eliminated. In FIG. 2B, two primary (adjacent) edges along the oligonucleotide backbone are shown in black, and two secondary (base pair) edges are shown in red. In this Example, nodes 1 and 2 were conjugate to nodes 3 and 4, respectively. In the following discussion, the secondary edges that remained after Step 2 filtering are referred to as “admissible” edges. FIG. 4 (right panel) illustrates the secondary edges that remained following execution of Step 2 for a serotonin aptamer. As shown in FIG. 5, the filtering step removed a little over half of the edges from the graph.

[0230] Step 3: This step involved eliminating additional matching candidate nodes to achieve a minimum free energy configuration, with at most one candidate pairing per nucleotide. This process required an understanding of face configurations and their associated energies. Five face classes were considered: hairpin loops, stacking regions, bulge loops, interior loops, and multibranch loops, as illustrated in FIG. 2C. This step could also be framed within the context of subgraph matching, wherein the face classes provided a natural subgraph grouping of matched nodes to simplify the search for the lowest energy configuration. Each face was assigned a scalar energy value, and the total energy was computed as the sum of all face energies. The objective was to identify the configuration with the minimal total energy by selecting among the remaining candidate edges between different matched nodes (nucleotide pairs). It should be noted that the energy depended on both the face configuration and the specific nucleotides composing each face.

[0231] Step 3 consisted of a comprehensive search over substrings derived from the DNA strand. For each candidate interaction (i, j), all possible face energies associated with any admissible face containing the pair (i, j) were computed. The search was initiated with smaller distances between i and j and was extended progressively to evaluate larger distances. The algorithm concluded by selecting the edge (i, j) with the minimal internal energy and designating it as the “last pair,” defined such that no faces existed outside of the corresponding substring. All nucleotide pairs included in the optimal configuration were then selected as part of the final predicted structure.

[0232] Optimization involving hairpin loops, interior loops, and bulge face energies was relatively straightforward compared to optimization required for multibranch structures. For multibranch configurations, the complexity arose from the need to select among all possiblesubsets of allowable edges, which constituted a combinatorially large search space. Further details regarding this optimization process are provided in the Examples 10-12.

[0233] The recursive process used to identify the minimum free energy configuration, as originally proposed by Zuker and Stiegler, Nucleic Acids Res. 9(1): 133- 148 (1981), is outlined herein. This process forms the basis of the algorithms employed by tools such as mfold and SeqFold. In this approach, structural data were stored in structure caches, Svand S which functioned as mappings from index pairs (i, j) to corresponding energy values and associated structural configurations. These caches were implemented as lists of lists, where each element contained an energy value and an interaction configuration. The cache Sv(i, j) represented the optimal energy and structure of the substring spanning positions [i, j], under the assumption that i and j formed an internal base-pair interaction. In contrast, Sw(i, j) represented the optimal energy and structure of the same substring without any interaction assumptions. In cases where i and j interacted in the optimal configuration, it followed that Sv(i, j) = Sw(i, j).

[0234] The recursive process was defined as follows.Formula (1) where Sv(i, j) must be computed explicitly using the values of the cache for all subsequences of [i, j], The method used to fill the Svcache was by computing the minimal energy that the string [i, j] can have if it ends in a face, as described in Example 2.EXAMPLE 2. Energy Optimization and Algorithmic Complexity

[0235] This Example describes the recursive method used to compute the minimum free energy configuration of DNA secondary structures by systematically evaluating face classes, such as hairpin loops, inner loops, bulges, stacks, and multibranches (FIG. 2C), within substrings of the DNA sequence.

[0236] In each step of the energy minimization algorithm, a decision was made regarding the type of face to be formed and its corresponding size. It was assumed that in the substring [i, j], a face would terminate at (i, j). If this face included no internal interactions, it was classified as a hairpin loop. If it included a single internal interaction, it was considered an inner loop, bulge loop, or stack. If two or more internal interactions were present, the structure was identified as a multibranch loop. To determine the optimal face configuration, it was generally necessary to identify the optimal arrangement of internal base pairs. The following describes the systematic process used to compute these structures and their associated energies.

[0237] Under the assumption that i interacted with j, the face structure with the lowest energy was determined by computing optimal energies for each of three face types: the hairpin loop, allpossible inner loop / bulge / stack combinations, and all feasible multibranch structures. For the hairpin structure, only one configuration existed, and the complexity of this calculation was sublinear in time.

[0238] The inner loop, bulge, and stack structures were grouped together, as each was defined by exactly two interior base pairs. A search was conducted over all admissible edges - identified in Step 2 - within the substring (i, j) to determine the face with the lowest associated energy. The computational complexity of this step was generally 풪(|j - i|2); however, in theimproved SeqFold and the presently described improved algorithm implementations, the number of evaluations required was significantly reduced, fewer than O(n²) checks, due to prior edge pruning via graph matching, as illustrated in FIG. 4.

[0239] Multibranch optimization involved evaluating all combinations of admissible base pairs within the substring [i, j] that could potentially form a multibranch structure. This step represented a combinatorially large search space. An upper bound on the search space was defined as 2N, where N represented the number of admissible base pairs within [i, j], In earlier work (Zuker and Stiegler, Nucleic Acids Res. 9(1): 133-148 (1981)), this was avoided by assuming the existence of a k E [i, j] such that the optimal multibranch configuration could be formed by joining the optimal structures of [i+1, k] and [k, j-1], Using this assumption and storing both the interaction and overall optimal structures, algorithms such as those implemented in SeqFold, mfold, and the presently improved SeqFold identified optimal multibranch configurations in O(\j - i\) steps.

[0240] In this Example, to gain greater control over the search process, a parameter m was introduced to limit the number of branches permitted in a multibranch structure. This reduced the complexity from evaluating all combinations of edges between i and j to evaluating only those combinations involving at most m edges. The resulting computational requirement was at most 풪 which was further reduced following the edge filtering conductedduring the graph matching step (Step 2).

[0241] Combining all of these complexities, the complexity of computing the optimal face forfor the algorithms based on the Zuker-Steigler approach, andfor the presently described improved algorithm.

[0242] This free energy optimization described herein was carried out for all admissible edges (z, j). An estimate for the number of admissible edges was not available; therefore, a naive upperbound based on all possible pairs (i, j) was used. This gives upper bounds for the overall complexity offor the Zuker-Steigler algorithms, andfor the presently described improved algorithm.

[0243] These estimates included all possible terms, many of which may not have been necessary. In the calculations, the parameter m was set to 4. For experimental runtime comparisons, see FIG.6, which illustrates the performance of SeqFold and the presently improved SeqFold. Notably, the observed computational scaling for SeqFold was 풪(n3.696), whereas the presently improved SeqFold exhibited a scaling of 풪(n3.477).EXAMPLE 3. Numerical Results for DNA Folding

[0244] This Example describes computational benchmarks comparing the performance of the present approach with existing algorithms for predicting minimum free energy (MFE) secondary structures of aptamers.

[0245] Computational examples were presented to illustrate the advantages of the proposed approach. Comparisons were made with two benchmark algorithms for MFE secondary structure prediction: SeqFold and mfold. Several benchmark sequences, including those from SeqFold’ s comparison library and actual aptamers with known targets, were used to evaluate differences in algorithm performance, as shown in this section and in Example 11.

[0246] The analysis began with five aptamer sequences from the SeqFold library (Table 1), previously used to validate SeqFold against mfold. From the original set of seven DNA sequences, only those where SeqFold and mfold disagreed were considered. For each sequence, secondary structures were computed using the present approach, the presently improved SeqFold, and the two baseline methods, mfold was used to determine the ground truth energies of the resulting folded structures. For sequences SI, S2, and S3, the MFE structures computed by the present method matched those of mfold in both energy and structure, suggesting that subgraph matching could be effectively employed.Table 1. Results from experiments conducted on five single-stranded DNA sequences.* These sequences were among those used by the authors of SeqFold to validate their approach against mfold and can be found in the SeqFold GitHub repository. The table presented the structure energies, computed using mfold, that were associated with the MFE structures predicted by the present approach, the presently improved SeqFold, SeqFold, and mfold, respectively. The lowest energy values were indicated in bold.

[0247] In contrast, sequences S4 and S5 provided examples where discrepancies were observed. The structures generated by the present approach differed from those predicted by mfold, as evidenced by higher associated energies. Despite these differences, the present method consistently outperformed SeqFold across all tested sequences. For SI, S2, and S3, SeqFold failed to compute the MFE structures identified by mfold.

[0248] Structural differences were further illustrated through additional examples. FIG. 7 depicted a longer aptamer from the SeqFold catalog, where the present method achieved a closer approximation to mfold than either SeqFold or the presently improved SeqFold, likely due to differences in multibranch energy computation. Specifically, the presently described machine learning algorithm did not find an interaction between positions 34 and 42 in one of the multibranches and found an interaction between positions 71 and 79 in the other multibranch. The presently improved SeqFold failed to find the additional branch on the external loop, and found different lengths of the branches in the multibranch loop, again due to a slight difference inhow the presently described machine learning algorithm and the presently improved SeqFold calculate multibranch energies.

[0249] FIG. 8 demonstrated that, unlike SeqFold and the presently improved SeqFold, the present approach identified a coaxial stacked multibranch structure predicted by mfold.

[0250] FIG. 9 depicts an aptamer for Sgc-3b, a particular membrane protein. Similar to the example in FIG. 7, the presently described machine learning algorithm found an aptamer most similar to what mfold predicts. The only discrepancy was the initial stacking region that the presently described machine learning algorithm found was one shorter than what mfold found. Both SeqFold and the presently improved SeqFold predicted entirely different multibranches. Using the mfold energies, SeqFold and the presently improved SeqFold showed much less favorable energies.

[0251] FIG. 10 depicts an aptamer for theophylline. This example illustrates noncanonical pairing, which was not found by the presently described machine learning algorithm or the presently improved SeqFold, due to Step 2. The interaction (12, 32) is a T-G non-canonical pairing found by mfold. This caused the G - C pair right before it to be deleted by Step 2. This caused the presently improved SeqFold and the presently described machine learning algorithm to neglect the central stacking region, as neither algorithm finds having the G-C, G-C stack to be favorable without the rest of that stretch. More examples comparing the secondary structures of various exemplary aptamers are shown in Example 11.EXAMPLE 4. Dataset and Data Pre-processing for Machine Learning

[0252] This Example describes the implementation of a solution-phase selection procedure and subsequent high-throughput sequence and structure analysis for aptamers targeting a small-molecule neurotransmitter.

[0253] Solution phase systematic evolution of ligands by exponential enrichment (SELEX) was performed according to previously published protocols to select new aptamers for the small molecule neurotransmitter target. Nakatsuka et al., Science 362(6412):319-324 (2018); Yang et al., Science 380(6648):942-948 (2023); Yang et al., J. Am. Chem. Soc. 134(3): 1642-1647 (2012); Yang et al., ACS Chem. Biology 12(12):3103— 3112 (2017); Yang et al., Methods 106:58-65 (2016); Cheung et al., ACS Sensors 4(12):3308-3317 (2019; Wang et al., Sci. Adv.8(1):eabk0967 (2022). Two libraries were used. One library had a 48-nucleotide randomized region and the other had a 58-nucleotide randomized region. Standard desalted oligonucleotides were used for both libraries, as well as their primers. Both libraries contained identical sequences flanking the randomized regions. The five nucleotides on the 5’ and 3’ ends werecomplementary to one another. As previously reported in Yang et al., J. Am. Chem. Soc. 134(3), this library design preferentially selected for aptamers that underwent stem closure upon target binding, a design approach advantageous for biosensing applications.

[0254] Iterative selection rounds were followed by semi-quantitative PCR to determine when to increase the stringency of the selection conditions. Selections were carried out in phosphate-buffered saline (PBS) with 2 mM MgCl2at pH 7.4. The PCR started at 95 °C for 2 minutes, followed by N cycles of [95 °C for 15 seconds, 60 °C for 20 seconds, 72 °C for 30 seconds] and a final cycle of 72 °C for 3 minutes. Each PCR run lasted 13 + 2 cycles.Generally, target concentrations were decreased when the band densities of the prewash in the previous step and the first wash in the next step were similar, and the bands remained visible. Negative selection was performed using dopamine, serotonin, and epinephrine, all of which were structurally similar to the target molecule.

[0255] Once the prewash and first target wash repeatedly showed no increase in band densities, the oligonucleotide samples were sent for amplicon-EZ next-generation sequencing (NGS). Samples from the 9th and 13th SELEX cycles were sequenced for the 48-nucleotide randomized region library. Samples from the 12th and 16th cycles were sequenced for the 58-nucleotide randomized region library. Longer overhang NGS primers were added to increase sequencing speed.

[0256] The NGS data were preliminarily cleaned by removing the primer sequences from both the 3’ and 5’ ends, as they did not contribute to aptamer function, tertiary structure, or sensing performance. The data were filtered to focus on the randomized regions by selecting nucleotides within the five outermost complementary nucleotides on the 5’ and 3’ ends, which formed the stem. Although all sequences were expected to have lengths corresponding to the 48 or 58 random region libraries, other lengths were observed, likely due to insertions or PCR errors. These variable-length sequences were retained because they persisted in the screening pool, indicating potential target recognition. All sequences were ranked by the number of reads (“counts”) in the NGS files. Typically only high-count sequences would be selected for empirical analysis of target recognition and sensing performance. However, this practice risked excluding high-affinity, high- specificity aptamers with low counts. New algorithms were therefore developed to mine low-count sequences based on structural similarity with high-count sequences via high-throughput secondary structure analysis and clustering.

[0257] The data cleaning process resulted in four files containing single-stranded DNA sequences and associated counts. Two files contained sequences from the 48-nucleotide randomized region library (after 9 and 13 PCR cycles) and two files contained sequences fromthe 58-nucleotide randomized region library (after 12 and 16 PCR cycles). All sequences arising from the data generation and cleaning procedures were considered, totaling 16,631single-stranded DNA sequences with associated counts. These sequences included the five complementary nucleotides on the 5’ and 3’ ends, which formed the terminal stem. Duplicate sequences (ones reported in files associated with different cycles) were removed, retaining the version associated with the highest count. Sequences shorter than 20 nucleotides or those containing an invalid nucleotide (i.e., not A, T, C, or G) were also removed. After filtering, a dataset of 4,933 sequences remained. Additional outlier sequences based on length were removed by retaining only those between the 5th and 95th percentiles of length. The final dataset contained 4,450 sequences with lengths ranging from 33 to 83 nucleotides. Secondary structures and face energies of all sequences in the processed dataset were computed using the described pipeline. The five complementary nucleotides on the 5’ and 3’ ends were forced to base pair. Folding of all 4,450 sequences took approximately three minutes on a commercial laptop equipped with an 11th Gen Intel Core i5- 1135G7 CPU at 2.42 GHz with 8 GB RAM.EXAMPLE 5. Mathematical Representation of Secondary Structures

[0258] This Example describes the use of multiple vector-valued descriptors to represent the secondary structure of single-stranded DNA sequences.

[0259] Three vector-valued descriptors were proposed to represent the secondary structure of single-stranded DNA sequences. The first descriptor was based on the adjacency matrix of the secondary structure graph, which was termed the structural matrix. The second descriptor was based on Motzkin paths, combinatorial objects used in mathematics, particularly in the study of lattice paths. The third descriptor type was inspired by the “Bag of Words” approach, commonly used in text-mining and natural language processing.Structural Matrix Descriptors

[0260] The structural matrix descriptors consisted of the adjacency matrix of the secondary structure graph, where all matrix entries corresponding to the primary structural bonds between nucleotides were set to zero (FIG.2D, lower left). These descriptors provided information solely about the topology of the secondary structures. They did not represent information related to the nucleotides composing the single-stranded DNA sequence or the energies of the folded structures.Motzkin Paths Descriptors

[0261] A Motzkin Path of length N was defined in NxN starting at (0, 0) and ending at (0, N) with each step being one of (+1, +1), (+1, 0), or (+1, -1). Such a path could be recorded as a balanced parenthetical sequence with characters, that is, a sequence 5 G {(,.,) }Nthat was balanced in parentheses. The mapping (■-> (+1, +1),. ■-> (+1, 0), ) ■-> (+1, -1) was used to translate between forms. The balanced condition ensured that the path never crossed below the x axis. Since each path step moved one to the right, the sequence could also be recorded as G {(-1, 0, 1)}N, where all partial sums were non-negative. A sequence 5 G I l'"vwas admissible as a partial sum sequence if= sN= 0, and |Sj+1— st< 1| for all i. These objects were known to correspond to the secondary structures of RNA sequences without pseudoknots.

[0262] Motzkin numbers are a generalization of Catalan numbers and have been studied thoroughly, with several identities known (Donaghey et al., Theory Ser. A 23(3):291— 301 (1977)). At the level of two-dimensional structures, DNA and RNA sequences are in bijection under the map U ■-> T,. However RNA also exhibits pseudoknots, a feature not commonly present in DNA. While single- stranded DNA sequences can form pseudo knots, the conditions under which they form are specialized. This means that modified versions of older RNA methods could be used to predict DNA structure efficiently, as described in the algorithms above.

[0263] The quickest heuristic was to look at the balanced parenthetical sequence with spaces, and note that when a pair of parentheses balances, those two indices would interact in the DNA strand. Notably, assuming a lack of pseudoknots, the pairings can be drawn as non-crossing paths in the plane, and thus the secondary structure of a sequence with length N was a noncrossing chording of an IV-gon. For what follows, illustrated herein is that this structure was completely determined by a balanced parenthetical sequence with “,”s inserted, or equivalently, a Motzkin path.

[0264] To go from a collection of chords of an IV-gon to a parenthetical sequence, a starting point of the IV-gon was chosen to hold constant (say the top) and move clockwise, replacing each start of a chord with ‘(’, each position without a chord by and each end of a chord by ‘)’. This results in a sequence with an equal amount of ‘(’ and ‘)’ parentheses. In terms of the DNA sequence, the sequence was laid along the IV-gon in order. When a base pairing was encountered, open and closed parentheses were placed on each side of the pair, and a chord was drawn from one nucleotide to the other. It was noted that the open parenthesis was placed before the close parenthesis in the sequence, and an equal number of both was maintained, as each pair corresponded to a chord, resulting in a balanced parenthetical sequence. As the sequence was balanced, two different chordings would have two different words, as they must differ by at least one chord position. One could also verify that putting the parentheses along the circle andconnecting corresponding open and closed parentheses would give a non-crossing chording of the circle, showing that this mapping was, in fact, a bijection. Note, that a balanced parenthetical sequence with periods could be replaced with a { 1, 0, -1 }wvector whose cumulative sum was non-negative via the substitution ( — 1,. — 0 ) — - 1.

[0265] To visualize how this description provided the secondary sequence of the DNA strand, the backbone was placed clockwise along the outside of the circle, and a chord was drawn between any two nucleotides that interacted in the secondary structure. The non-crossing condition was derived from the observation that DNA sequences tended not to form so-called pseudoknots. An example sequence is shown below FIG. 11.

[0266] Note that all information about the topology of the secondary structure was stored in this Motzkin path formalism. Thus, the Motzkin path vectors (both v E {— 1, 0, 1}Wand w, wt= jl=ovj) could be sued to represent the sequence secondary structures. By comparing the Motzkin path vectors, structural similarities between different folded sequences could be identified. The cumulative sum vector was used to describe the secondary structure of the sequences in our numerical experiments. Moreover, dimensions where the Motzkin path vectors were identical across all data points were excluded.

[0267] This additional processing step eliminated redundant information from the Motzkin path descriptors, resulting in a more compact representation. Six dimensions were removed: five related to the initial stack of five base pairs present in all the secondary structures in the dataset, and one from the fact that all Motzkin path descriptors considered herein have value zero in the last dimension. The processed Motzkin path descriptors had 77 dimensions.

[0268] The Motzkin path and structural matrix descriptors only described the secondary structure’s topology. However, the chemical properties of a folded single-stranded DNA sequence also depended on the secondary structure’s energetic configuration. Next, the Bag of Faces representation is used to characterize the energetic configurations of secondary structures.Bag of Faces Descriptors

[0269] The folded secondary structure of a single-stranded DNA sequence was characterized by the faces of the secondary structure graph and their associated energies. Based on this observation, a descriptor was proposed that counted the occurrences of face / energy configurations in a given folded sequence. This descriptor was referred to as Bag of Faces (BoF).

[0270] The concept of BoF was inspired by the Bag-of-Words descriptor commonly used in natural language processing to provide vectorized representations of text. The Bag-of-Wordsdescriptor encoded the frequency of word occurrences in a given text. In this work, faces and their corresponding energies in the secondary structures were used in an analogous manner to generate Bag of Faces descriptors.

[0271] In the BoF model, a secondary structure was mapped onto a vector consisting of bags, where each bag, or entry of the vector, counted the occurrences of a particular face / energy type, e.g., [stack, -1.5], [hairpin, 2.2], and so on. A visual representation of a BoF descriptor representing the secondary structure of an exemplar single-stranded DNA sequence was provided in FIG. 2D (upper right). The sequences in the datasets were characterized by a total of 277 face / energy configurations.EXAMPLE 6. Visualizing the Space of Secondary Structure

[0272] This Example describes the use of t-distributed stochastic neighbor embedding (t- SNE) to visualize and analyze structural similarities among thousands of single-stranded DNA sequences based on their secondary structures derived from SELEX-NGS data.

[0273] Data visualization was used to assist in exploring the configuration space by transforming large amounts of secondary structure data into compact visual formats. This transformation facilitated the identification of patterns, relationships, and outliers that might not have been evident if only a few data samples had been analyzed at a time. The data visualization process was used to provide two- or three-dimensional representations of vector- valued descriptors that represented sequences in the dataset. Such lower-dimensional representations were obtained using dimensionality reduction techniques.

[0274] A key aspect of dimensionality reduction was the preservation of the intrinsic geometry of the high-dimensional data manifold, e.g., the preservation of pairwise distances. Recent advancements in ML had produced efficient dimensionality reduction algorithms, such as t-distributed stochastic neighbor embedding (t-SNE) (van der Maaten et al., J. Mach. Learn. Res.9(86):2579-2605 (2008)) (which was used here) and principal component analysis (PCA) (Jolliffe et al., Philos. Trans. R. Soc. A: Math. Phys. Eng. Sci. 374(2065) (2016)). Conceptually, a t-SNE was constructed by minimizing the divergence between two probability distributions: one that measured pairwise similarities of the points in the high-dimensional space and another that measured pairwise similarities of the points in the lower-dimensional space. The goal was to ensure that similar points in high-dimensional space remained close in the lower-dimensional representation.

[0275] A t-SNE was used to perform non-linear dimensionality reduction. The t-SNE projection into the lower-dimensional space was constructed to preserve complex, non-linearrelationships, in contrast to a PCA projection, which was linear. The t-SNE algorithm also depended on a hyperparameter called perplexity, which balanced the focus between local and global data structures during dimensionality reduction. Perplexity influenced how many neighbors each point considered when the data were embedded into lower dimensions. Smaller perplexity values emphasized local relationships by focusing on nearby points, while larger values considered a broader neighborhood, capturing more global patterns.

[0276] Choosing the appropriate perplexity, typically through trial and error, was found to be crucial, as it could significantly affect the visualization outcome. The perplexity value was set to 100. It was noted that the computational cost of t-SNE was quadratic in the number of data points, as similarities for all pairs of data points needed to be computed. Therefore, this method was considered potentially unsuitable for analyzing extremely large datasets. The t-SNE algorithm from the scikit-learn Python library was employed. Pedregosa et al., J. Mach. Learn. Res. 12:2825-2830 (2011). In the experiments, t-SNE computations were completed in less than 32 seconds on the dataset consisting of 4450 sequences represented by the Motzkin-path descriptors.

[0277] FIG. 12 was used to display the t-SNE two-dimensional representations of the Motzkin path descriptors of the sequences in the dataset introduced in Example 4. Each data point corresponded to a specific sequence, with a color corresponding to its count number from the SELEX-NGS output; brighter colors indicated higher counts. FIG. 12 clearly illustrated that data points tended to cluster, indicating that sets of sequences whose secondary structures shared high similarities were present. Additionally, the plot indicated that higher-count structures with brighter colors were generally located in specific clusters or regions of the visualized space. Conversely, regions with solely low-count structures were also observed, such as the lower right corner. It was noted that the high-count structures were not isolated, so that low-count sequences with secondary structures similar to the high-count sequences could be identified.

[0278] For example, the four sequences with the highest SELEX-NGS counts were considered. For each of these sequences, another sequence was identified that possessed a similar secondary structure and a different primary sequence but had the lowest possible count of 1. In FIG. 12, exemplary folded sequences were illustrated, the associated counts for each sequence were provided, and the regions of the two-dimensional space where those structures were represented were indicated. FIG. 13 was used to show the foldings of the top count aptamer from this screen.EXAMPLE 7. Machine Learning Similarity Search

[0279] This example describes a structural and energetic similarity search conducted on a dataset of single- stranded DNA sequences to identify sequences with secondary structures that matched the topological and energetic configurations of high-count aptamer candidates, with the goal of discovering potentially high-affinity sequences beyond those with the highest NGS counts.

[0280] One important task for analyzing large datasets of single- stranded DNA sequences representing aptamer candidates was to identify sequences with folded structures similar to those hypothesized to have favorable binding properties, e.g., those associated with high NGS counts. When the similarities between secondary structures were investigated, two aspects were considered: structural and energetic. Structural similarity was defined based on the topology of folded sequences, while energetic similarity was defined based on the energetic configurations of the secondary structures, such as the energies associated with the faces composing the secondary structures.

[0281] Using Motzkin path (or structural matrix) descriptors, sequences with secondary structures topologically equivalent to that of a given sequence of interest were identified. Using the BoF, sequences with secondary structures having similar energetic configurations were identified.

[0282] Consider D =dataset of single-stranded DNA sequences, Mof associated Motzkin path descriptors, and F ••=which is a set of related BoF representations, where mtand ftare the Motzkin path and BoF representations of the sequences Xi, respectively. Moreover, consider a sequence of interest x, e.g., a high count sequence, with Motzkin path and BoF descriptors m and f, respectively. Given em, £f E IK+, the sets8^£mc B of sequences with secondary structures could be defined topologically em-similar and energetically 6y-similar to x as followsTx,em■= Xi e D such that^x,£f-= {Xi E D such that || i - f||2< efwhere || ■ ||2is the Euclidean norm and emare parameters quantifying structural and energetic similarities between the representations. In particular, T20and 8^0are the set of sequences with secondary structures topologically and energetically equivalent to that of the sequence x, respectively.

[0283] A structural and energetic similarity search was performed on the dataset introduced in Example 4, which consisted of 4450 unique single- stranded DNA sequences. One of the goals of the similarity search was to identify sequences with secondary structures that possessed thesame topological and energetic configurations as those with the highest-count sequences, which were then predicted (but not known) to have good binding properties. In particular, sequences with properties similar to those whose counts were in the top 0.1 percentile were targeted. These were sequences whose count values exceeded those of 99.9%

[0284] Notably, the third-highest-count sequence (with count 48,708) was found to have a secondary structure with the same topology as the top high-count sequence (with count 77,352). Thus, to restrict the analysis to a diverse set of topological structures, only the first, second, fourth, and fifth highest-count structures were studied. The four sequences analyzed were referred to as x with count 77,352, x2with count 67,049, x3with count 28,432, and x4with count 12,126. The secondary structures of x, x2x3, and x4are illustrated in FIG. 12.

[0285] Structural similarity searches were performed on the dataset, and it was found that it contained 103 sequences with secondary structures sharing the same topology as x, 59 sequences similar to x2, 34 sequences similar to x3, and 36 sequences similar to x4. That is, fe,o| = 103, |^2j0| = 59, |^3j0| = 34, and |^4>0| = 36.

[0286] The bar plots in FIG. 14 were used to illustrate the number of different energetic r iprofiles associated with sequences in j.0J4_iand the number of sequences for each of the energetic profiles. It was clearly shown that sequences could have shared the same topology in their folded configuration but were associated with different energetic profiles. In the dataset, the topology of the secondary structure of x was associated with nine distinct energetic profiles, x2with three, x3with seven, and x4with eleven energy profiles. In the bar plots in FIG. 14, the striped bars were associated with the energetic configurations of at least one high-count sequence. In particular, FIG. 14 indicated that the dataset contained 80 sequences with secondary structures sharing the same topology and energetic configuration as x, 56 as x2, 26 as56, | T^3I0A £^>01 = 26,

[0287] Recalled that the bar plots in FIG. 14 related to the first, second, fourth, and fifth highest-count sequences, which were referred to as Xi, Xi, X3, and x4, respectively. The sequence with the third-highest count (48,708) was found to have had a secondary structure with the same topology as Xi, the highest-count sequence.

[0288] In (a) of FIG. 14, only one bar with stripes was shown. That is, only one configuration (Configuration I) was found to correspond to a high-count sequence. This occurred because the third and the first highest-count sequences were determined to share the same energetic configuration.EXAMPLE 8. Machine Learning Topic Modeling and Spectral Clustering

[0289] This Example describes how latent topic modeling and spectral clustering were applied to BoF descriptors of single-stranded DNA sequences to identify clusters of sequences with structural and energetic configurations similar to high-count aptamers, thereby facilitating the targeted selection of promising candidates from large datasets.

[0290] Sequences were also clustered based on the energetic configurations of their secondary structures. To achieve this, latent topic models were built using BoF descriptors. Next, the identified topic distributions were used to describe and cluster the sequences.

[0291] Topic modeling was used as a machine learning technique to discover hidden semantic structures within a corpus of documents. By analyzing word co-occurrences across documents, underlying topics were identified, where each topic was represented as a distribution of words, and each document was considered a mixture of these topics. Common unsupervised ME methods, such as Eatent Dirichlet Allocation (LDA) and Nonnegative Matrix Factorization (NMF), enabled the automatic association of each text with a distribution over topics. Given a text and a number of topics, these methods quantified how relatable the given text was to each identified topic.

[0292] To perform topic modeling, Nonnegative Matrix Factorization (NMF) was used, which was a linear algebraic method. Einear approaches had the advantage of being computationally efficient and scalable to massive datasets. The NMF model was used to discover interpretable latent components in high-dimensional unlabeled data. It analyzed a set of documents described by the counts of unique words. In particular, NMF took as input the so-called term-document matrix X E Blvxt. Each of the w E H rows of X corresponded to a unique word in the vocabulary, and each of the t E H columns corresponded to a text. The (i, j)-th entry of X represented the number of occurrences of the ith word in the j th document.

[0293] Next, NMF was used to decompose X E Blvxtinto two lower-dimensional nonnegative matrices, W E Blv,kand H E IKk,t, such that: X ~ W • H. W was used to capture the basis vectors (e.g., topics), and H was used to encode the coefficients or contributions of each basis vector (e.g., the importance of each topic in each document). Both W and H were constrained to be non-negative, which aligned with many real- world scenarios, such as text data where counts could not be negative. The hyperparameter k was defined a priori. It determined the number of latent components (or topics) identified by the NMF algorithm. The NMF model computed W and H by solving the following optimization problem:W. H > 0 "X-W H"i ’where || ■ ||Fis the Frobenius norm.

[0294] Here, single-stranded DNA sequences were considered instead of text. Occurrences of face-energy configurations across secondary structures were analyzed instead of words.Consequently, in the matrix X that was given as input to NMF, each column was associated with a sequence, and each row corresponded to a face-energy configuration. The (i, j)-th entry of X counted the number of occurrences of the i-th face-energy configuration in the j -th document. Here, each column of the matrix X was the BoF descriptor of the related sequence.

[0295] The rows of W were used to characterize each topic as a distribution of face-energy configurations, while the columns of H were used to describe the secondary structure of each sequence by quantifying how relatable it was to each one of the identified topics. In particular, the columns of the matrix H were used to associate the secondary structure of each sequence with a mixture of different topics. Hence, secondary structures were described with their associated topic mixture distributions. An additional concept of similarity was introduced. The distance between any two secondary structures was measured by comparing how dissimilar their topic mixture distributions were. This meant determining how different the columns of matrix H were, which were associated with the two sequences of interest. The NMF algorithm from the scikit-leam Python library was employed, considering 25 topics. The number of topics for topic modeling was selected by the user. In the experiments, NMF on the 4450 BoF descriptors was completed in less than two seconds.

[0296] Next, clusters of sequences described by similar topic mixture distributions were identified. Spectral clustering, a well-known graph-based clustering machine learning technique, was employed. This method leveraged the eigenvalues of a similarity matrix to partition data into distinct groups. A similarity graph was constructed from the input data points, and the eigenvectors of the graph Laplacian were used to identify clusters. Spectral clustering effectively captured non-linear structures in the data, making it particularly useful for complex clustering tasks in machine learning. The spectral clustering algorithm from the scikit-learn Python library was used. The number of clusters to be identified was defined a priori and provided as input to the clustering algorithm. In line with the choice to consider 25 topics, 25 distinct clusters were targeted. In the experiments, spectral clustering was completed in less than three seconds on the dataset consisting of 4450 sequences represented by the topic mixture distributions.

[0297] Recall that the primary goal of the data analysis was to identify sequences that had similar structural and energetic configurations as high-count sequences. A sequence was defined as high-count if the associated count was in the top 0.1 percentile. There were five suchsequences in the dataset, four of which were x15x2, x3, and x4, which had been analyzed in the previous section. The fifth sequence was the third highest count sequence, which was found to be topologically and energetically equivalent to x. The bar plot in FIG. 15 was used to illustrate the number of sequences per cluster. The bars with stripes were those associated with a cluster containing at least one high-count sequence. Only four of the twenty-five clusters were found to contain one of the five high-count sequences, as the first and third highest-count sequences shared the same structural and energetic configuration and were associated with the same cluster. Cluster M was found to contain the highest count sequence x and the third highest count sequence, cluster E contained x2, cluster W contained x3, and cluster X contained x4. According to the clustering results, the clusters with at least one high-count sequence contained 457 elements. That is, fewer than 11% of the sequences in the dataset were found to share similarities with a high-count sequence, suggesting that further research efforts should be prioritized on the analysis of this subset of sequences, while screening out the remaining 89%.

[0298] FIG. 15 displays a t-SNE two-dimensional representation of the descriptors that were obtained via topic modeling. Each data point was associated with a specific color, with each color corresponding to a specific cluster. Points with the same color were grouped in the same cluster. FIG. 15 enabled visualization of the entire space of folded structures and the cluster associated with each sequence. In the figure, the clusters that included at least one high-count sequence were also highlighted with the corresponding alphabetic character. For each of these clusters, an example of a low-count secondary structure that shared structural and energetic similarities with the high-count sequence in the cluster was proposed. The similarities were determined based on the descriptors that were obtained through topic modeling.

[0299] Interestingly, even though the clusters had been computed on the high-dimensional descriptors, the two-dimensional representations within each cluster were observed to aggregate in specific regions of the Euclidean space. This suggested that the similarity of the highdimensional descriptors was well represented in the lower-dimensional space. Moreover, based on the geometric distribution of the two-dimensional representations, it was observed that points within the same cluster could aggregate in more than one region of the Euclidean space. This suggested that the dimensionality reduction highlighted additional hidden patterns between sequences within the same cluster that could potentially have been exploited to further segment the larger clusters.EXAMPLE 9. Corrections to the SeqFold Energy Functions

[0300] This example describes improvements to energy computations in the SeqFold algorithm by correcting the assessment of terminal mismatches in single base pair loops and modifying the energy model for junctions to better align with mfold behavior.

[0301] The presently improved SeqFold and the presently described machine learning algorithm were primarily based on the functions used in SeqFold to compute energy values of the various faces, such as hairpins, bulges, and inner loops. However, two straightforward improvements were made to these energy functions. Specifically, the computation of the energies associated with internal loops generated by single base-pair mismatches and with junctions, which were multi-branches with no unpaired nucleotide between the branches, was modified.

[0302] The energy associated with internal loops from single base pair mismatches was determined based on the energy of the terminal mismatches of the loop. Terminal mismatches, which occurred at the ends of a loop’s double-stranded regions, consisted of four nucleotides: two that were base pairing and two that were inside the loop and not matching. A left and right terminal mismatch was associated with each loop depending on the 5’ to 3’ direction. FIG. 16 illustrates a folded sequence with an internal loop resulting from a base pair mismatch and highlighted the left and right terminal mismatches. The right terminal mismatch, enclosed in the dashed box, consisted of nucleotides CC / AG, while the left one, enclosed in the continuous line box, consisted of nucleotides GC / CA. Following along Santalucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004), the free energy associated with a base pair mismatch loop was computed using the following formula:GT° (single mismatch loop)= GT° (left terminal mismatch) + GT° (right terminal mismatch), where GT(left (right) terminal mismatch) was calculated using Formula (2) from the enthalpy and entropy associated with the terminal mismatches. Table 2 and Table 3 report the energies SeqFold assigns to terminal mismatches of inner loops. Table 2 was used when the terminal mismatch did not involve any nucleotide that was either the first or the last in the sequence. Such mismatches are referred to as ‘internal.’ Table 3 was used when the terminal mismatch includes either the first or the last nucleotide in the sequence.Table 2. Enthalpies and Entropies associated with DN internal terminal mismatches*.* The reported values are taken from SeqFold GitHub repository and are valid for temperature T=37°. Internal terminal mismatches are sets of four nucleotides WX / YZ where W is not the first nucleotide of the DNA string and Z is not the last. The energies are the same for each terminal mismatch in the reverse direction. That is, WX / YZ has the same energy as ZY / XW.Table 3. Enthalpies and Entropies associated with DNA terminal mismatches.* The reported values are taken from SeqFold GitHub repository and are valid for temperature T=37°. The energies are the same for each terminal mismatch in the reverse direction. That is, WX / YZ has the same energy as ZY / XW.

[0303] The energy associated with single base pair mismatch loops was computed inaccurately by SeqFold due to an incorrect assessment of the left terminal mismatch.Specifically, according to SeqFold, the left terminal mismatch was considered to consist of the base pairs at the end of the two stacks that determined the loop. For instance, in the example shown in FIG. 16, SeqFold identified GC / CG as the left terminal mismatch. This issue was corrected, and terminal mismatches were assessed accurately.

[0304] Regarding the computation of energy associated with junctions, junctions were associated by SeqFold with energy < Em= 4.6 + Est, where Est was the sum of the energies associated with each of the branches (or stacks) of the junction. Alternatively, junctions were associated with lower energy by setting < Em= 0.6 + Est. By associating junctions with lower energy, these structures were considered to be more stable than as assessed by SeqFold. This modification of the energy computation was not supported by empirical laboratory experiments. Instead, the constant in the junction energy computation was tuned using a heuristic approach to emulate mfold folding behavior. It was noted that, based on experience, both SeqFold and the presently improved SeqFold consistently failed to predict junctions where both the presently described algorithm and mfold succeeded. Examples were provided in FIG. 8, which illustrated a comparison of the secondary structures of two exemplary sequences computed using SeqFold, the presently improved SeqFold, the presently described machine learning algorithm, and mfold.FIG. 8 shows that SeqFold and the presently improved SeqFold failed to predict the junction in both examples. After an initial empirical investigation, it was speculated that this limitation of the SeqFold approaches was not only due to the junction energy computation but also due to the implementation of the dynamic programming approach. In particular, it was believed that SeqFold and the presently improved SeqFold did not consider junctions as possible faces.EXAMPLE 10. Detailed Energy Functions

[0305] This Example describes how energy values for different RNA secondary structure motifs, such as hairpins, bulges, internal loops, and multibranch loops, were computed usingthermodynamic parameters and formulas from SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004), and implemented in SeqFold, the presently improved SeqFold, and the presently described machine learning algorithm.

[0306] SeqFold, the presently improved SeqFold, and the presently described machine learning algorithm were used with the same energy functions for hairpin loops, internal loops, bulges, and stacks. As described in SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct.33(l):415-440 (2004), the energies for 37 °C with a salt concentration of 1 M NaCl were employed. The code was configured to allow modification of the temperature. The free energy values were calculated from enthalpy and entropy values provided in a table from SantaLucia et al., and were adjusted for the appropriate temperature using the following formula (noting that the enthalpy increment was given in kilocalories / mol and the entropy increment was given in calories / mol, hence the factor of 1000):AS AGT= AH — T ■ — —T1000Formula (2)

[0307] For ease of writing, everything was expressed in terms of AG° from that point onward, although it was noted that all of the tables stored values of entropy and enthalpy.

[0308] As previously discussed, the energy of the structure associated with the substring [i, j] was computed as the sum of the face energy having the pair (i, j), and for each internal pair in this face, (i',j'), the energy of the structure associated with the substring [i', j'] was added. This was implemented by finding the minimum of three energies, corresponding to the three types of faces that could be present.

[0309] For hairpins, Table 4 was used for all loops of length 3 and 4, which depended on the base pairs present (SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004)). For hairpin lengths 5 through 30, Table 5 depending only on the length of the loop was used. For anything longer, the Jacobson-Stockmayer energy extrapolation formula was applied:Formula (3)

[0310] Here n is the length of the sequence one wants to calculate, nmax is the maximum length with an experimentally known value, R is the ideal gas constant, and 2.44 is an experimentally determined constant. The terminal mismatch penalty was added to this energy value.Table 4. Enthalpies and Entropies associated with all hairpin loops of length 3 or 4*.* AH is in units of kcal / mol, and AS is in units of cal / mol. The reported values are taken from SeqFold GitHub repository and are valid for temperature T=37°. literature values given in SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004).Table 5. Loop* These values are used to compute the energy of inner loop, bulge loops, and hairpin loops. The energy of larger loops is calculated using the Jacobson-Stockmayer energy extrapolation formula (SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004)).

[0311] The formula for an inner loop was a bit more complicated. From SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004), the following formula was used:Formula (4) where the components are AG°, the energy given the length of the inner loop (Table 5).

[0312] AG°sym= \nt— nr\ x 0.3 kcal / mol is the asymmetric penalty (nt, nrare the lengths of each side of the loop). AGT°N Nis the terminal mismatch penalty, an d is the sum of the penalty from the inner end of the loop and the outer end of the loop.

[0313] If the interior loop was instead a stack, the stacking energy was looked up in the Tables of the present disclosure. If the bulge had size 1, then the stacking energies on each side were added, followed by the penalty for the internal nucleotide and the penalty due to the bulge strain; if one of the closing ends was an A, an AT penalty was applied:

[0314] For longer bulges, the penalty for the loop was no longer made dependent on which nucleotides were contained in the bulge, but only on the length; the values for this were provided by Table 5 for lengths 1 through 30, and the Jacobson-Stockmayer formula was used to extend to loops longer than 30.

[0315] For multi-branches, the face energy is given bywhere a, b, c, d are energy parameters defined in SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004) (values (2.6, 0.2, 0.2, 2.0) were used herein), Nbr is the number of branches in the multibranch structure, Nunpaired is the number of unpaired nucleotides in the multibranch, and 6fuUy stacked iszeroif Nunpaired0, and one otherwise, and corresponds to the stabilization that is gained by the multibranch being fully stacked, sometimes called the coaxial stacking stabilization. Note that this particular energy function is linear with the length of the multibranch. SantaLucia et al., Annu. Rev. Biophys. Biomol. Struct. 33(l):415-440 (2004) suggest that a logarithmic dependence may be more optimal, i.e. an energy function similar to the Jacobson-Stockmayer equation. If given appropriate energy values, such a function could beimplemented in either the presently improved SeqFold or the presently described machine learning algorithm.EXAMPLE 11. Illustration of the Differences in the Algorithms with Aptamers

[0316] This Example describes differences in folding outcomes between the presently described machine learning algorithm and mfold, highlighting variations in coaxial stacking and multibranch stability due to distinct optimization methods and energy computations.

[0317] As described in Section II, C supra, the presently described machine learning algorithm was noted to have a few key differences from mfold not only in the optimization technique but also in energy computation. Several examples were provided to illustrate these differences. FIG. 17 was used to exemplify the difference in coaxial stacking computation, as the mfold structure appeared to gain additional stability in the stacked configuration when compared to the configuration found by the other software. FIGS. 18 and 19 displayed aptamers for Cocaine and a membrane protein of a glioma cell line SHG44, respectively. These aptamers served as examples of an open multibranch structure (z.e., the “open face” had multiple branches), where it appeared that methods other than mfold did not identify the same stability for this type of structure. The examples provided herein also indicated that there may have been subtle changes in the energy values of particular structures used in mfold’ s code that were not found in the literature.

Claims

CLAIMSWe claim:

1. A computer-implemented method for high throughput classification and functionality assessment of aptamer structures, the method comprising:(a) receiving, by a processor, a plurality of aptamer sequences;(b) generating, by the processor, a structure for each aptamer sequence of the plurality of aptamer sequences to generate a plurality of structures;(c) identifying, for each structure of the plurality of structures, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures; and(d) classifying, by applying the set of parameters of the each structure into a trained machine learning model, the plurality of structures into one or more aptamer clusters.

2. The method of claim 1, further comprising (e) determining a target binding characteristic for an aptamer of the one or more aptamer clusters.

3. The method of claim 1 or 2, further comprising generating a visualization of the one or more aptamer clusters.

4. The method of any one of claims 1-3, wherein the structure, for the each aptamer sequence, is a secondary structure or a tertiary structure based on the each aptamer sequence.

5. The method of any one of claims 1-4, wherein the generating the structure for each aptamer sequence comprises:optimizing a selection of a set of bonds between pairs of nodes of the each aptamer sequence to achieve a minimum free energy.

6. The method of any one of claims 1-5, wherein the generating the structure for each aptamer sequences comprises:(i) identifying a plurality of nodes in the each aptamer sequence;(ii) determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the each aptamer sequence; and(iii)selecting, among the plurality of admissible bonds, a set of bonds between the each pair of nodes to form the structure by optimizing to achieve the minimum free energy for the structure.

7. The method of claim 6, wherein the selecting the set of bonds comprises:(i) determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to a face class comprising one or more template graphs of the structure; and(ii) filtering, from the selecting the set of bonds, the admissible bond that fails to conform to the one or more template graphs via a two-level stacking region template.

8. The method of claim 7, wherein the face class comprises a hairpin loop, a stacking region, a bulge loop, an interior loop, or a multi-branch loop.

9. A system for high throughput classification and functionality assessment of aptamer structures, the system comprising:(a) a processor; and(b) memory storing instructions that, when executed by the processor, cause the processor to:(i) receive a plurality of aptamer sequences;(ii) generate a structure for each aptamer sequence of the plurality of aptamer sequences to generate a plurality of structures;(iii)identify, for each structure, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures; (iv) classify, by applying the set of parameters of the each structure into a trained machine learning model, the plurality of structures into one or more aptamer clusters; and(v) determine a target binding characteristic for an aptamer of the one or more aptamer clusters.

10. The system of claim 9, wherein the memory storing instructions, when executed, further cause the processor to: generate a visualization of the one or more aptamer clusters.

11. The system of claim 9 or 10, wherein the structure, for the each aptamer sequence, is a secondary structure or a tertiary structure based on the each aptamer sequence.

12. The system of any one of claims 9-11, wherein the instructions, when executed, cause the processor to generate the structure for the each aptamer sequence by:optimizing a selection of a set of bonds between pairs of nodes of the each aptamer sequence to achieve a minimum free energy.

13. The system of any one of claims 9-12, wherein the instructions, when executed, cause the processor to generate the structure for the each aptamer sequence by:(i) identifying a plurality of nodes in the each aptamer sequence;(ii) determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the each aptamer sequence; and(iii)selecting, among the plurality of admissible bonds, a set of bonds between nodes to form the structure by optimizing to achieve the minimum free energy for the structure.

14. The system of claim 13, wherein the instructions, when executed, cause the processor to select the set of bonds by:(i) determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to a face class comprising one or more template graphs of the structure; and(ii) filtering, from the selection of the set of bonds, the admissible bond that fails to conform to the one or more template graphs via a two-level stacking region template.

15. The system of claim 14, wherein the wherein the face class comprises a hairpin loop, a stacking region, a bulge loop, an interior loop, or a multi-branch loop.

16. A non-transitory computer-readable medium (CRM) having stored thereon computer- readable instructions executable to cause performance of operations comprising:(a) receiving, by a processor, a plurality of aptamer sequences;(b) generating, by the processor, a structure for each aptamer sequence of the plurality of aptamer sequences to determine a plurality of structures;(c) identifying, for each structure, a respective set of parameters assessing a similarity criteria between the structure and another structure of the plurality of structures;(d) classifying, by applying the set of parameters of the each structure into a trained machine learning model, the plurality of structures into one or more aptamer clusters; and(e) determining a target binding characteristic for each aptamer cluster of the one or more aptamer clusters.

17. The non-transitory CRM of claim 16, wherein the computer-readable instructions are executable to cause performance of operations further comprising: generating a visualization of the one or more aptamer clusters.

18. The non-transitory CRM of claim 16 or 17, wherein the structure, for the each aptamer sequence, is a secondary structure or a tertiary structure based on the each aptamer sequence.

19. The non-transitory CRM of any one of claims 16-18, wherein the generating the structure for each aptamer sequence comprises:optimizing a selection of a set of bonds between pairs of nodes of the each aptamer sequence to achieve a minimum free energy for the structure.

20. The non-transitory CRM of any one of claims 16-19, wherein the generating the structure for each aptamer sequences comprises:(i) identifying a plurality of nodes in the aptamer sequence;(ii) determining a plurality of admissible bonds between each pair of nodes of the plurality of nodes in the aptamer sequence; and(iii) selecting, among the plurality of admissible bonds, a set of bonds between nodes to form the structure by optimizing to achieve the minimum free energy for the structure.

21. The non-transitory CRM of claim 20, wherein the selecting the set of bonds comprises:(i) determining, for an admissible bond of the plurality of admissible bonds, that the admissible bond fails to conform to a face class comprising one or more template graphs of the structure; and(ii) filtering, from the selection of the set of bonds, the admissible bond that fails to conform to the one or more template graphs via a two-level stacking region template.

Images