Aptamers are defined as synthetic ligands composed of single stranded RNA or DNA which have unique structures that can be designed to target specific biomolecules with high affinity [9]. The production of aptamers is achieved through a process of systematic evolution of ligands by exponential enrichment (SELEX) and was discovered by three groups independently at the beginning of the 1990s, beginning with RNA aptamers then expanding into DNA and peptide-based aptamers (referred to as affimers) developed today [10, 11, 12]. Nucleic acid-based aptamers are single stranded nucleic acid structures which can form secondary structures by folding in on themselves into thermodynamically stable shapes [13]. These structures consist of double helices, and hair-pin structures of varying sizes governed by the number of mis-match pairing of bases which produce internal loops and three-way helical regions [14]. While SELEX has been used previously for the development of aptamers, there have been movements towards the use of in silico computer-based modelling to streamline the production and design of aptamers [15]. Techniques like M.A.W.S (making aptamers without SELEX) are the most notable recent examples of in silico aptamer development and has been seen as a gold standard for the development of aptamers in silico [16]. The focus of aptamer generation in silico looks to remove some of the initial high costs associated with aptamer selection, as well as looking to produce high affinity aptamers prior to the purchase of reagents needed for in vitro validation of aptamer binding [17].

One of the main advantages of aptamers over their monoclonal antibody counterparts is tissue penetration, as well as a capability to be raised against targets that do not elicit an immune response, such as charged ions, intracellular proteins and vitamins [18]. The use of aptamers has emerged as a potential detection strategy for a wide number of biomarkers, with particular interest in the detection of markers for neurological disease [19]. In particular, alpha synuclein ($\alpha$-syn) has emerged as a candidate in the diagnosis of neurological diseases such as Parkinson’s (PD) and Alzheimer’s (AD).

$\alpha$-syn is a 140 amino acid protein encoded by the SNCA gene located on chromosome 4 and has been implicated in the progression of Lewy body, Alzheimer’s and Parkinson’s disease associated dementia [20]. Under physiological conditions, $\alpha$-syn is maintained in an unfolded state but in disease states it is capable of forming toxic oligomers and aggregates with local proteins and has been linked to toxicity via mitochondrial damage, inflammation and axonal transport dysfunction (Fig. 1, Ref. [21]) [22].

Fig. 1.

Summary of the interactions in AD. Notably $\alpha$-syn is associated directly with synapse loss and neurodegeneration in AD, however this is separate from the formation of amyloid plaques which are typically associated with the condition as described in [21].

The link between AD and $\alpha$-syn has been established, however the link with PD has not been fully confirmed but has the potential to be used to establish differences in diagnosis with lower cerebrospinal fluid (CSF) levels being detected in PD compared to AD [23]. Efforts to prevent the aggregation of $\alpha$-syn, and by extension alleviate symptoms of AD, have focused on the use of aptamers due to their inherent ability to penetrate the blood brain barrier. Several aptamers have been identified which have successfully targeted $\alpha$-syn and have demonstrated mitochondrial protective effects as well promoting $\alpha$-syn degradation in vitro [22, 24, 25]. The approaches taken when producing anti $\alpha$-syn aptamers has been the traditional SELEX process, with minor modifications that have been incorporated over time to improve the methodology. With the advances in computer simulations and molecular docking experimentation, interest in improving the experimentally derived aptamers has increased exponentially.

The concept of computer aided drug design, and by extension in silico aptamer design, has been gaining traction in recent years, specifically due to its ability to streamline the design of innovative treatment and detection strategies for a wide variety of applications. Advances in computer algorithms capable of predicting the folding characteristics of single stranded DNA, as well as increasing availability of these software, have been critical to the development of in silico methods [26]. The standardised approach to in silico aptamer design is summarised in Fig. 2 (Ref. [27]).

Fig. 2.

Summary of steps within in silico aptamer design. RNA folding algorithms have proven to be effective at predicting DNA aptamer structure despite minor chemical differences. It is common practice to correct these structures in post. Image based on information contained within [27].

Sequence generation can be upscaled to be comparable with in vitro SELEX, allowing for the screening of an exceptionally large number of sequences. The difficulties come when establishing selection criteria for the computer-generated sequences, and which ones to proceed with in a docking experiment as the large number of sequences make it impractical to dock every aptamer produced. A common approach is to look at the free energy potential of the aptamers secondary structure [26] and secondary structural motifs [28] to screen potential aptamer candidates. The second option is particularly useful when looking at aptamer optimization studies, as the original sequence can be mutated at specific points to see what impact this will have on overall binding.

Following sequence generation, aptamers are then docked to the target molecule in question. There are a number of programmes which are available for use in this field, some which are freely available and some requiring a subscription [28]. Irrespective of the programme used, the objective is to create a pose (two molecules which have formed a complex) with the lowest Gibbs free energy ($\mathrm{\Delta }$G) binding energy [27]. Experiments themselves can take two main approaches, fixed or flexible docking which can be further sub-divided into 5 categories. The approaches are summarised in Table 1 [29].

In silico docking is not without its limitations however, despite the advances in computer modelling and DNA/RNA folding algorithms. Issues with computer processing load [27], particularly with more flexible approaches to docking, as well as the difficulties in accurately mimicking the electrochemical environment these molecules would experience in vivo [30] do mean that at this stage, predictions need to be validated in vitro before definitive conclusions can be made. However, this does not detract from the usefulness of in silico methods for they have the potential to both streamline and improve the aptamer generation process.

In this study, an experimentally produced aptamer was docked alongside an in silico developed aptamer against human $\alpha$-syn. The intention of this study was to demonstrate a novel method by which aptamers can be developed in silico and validated against previously observed interactions, with the hope of refining that aptamer to improve the binding efficiency to its target. To our knowledge, this is the first investigation looking into the refinement of an anti-human $\alpha$-syn DNA aptamer to produce a structure with a predicted higher binding affinity to the experimentally derived sequence.

Aptamer sequences were generated using Python script to create a library of 1000 sequences with an N30 wobble region flanked by constant primer binding regions (ATAGTCCCATCATTCATT-N30-AGATATTAGCAAGTGTCA). Sequences were retained following prediction of thermodynamic stability using the DINAMelt server-QuickFold web server [31]. Retained sequences were required to have negative $\mathrm{\Delta }$G values, using the parameters listed in Table 2.

Secondary structural information for 50 structures, selected based on their $\mathrm{\Delta }$G, was obtained using MFold [32]. Structure prediction in MFold used the same parameters as were used within the DINAMelt server-QuickFold.

Vienna sequence outputs from MFold were used to assist in the production of 3 dimensional models using RNAcomposer [33]. To account for this being a programme which folds RNA, all thymine residues within the sequences were replaced with uracil groups for the purpose of folding. Vienna sequences were input along with the base sequence data, with the top 50 sequences being run in batch format within RNAcomposer. Structure prediction used default parameters for the purposes of folding.

Protein data bank (PDB) output files of the predicted sequences were converted into mol2 format using CSD-Discovery-Mercury [34]. Structures were then modified to reintroduce the deoxyribose groups, replacing the ribose sugar groups within their structures using PyMol [35] via command line editing within the software (Schrödinger, New York, USA). Commands were set to remove O2’ and HO2’ atoms from all residues and to add a hydrogen to the C2’ atom of each residue. Uracil groups were replaced with thymine groups using the mutagenesis wizard also contained within the PyMol software. Files were then exported back into CSD-Discovery-Mercury to add any hydrogens to the structure. In addition, bond angles and types were corrected as some had been incorrectly modelled during the file conversion process. Of the 50 structures that were originally generated, the top 10 were selected to dock with the target molecule. Primary structural similarity to M5-15 was used as the primary criteria when selecting the aptamer, followed by thermodynamic stability. This is due in part to the nature and composition of free hydrogen linkage points within the short structure determining how the single stranded nucleotide will fold in on itself [36].

Docking experiments were analysed determining the predicted interactions between the components of the final pose. Bond distance, nature of interactions and type of interaction were all documented and analysed. Overall fitness of poses was determined using ChemPLP fitness scoring within CSD-Discovery. Structures were virtually screened in blocks of 10 using LYS 96, O as a central atom and a radius of 81Å, encompassing the target protein completely. All rotatable bonds were fixed to create a rigid dock and search efficiency was set to low/virtual screening. Fitness scoring function was set to ChemScore with no re-scoring function selected.

Binding energies were determined using follow up confirmation docking experiments using AutoDock Vina. Docking parameters between the programmes were kept constant. The docking positions generated by CSD-Discovery were used to inform the choice of grid location for docking in AutoDock Vina. The top 3 aptamers selected from CSD-Discovery were docked against $\alpha$-syn using the following settings: grid size: 110, 110, 110; coordinates: center_x = 23.673, center_y = 56.090, center z = –15.985; energy range = 4; exhaustiveness = 8; ligand torsions set to zero during setup. A summary of data files used within the work can be found in Table 3 [37, 38]. A spinach protein aptamer, 4ts2, was chosen as a control as it has no specificity to the target protein. Comparisons with a sequence which binds to a protein which is in an unrelated organism gives an indication of the baseline interaction between the target and a nucleotide sequence.

Top poses from CSD GOLD were used to determine the best aptamers to move forward to docking in Autodock Vina based on a combination of Fitness score, number of hydrogens and average bond length of interactions. Analysis of variance was conducted on the data produced by Autodock Vina to determine any significant changes in binding energies in comparison to M5-15.

The in silico selection of ssDNA aptamers against $\alpha$-syn using simulated SELEX produced numerous potential candidates which were shown to be comparable to the reference aptamer made by [37] within the limitations of the model. Virtual screening of the top 50 aptamers based on $\mathrm{\Delta }$G produced a pool of 10 aptamers which are shown to bind to $\alpha$-syn (see Table 4). The pool of candidates displays a large variation in $\mathrm{\Delta }$G due to differences in GC content and the presences of Watson-Crick base pairing.

The known sequence which binds to $\alpha$-syn was shown to form three hairpin/stem loop structures when folded at 25 ${}^{\circ }$C under the ionic conditions [37]. Many of the in silico sequences selected based on $\mathrm{\Delta }$G have shown the formation of similar binding motifs which may be effective at targeting $\alpha$-syn due to structural homology (see Fig. 2). Lower $\mathrm{\Delta }$G values were also shown to be heavily correlated to the formation of Watson-Crick base pairing as opposed to higher $\mathrm{\Delta }$G values which are the result of linear regions within secondary structures.

Pairwise alignment (Fig. 3) was conducted to determine the sequence homology between each candidate and M5-15. When aligned, TMG-79 was shown to have the largest percentage similarity (51.4%), suggesting a degree of similarity between the structural motifs of TMG-79 and the M5-15 (see Figs. 2,4). Aptamers 123 and 125 both show lower similarities than TMG-79 which would explain structural differences within aptamer 123 and 125 such as the appearance of structural motifs not present in M5-15 and the change in relative position of stem loops.

Fig. 3.

Pairwise alignment of N30 wobble region for top three selected aptamers based on PLP Fitness score against experimentally created aptamer. Similarity: TMG-79 = 51.4%; Aptamer 123 = 20.8%; Aptamer 125 = 9.3%.

Fig. 4.

Predicted secondary structures for top 10 aptamers based on $\mathrm{\Delta }$G predicted by mfold for sequences in digitally created library. M5-15, referred to here as the reference aptamer [37] raised against $\alpha$-syn, is included as a control within the process. (Energy rules: DNA, temp = 24 ${}^{\circ }$C, 125 mM [Mg], 0 mM [Na]).

GOLD docking results show a large improvement in fitness scores when comparing positive and negative controls to in silico aptamers. TMG-79 PLP fitness (Fig. 5) indicates it could be a candidate for greater binding potential. Further analysis of docking data shows that aptamer complexes for candidates 123 and 125 have a larger number of predicted H bonds with a comparable average bond length (Table 5).

Fig. 5.

Data showing the ChemPLP fitness score for the top poses of the previously selected aptamer sequences. The reference aptamer and spinach aptamer serve as positive and negative controls, respectively. Each selected aptamer shows a higher docking score for their top pose when compared to the reference aptamer.

Despite having a lower fitness score than that of TMG-79, aptamer 123 was predicted to form a larger number of H-bonds making it a promising candidate. This aptamer also appears to show some minor structural similarities to the reference aptamer, particularly when examining its hairpin loop structures. Large variations in ChemPLP score were observed across the poses for each aptamer due to the low search efficiency used for virtual screening.

In addition, follow up experiments using AutoDock Vina demonstrating top pose and average binding affinity for the reference, spinach and the top 3 in silico generated aptamer sequences. Binding free energy was lower for all generated sequences compared to the M5-15, these changes were not statistically significant (TMG-79, p = 0.18; aptamer 123, p = 0.17; aptamer 125, p = 0.42; Spinach, p = 0.18).

Further analysis of the reference aptamer and TMG-79 showed that interactions occurred mainly in or around the stem loop structural motifs of each aptamer. The reference aptamer was predicted to form H-bonds at residues LYS97 and LYS21. In contrast, TMG-79 was predicted to form H-bonds at protein residues LYS80 and GLN99 suggesting that it may interact in a similar region to the reference aptamer. Differences in the overall structure of TMG-79 compared to the reference aptamer support the observation that the interactions are similar, but at different points within the protein. Aptamer 123, whilst having a lower overall ChemPLP fitness score, had a greater number of interaction points in these simulations, importantly with polar side chains of amino acids such as THR81, suggesting stability of these bonds within the complex. However, given that ChemPLP score is a summary of critical factors, including binding energy, this suggests that while there are greater interaction points, these interactions are less likely to form spontaneously and be stable. Bond distances observed from the reference aptamer are also greater than in 2 of the 3 aptamers generated, and mean binding energy predicted greater affinity across all aptamers screened. This indicates a tighter binding potential with the generated sequences compared to M5-15. Structural comparisons of the poses suggests that all aptamers bind to a similar region within $\alpha$-syn, primarily to residues contained within the 2 alpha helix motif containing LYS21 LYS43 LYS45 THR59 LYS80 THR81 SER87 LYS96 LYS97 GLN99 ASN103. These residues appear to be critical to the stabilisation of aptamer binding within the protein; the poses containing the in silico aptamers (Fig. 6) share remarkable similarity to the reference aptamer. Prediction of the final location of the aptamer within the central cleft within the target molecule also suggests that this may be a preferential location for a DNA molecule within $\alpha$-syn.

Fig. 6.

Showing the binding position and orientation of each aptamer in their top pose when docked against alpha synuclein. Each aptamer is predicted to bind within the alpha helices of $\alpha$-syn similarly to the reference aptamer.

When comparing the reference aptamer interaction with the generated sequences, there are similarities between the way they interact with $\alpha$-syn, specifically around the preference to bind within the alpha helix flanking structure within the protein. However, there are differences in the predicted binding between the reference and generated aptamers. Molecular docking relies on predicting the presence of intermolecular forces between two molecules, including van der Waals forces, hydrogen bonds, electrostatic interactions, side chain interactions, hydrophobic interactions and interaction of the solvent on docking (usually water) [39]. Scoring of a docking pose is a combination of all these factors, based on the assumptions and conditions that were imposed during the docking experiment [40]. The likelihood of a docking pose being true can be argued based on known interactions between amino acid residues, the involvement of any side chains and how they interact with DNA molecules. Within the docking poses proposed here, there are a number of interaction points that are observed particularly within amine groups and aromatic side chains between the DNA and the amino acids. These observations have been documented previously. Within all of the structures proposed, there were 5 amino acid residues which were consistently involved in some capacity. These amino acids were, threonine (Thr), lysine (Lys), glutamine (Gln), asparagine (Asn) and serine (Ser).

Thr and Lys interactions with DNA are well documented, with Thr interacting predominantly with phosphate backbone groups [41], and Lys predominantly interacting with guanine groups [42]. Specifically, when looking at Lys, it has been reported to interact with oxidized guanine residues (8-Oxo-G). Binding induces conformational change in the 8-oxoguanine DNA glycosylase (Ogg) required for the removal of 8-Oxo-G, an important DNA repair mechanism [42]. While there are no 8-Oxo-G residues present within the structure of the aptamer proposed here, it appears that there is a natural interaction between Lys and guanine. Thr based hydrogen bonding has the highest propensity in highly hydrophilic environments [43], and since the target protein is water-soluble and primarily found in neuronal tissue, there is a suggestion that this interaction may occur in vivo. Following on from this, we should consider the possibility of the predicted Gln, Asn and Ser interactions between the generated aptamer and the target. Gln has been reported previously to have strong association with polar side chains and phosphates contained within the backbone, as well as interacting with the nitrogenous base component of adenosine [41]. Asn is generally associated with base recognition binding, predominantly to adenosine, and has been demonstrated in synthetic zinc finger constructs like 1meyC [44]. These interactions have also been observed in B-DNA binding of transcription factor-like elements [45]. Ser rich proteins are essential pre mRNA splicing factors. They regulate constitutive splicing and influence alternative splicing processes through an extensive and controlled phosphorylation of the Ser residues [46]. Ser residues are highly conserved within the DNA binding domains of bZip transcription factors of Arabidopsis thaliana, and were reported to be in physical contact with DNA [47] specifically providing a binding point through the hydroxymethyl side chain, something that was observed here in our model.

Interestingly, the non-species aptamer 4ts2 demonstrated a greater binding affinity for the target M5-15 which was an unexpected outcome of the work. This may be in part due to the very nature of $\alpha$-syn and its ability to associate with other “sticky” proteins. The unfolded nature of $\alpha$-syn allows it to function as a molecular chaperone assisting in the transport of vesicles and facilitating SNARE formation [48]. This function could extend to the adhesion of nucleotides to $\alpha$-syn, and could also suggest that RNA has a naturally greater affinity for $\alpha$-syn and might be a viable alternative to the DNA aptamers raised here and experimentally. This has been observed experimentally, with $\alpha$-syn being capable of binding its own mRNA to prevent translation of further protein [49]. Using catRAPID, an algorithm used to predict protein/RNA interactions, demonstrated $\alpha$-syn easily binds to its own mRNA without the need for an intrinsic riboswitch sequence within the mRNA [50].

When considering the structures of the aptamers, it is important to consider the presence and size of key structural motifs within the aptamer. Loops, both external and internal, tend to be much more flexible in nature when compared to standard Watson and Crick paired sections or other loops such as stems and bulges [51]. This inherent flexibility, combined with the relatively small size of aptamers in three-dimensional space, allows for the aptamer to integrate with regions of proteins that larger more complicated structures, like antibodies, are simply incapable of achieving [14]. This has been demonstrated previously by Autiero et al. [52] when looking into RNA aptamer binding to S8 ribosomal proteins. Here the ability of the aptamer to adjust its structure to bind to the target not only improves the binding efficiency of the aptamer but is also a critical property of aptamers that warrants further investigation.

When looking at the reference aptamer and the top aptamer generated as part of this study, displayed in Fig. 7, the most notable differences between them are the extension of the stem loop between C7 and G25 and the incorporation of a greater number of Watson and Crick pairing between G27 and C48 in TMG-79. This increase in bulge size, and the incorporation of unpaired G44 and A45, gives greater flexibility to this extended region in comparison to the reference aptamer’s similar structure between C26 and G39. An additional consideration is the increase in the number of potential binding points which have been predicted within TMG-79. This may be due to the increased flexibility of the aptamer compared to the reference aptamer which, while not being a rigid structure, has a greater number of structural motifs limiting its flexibility [52]. It is also important to remember that, while the increase in flexibility is important, the structural motifs themselves provide unique binding properties which allow for aptamer/protein interactions. However, there is evidence to suggest that both stable and unstable secondary motifs are important in aptamer binding. Adachi et al. [53] demonstrated that their RNA aptamer’s inherent instability actually promoted tighter binding with higher specificity compared to aptamers which had more stable structures. This conflicts with Troisi et al. [54], where the structural motifs interact in a specific fashion to promote greater affinity for thrombin of NU172. The difference between these studies and their approach to the development of the aptamer may be down to either differences in complexity of the proteins or their specific locations. In both examples presented, they are looking at relatively complex proteins, either thrombin [54] or interleukin 17 receptor A/F (IL17A/F) inhibition [53]. Within the second example, the hope is to inhibit formation of IL17A/F, requiring only to be able to bind to one component to competitively inhibit this function. In contrast, Troisi et al. [54] identified sequence specific binding to the thrombin protein which could be significantly reduced based on single amino acid mutations. The conclusion that can be drawn from this is that location, structure and the complexity of the proteins involved dictates the requirement for either rigidity or flexibility within the aptamer’s structure which influences its ability to bind to its target. Considering the inherent properties of $\alpha$-syn being a mainly unstructured protein with the exception of the two alpha helix motifs, it is interesting to see that both the reference aptamer and TMG-79 both bind within this alpha helix “sandwich”. Incorporation into this region and the increased affinity for the target observed with TMG-79 suggests that, while the protein may be mainly in an unstructured state, this is a critical recognition element for biomarker detection. While this data is promising, the requirement to get a solved crystal structure for the interaction is necessary to confirm the interactions predicted as part of this work.

Fig. 7.

Binding interaction points highlighted in red. (A) M5-15. (B) TMG-79. (C) While there are similarities between the two aptamer structures (51.4% sequence homology), the notable difference is the lack of stem loop between G43 and C58 within the reference aptamer. This change results in a larger bulge region in TMG-79, allowing for increased flexibility of the overall aptamer structure. Top pose predictions of 4ts2 (left) and M5-15 (right) against $\alpha$-syn are shown.