Computer-aided design of aptamers for SMMC-7721 liver carcinoma cells

Introduction

Aptamers generally have 15–60 nucleotide bases in length and can fold into complex three-dimensional conformations and bind targets, including molecules, cancer markers and cells, with high affinity and high specificity [1], [2]. This binding reaction enables their application not only as receptors in the analysis of food and environment, but also as new molecular medicines, or as molecular probes for clinical diagnostics and related applications. Compared to other molecular probes, such as antibodies, ligands, enzyme substrates, nucleic acid aptamers are relatively new and promising. They not only have many inherent advantages in stability and facility of generation and synthesis, but also exhibit fast blood clearance and rapid tissue and tumor penetration [3], [4], [5]. Aptamers are typically selected from libraries of random DNA (or RNA) sequences with systematic evolution of ligands by exponential enrichment (SELEX) [6], [7]. SELEX process partitions bound candidate aptamers from non-binding sequences by means of an affinity method. This technique amplifies the bound sequences by polymerase chain reaction (PCR). In cell-SELEX, the whole cells are used as targets [8], [9]. This method has been used to obtain aptamers identifying molecular differences among cancer cells [10], [11].

Generally, the SELEX experiments need 10–15 selection rounds to be completed. Thus this technique is relatively labor-intensive, time-consuming and expensive. Structure-activity relationships (SAR) models including the pattern recognition (classification) can be used to reveal the relationships between chemical structures and properties. According to classical chemical theory, it is the molecular structure that predetermines the chemical, biological, and physical properties. Therefore, these properties may be elucidated from the structural information, which can be reflected with descriptive parameters (descriptors). SAR models can be used to screen a series of molecules, including those not yet synthesized, on the computer in order to select chemical structures with the properties desired. It is then possible to synthesize the most-promising candidates for laboratory testing [12]. Thus, computer-aided aptamer design based on pattern recognition approach can accelerate the process of selecting DNA (or RNA) sequences for use as aptamers or for any other purpose. However, as we know, there are only a few SAR models reported on aptamer sequences [13], [14], [15]. There may be two reasons for this phenomenon. One is that experimental data such as affinity for aptamers are insufficient. The other is that aptamer sequences have large molecular weight, which leads to difficult in calculating descriptors.

To design and synthesize aptamers for hepatocellular carcinoma (HCC) with high sensitivity and specificity would be of great importance for scientific research, and clinical diagnosis and treatment of cancer [11]. Schütze et al. [16] found that reverse-complement sequences can recognize targets. The aim of this work is to develop a pattern recognition model (or SAR model) for aptamers against SMMC-7721 liver carcinoma cells, using molecular descriptors from reverse-complement sequences as the input variables of the support vector machine (SVM) classification model. Subsequently, according to computer-aided aptamer design, eight DNA aptamer sequences are synthesized and tested.

Materials and methods

Candidate aptamer sequences selected for SMMC-7721 liver carcinoma cells through whole cell-SELEX were taken from our laboratory (see Table S1 in the Supplementary data). After deleting the sequences whose secondary structure containing no pair, 226 reverse-complement sequences (5′-GATGACGACCGACTGACTTC – (center sequence) – GGTAGTCAGTGGCAAGTCAA-3′) were obtained and split into two sets: a training set (60 reverse-complement sequences) and a test set (166 reverse-complement sequences).

The training set was randomly taken from the reverse-complement sequences of the 3rd, 11th and 13th round of selection. The class labels (i.e. target values) of aptamers from the 3rd, 11th and 13th rounds of SELEX were 1, 2, and 2, respectively. The class labels 1 and 2, respectively, denote the high and low affinity and specificity aptamers. The training set that has 35 sequences with the class labels being 1 and 25 sequences with the class labels being 2 was used to train and optimize pattern recognition models. The test set including reverse-complement sequences from the 3rd, 5th, 7th, 9th, 11th, and 13th rounds of experiments was used later to check the predictive power of the developed model.

During the first few SELEX rounds, the binding affinity K_d increases significantly [17], [18]. When a certain number of rounds are performed, the affinity reaches a plateau. Therefore, it is reasonable to set class labels of aptamer sequences from the 3rd, 11th and 13th rounds of SELEX as 1, 2, and 2, respectively.

RNA structure software version 5.3 [19] was used to calculate the free energy (E) and obtain the secondary structures of center sequences of reverse-complement sequences. The loop of each secondary structure that has smaller hairpin structures and larger stems was used to calculate molecular descriptors. Each loop structure was sketched using ChemBioDraw Ultra 11.0 [20], and optimized using ChemBio3D Ultra 11.0 [20] with the default convergence criteria in the soft. Each energy minimized molecule was then input for Dragon software [21] by which 1664 molecular descriptors categorized by 20 blocks were calculated. It should be noted that ChemBio3D could not be used for determining the absolute or precise 3D structures of such sophisticated molecules as these in this work. However, the errors due to the approximate nature of ChemBio3D method are largely transferable within structurally related series [12]. Thus, relative structures optimized with ChemBio3D can be meaningful even though their absolute and precise structures are not directly applicable.

Besides above the free energy (E) and the 1664 molecular descriptors, the percentages for adenine (A) (A%), guanine (G) (G%), cytosine (C) (C%), and thymine (T) (T%) in each center sequence were calculated. Totally, 1669 molecular descriptors were calculated in this work for each aptamer.

One of the 20 descriptor blocks, 3D-MoRSE (3D Molecule Representation of Structures based on Electron diffraction) descriptors consist of 160 molecular descriptors. Our pervious study suggests that 3D-MoRSE descriptors are important in reflecting the structure information [15]. Three dimensional-MoRSE descriptors are calculated with the following expression [20]:

where Morsw is the scattered electron intensity, w is an atomic property, i.e. the unweighted case (u), atomic mass (m), the van der Waals volume (v), the Sanderson atomic electronegativity (e) or the atomic polarizability (p). r_ij is the interatomic distances and nAT is the number of atoms. The term s represents the scattering in various directions by a collection of nAT atoms and takes integer values in the range 0–31.

Results and discussion

The rounds of SELEX as the dependant variable, 1669 molecular descriptors as the independent variables, stepwise multiple linear regression (MLR) was adopted to select the optimal subset of variables used to develop the SVM classification model from the training set. By investigating the correlation between the 1669 descriptors and the rounds of SELEX with SPSS 11.5, an optimal descriptor subset, E, Mor26e (3D-MoRSE - signal 26 / weighted by atomic Sanderson electronegativities), Mor32p (3D-MoRSE - signal 32 / weighted by atomic polarizabilities), A%, and G%, were obtained to describe the structure features of each aptamer, which are listed in Table S2 (see Supplementary data). The two 3D-MoRSE descriptors, Mor13m and Mor32p, retain important structural features of electronegativities and polarizabilities, respectively [15], [21]. Thus the 3D-MoRSE descriptors used denote the electronic distribution of loops, which is related to the induced fitting behavior and molecular recognition of aptamers. The free energy (E) is related to the stem length of the secondary structure of center reverse-complement sequences and is correlated with conformational stability of the secondary structure.

Table S2 indicates that an aptamer with a higher A% (or lower G%) value would possess high affinity and selectivity for SMMC-7721 liver carcinoma cells. According to classical theory, entropy is commonly associated with the degree of disorder in a system. During the SELEX process, the amount of disorder (or the entropy) in a center sequence gradually decreases. This means that the fraction of a certain base would increase or decrease. On the other hand, the bases A and G have different chemical composition and structure that may contribute to molecular recognition, such as shape, electrostatics, dynamics and entropy, and result in difference with respect to the affinity and selectivity of aptamers. Therefore, the descriptors A% and G% are related to the binding affinity of aptamers to targets.

The optimal descriptor subset (E, Mor26e, Mor32p, A%, and G%) was then used to develop a support vector classification (SVC) model, with Gaussian radial basis function (RBF) as the kernel function. The SVM parameters C and γ were optimized with the particle swarm optimization (PSO) method [22]. During the optimization process, both of the cognition learning factor c₁ and the social learning factor c₂ were set as 2. The number of particles and the maximum number of iterations were set to be 20 and 200, respectively. The searching range of SVM parameter C was between 0.1 and 1000. The searching range of γ was 0.1–1000. The leave-one-out (LOO) cross-validation procedure was carried out for the training set during the optimization of SVM parameters C and γ. LibSVM package (Version 3.12) [23] was used to develop the SVC model.

The optimization results show that the relatively optimal SVM parameters C and γ were 4.6 and 0.1, respectively. The SVC model with C=4.6 and γ=0.1 has classification accuracy of 98.3% for the training set. To evaluate the model, we calculated prediction values for the test set, which are listed in Table S2. For the classification parameters of the 3rd, 11th, and 13th round experiments, specificity and sensitivity are 90.91% and 84.91%, respectively. The fractions of winner aptamers with high affinity and specificity against targets for the 3rd, 5th, 7th, 9th, 11th, and 13th round experiments are 0.09, 0.17, 0.69, 0.84, 0.90 and 0.98, respectively. Thus, we can obtain the following fitting curve and the corresponding exponential equation (see Figure 1). Obviously, the results are in agreement with the aptamer evolutionary principles of SELEX based screening, as stated above.

Figure 1:

Evolutionary curve of aptamers for SMMC-7721 liver carcinoma cells.

Table S2 shows that the reverse-complement sequences of Nos. 53, 59, 60, 224, 225, and 226 have prediction class labels being 2, and belong to the class of sequences with high binding affinity (effective equilibrium dissociation constant, K_d). Their experimental K_d values are 7.28 nM, 16.44 nM, 14.28 nM, 1.47 nM, 29.79 nM, 5.34 nM, respectively. The theoretical predictions conform to experimental results.

To further evaluate the SVC models, we designed, synthesized and tested eight DNA candidate aptamer sequences with HCC cell line SMMC-7721 as target cells. The experimental conditions can be found in Reference [11]. Figure 2 is the experimental curves of K_d for new designed sequences measured by Becton Dickinson FACScalibur flow cytometer. According to the single binding site scheme [24]: B=(B_max×C)/(K_d+C), where C is the concentration of candidate aptamer sequences, B is the aptamer adsorption (i.e. the relative Geo Mean values), B_max is the maximal adsorption of candidate aptamer sequence at high aptamer concentrations. The experimental K_d values of eight new designed sequences were obtained and listed in Table 1. Table 1 shows that these new sequences belong to the class of sequences with high binding affinity and have experimental dissociation constants K_d in the nanomolar range. The results show that the eight new sequences can recognize the target SMMC-7721 cells with high binding affinity and may serve as effective tools for early diagnosis of hepatocellular carcinoma. The feasibility of applying SAR model based on reverse-complement sequences to design aptamers has been demonstrated.

Figure 2:

The experimental K_d curves of new designed sequences.

Conclusions

By means of the SVC technique, a pattern recognition model was constructed for two-class candidate aptamer sequences with high or low affinity against SMMC-7721 liver carcinoma cells. The five descriptors, E, Mor26e, Mor32p, A%, and G%, obtained from reverse-complement sequences were used to describe the structure features of each aptamer and reflect the properties of aptamers-targets interactions. According to the predictions from the SVC model, eight new DNA aptamer sequences with HCC cell line SMMC-7721 as target cells were designed, synthesized and tested. The experimental K_d values of eight new designed sequences are in the nanomolar range. The results show that the eight new sequences can bind the target SMMC-7721 cells with high binding affinity. These results encourage the further application of pattern recognition models based on reverse-complement sequences to other designs of aptamers.