Comparative Molecular Field Analysis (CoMFA) and Comparative Molecular Similarity Indices Analysis (CoMSIA) Studies on α1A-Adrenergic Receptor Antagonists Based on Pharmacophore Molecular Alignment

The α1A-adrenergic receptor (α1A-AR) antagonist is useful in treating benign prostatic hyperplasia, lower urinary tract symptoms, and cardiac arrhythmia. Three-dimensional quantitative structure-activity relationship (3D-QSAR) studies were performed on a set of α1A-AR antagonists of N-aryl and N-nitrogen class. Statistically significant models constructed from comparative molecular field analysis (CoMFA) and comparative molecular similarity indices analysis (CoMSIA) were established based on a training set of 32 ligands using pharmacophore-based molecular alignment. The leave-oneout cross-validation correlation coefficients were q2 CoMFA = 0.840 and q2 CoMSIA = 0.840. The high correlation between the cross-validated/predicted and experimental activities of a test set of 12 ligands revealed that the CoMFA and CoMSIA models were robust (r2 pred/CoMFA = 0.694; r2 pred/CoMSIA = 0.671). The generated models suggested that electrostatic, hydrophobic, and hydrogen bonding interactions play important roles between ligands and receptors in the active site. Our study serves as a guide for further experimental investigations on the synthesis of new compounds. Structural modifications based on the present 3D-QSAR results may lead to the discovery of other α1A-AR antagonists.


Introduction
Adrenergic receptors are members of the G-protein coupled receptor superfamily of membrane proteins that mediate effects of the sympathetic nervous system through the actions of epinephrine and norepinephrine and control the equilibrium of the cardiovascular system. The receptors have been classified into three classes (α 1 , α 2 , and β) with three members in the α 1 subfamily, namely, α 1A , α 1B , and α 1D [1]. Many physiological processes, including smooth muscle contraction, myocardial inotropy and chronotropy, and hepatic glucose metabolism are related to the α 1A -adrenergic receptors (α 1A -ARs); thus, clinical applications of these receptors are of interest to many scientists. The α 1A -AR antagonists are considered as the first-line of therapy for lower urinary tract symptoms associated with clinical benign prostatic hyperplasia (BPH), because they have proven efficacy in mediating smooth muscle contraction in the human prostate [2].
Reports on the application of QSAR analysis to α 1 -AR species are relatively limited [3][4][5][6]. The SAR of 39 α 1 adrenoceptor antagonists derived from the antipsychotic sertindole with respect to affinity was investigated using a 3D-QSAR approach based on the GRID/GOLPE methodology by Balle's group [7]. Montorsi reported quantitative size and shape affinity/subtype selectivity relationships derived from a large set of very heterogeneous α 1a -, α 1b -, and α 1d -adrenergic receptor antagonists [8]. Maciejewska et al. reported a structure-activity analysis of hexahydro and octahydropyrido-(1,2-c)-pyrimidine derivatives as α 1A -AR antagonists, using comparative molecular field analysis (CoMFA) [9]. Li's and Xia's groups provided insights into α 1A -adrenoceptor antagonists through self-organizing molecular field analysis in the same year [10]. Shakya et al. reported 3D-QSAR models based on a series of 25 dihydropyridine class compounds using the APEX-3D program, which can automatically identify biphore (pharmacophore) sites, and 3D-QSAR modeling [11]. Compared with SOMFA studies, the model has strong predictability, particularly for compounds with enantio-selectivity. Last year, Gupta reported on studies based on CoMFA and comparative molecular similarity indices analysis (CoMSIA) performed on a set of diverse α 1A -AR antagonists to understand the relationship between structure and antagonistic activity [12]. This study was based on common structural alignments. The generated models suggest that steric, electrostatic, and hydrophobic interactions play important roles in structure-activity analysis.
The superimposition of compounds in an alignment is crucial in 3D-QSAR studies. Many studies have reported that pharmacophore alignment is a useful tool that is superior to others [13,14]. Genetic algorithm with linear assignment of hypermolecular alignment of datasets (GALAHAD ® ) is regarded as a superior tool for molecular alignment compared with classical common structural alignment [15], especially for compounds that share few structural commonalities. GALAHAD uses the proprietary technology of Tripos ® and generates pharmacophore alignments and hypotheses from sets of ligand molecules using a genetic algorithm [16]. The relationship between the chemical structures and the biological function studies based on the pharmacophore alignments of α 1A -adrenergic receptor using GALAHAD has yet to be reported. Based on the alignment obtained, CoMFA and CoMSIA models were developed. The information derived from these models will be helpful in predicting the activity and guiding the design of other α 1A -AR antagonists.

Data Collection
All α 1A -AR antagonist structures were collected from recent reports ( Figure 1 and Table 1) [17][18][19][20][21][22][23][24][25]. The selected compounds have diverse structural features and cover a wide range of biological activity, spanning over four orders of magnitude (0.1-630 nM). The training and test sets were classified to ensure that both sets could completely cover the whole range of biological activity and physicochemical and structural diversity studied. The need for models to be tested on a sufficiently large test set (25%-33% of the total samples) to establish a reliable QSAR model had been previously suggested [26]. Therefore, the dataset, which was composed of 44 compounds, was divided into training (32 compounds) and test (12 compounds) sets. All data from various references were from the same laboratory. The affinity constants of compounds were evaluated by radio-receptor binding assays. [ 3 H] Prazosin was used to label cloned human α 1 -ARs expressed in CHO cells and WB4101 was used as a reference compound. K i values were derived from IC50 using the Cheng-Prusoff equation [27]. All

Structural Sketch and Alignment
Structural sketches and refinement of the entire set of α 1A -AR antagonists were accomplished using SYBYL ® 8.1 molecular modeling software (Tripos Associates Inc.) and their 3D structures were generated using CONCORD ® . All compounds were minimized under the Tripos standard (TS) force field [28] with Gasteiger-Hückel atomic partial charges [29]. Minimizations were done using the Powell method, in which calculations were set to terminate at an energy gradient value of 0.01 kcal/mol.
Pharmacophore results for N-Aryl and N-heteroaryl piperazine α 1A -AR antagonists using GALAHAD have been reported previously [30], and the optimized model from the previous report was used as a template in the present study. All compounds in the training and test sets were aligned with this template using the "Align Molecules to Template Individually" option in GALAHAD. The other parameters for calculation were set to default values.

3D-QSAR Studies
CoMFA, a method that reflects the non-bonding interaction between the receptor and the ligand, is widely used in drug design. The steric (Lennard-Jones) and electrostatic (Coulombic) potential energies of the TS force fields implemented in SYBYL were evaluated by CoMFA. For each pharmacophore alignment ligand, a 3D cubic lattice with a grid spacing of 1.0 Å in the x, y, and z directions was generated to enclose the molecule aggregate. A sp 3 carbon atom with a charge of +1.0 and a Van der Waals radius of 1.52 Å was used as a probe; this atom was placed at every lattice point to calculate various steric and electrostatic fields. An energy cut off value of 30 kcal/mol was imposed on all CoMFA calculations to avoid excessively high and unrealistic energy values within the molecule. Then, partial least-squares (PLS) analysis was applied to obtain the final model [31].
During calculation of the steric and electrostatic fields in CoMFA, many grid points on the molecular surface were ignored due to the rapid increase in Van der Waals repulsion. To avoid a drastic change in the potential energy of the grid points near the molecular surface, CoMSIA employed a Gaussian-type function based on distance. Thus, CoMSIA may be capable of obtaining more stable models than CoMFA in 3D-QSAR studies [31][32][33]. The constructed CoMSIA model provided information on steric, electrostatic, hydrophobic, hydrogen bond donor, and hydrogen bond acceptor fields.
The grid constructed for the CoMFA field calculation was also used for the CoMSIA field calculation [32]. Five physico-chemical properties (electrostatic, steric, hydrophobic, and hydrogen bond donor and acceptor) were evaluated using a common probe atom placed within a 3D grid. A probe atom sp 3 carbon with a charge, hydrophobic interaction, and hydrogen-bond donor and acceptor properties of +1.0 was placed at every grid point to measure the electrostatic, steric, hydrophobic, and hydrogen bond donor or acceptor field. Similar to CoMFA, the grid was extended beyond the molecular dimensions by 1.0 Å in three dimensions and the spacing between probe points within the grid was set to 1.0 Å. Different from the CoMFA, a Gaussian-type distance dependence of physicochemical properties (attenuation factor of 0.3) was assumed in the CoMSIA calculation.
The partial least squares (PLS) method was used to explore a linear correlation between the CoMFA and CoMSIA fields and the biological activity values [34]. It was performed in two stages. First, cross-validation analysis was done to determine the number of components to be used. This was performed using the leave-one-out (LOO) method to obtain the optimum number of components and the corresponding cross-validation coefficient, q 2 [35]. The value of q 2 that resulted in a minimal number of components and the lowest cross-validated standard error of estimate (S cv ) was accepted. The column filtering values (σ min ) was set to 2.0 kcal/mol in order to speed up the analytical process and reduce noise. Second, the optimum number of components were used to derive the final PLS model, with no validation method [36,37]. The CoMFA and CoMSIA results were then graphically interpreted by field contribution maps.

Predictive Power of the Models
The predictive power of the 3D-QSAR model was determined from a set of 12 compounds in the test set (Table 2). These molecules were superimposed using the same pharmacophore molecular alignment method described above, and their activities were predicted using the CoMFA and CoMSIA models generated by the training set. The predictive correlation coefficient, r 2 pred , of the CoMFA and CoMSIA models were calculated using the test set, according to the formula, where SD is the sum of the squared deviations of each experimental value from the mean, and PRESS is the sum of the squared differences between the predicted and actual affinity values for every molecule [31].

Alignment
The pharmacophore model generated by GALAHAD contained one donor center, one positive charged nitrogen atom, two acceptor centers, and two hydrophobic centers, more than other reported models, which is shown in [Figure 2(a)]. The internal distance of the pharmacophore feature is consistent with previously reported models. The pharmacophore model has reasonable assessment parameters, such as a high specificity value (4.212), high steric score (4903.5), and low energy value (92.9 Kcal/mol) [30,38]. All of the compounds were superimposed based on the template [ Figure 2(b)]. Ligands have the same binding site and share common interactions. However, they sometimes do not share common structures. GALAHAD, as an alignment tool, considers the flexibility of each compound and performs a compact alignment based on essential structural elements related to activity by variation of the torsional degrees of freedom. Compared to the QSAR model previously reported, the classical approach based on a rigid alignment of minimized structures, does not show any uniformity with the active conformation. Pharmacophore-based molecular alignment using GALAHAD in the current manuscript will generate a rational alignment for 3D-QSAR studies.

CoMFA and CoMSIA Models
A statistically significant 3D-QSAR model was obtained using the properly selected training set of 32 ligands. Results of the statistical analysis are presented in Table 3.
In the CoMFA model, initial PLS analysis of the aligned training set was done using a default σ min data filter of 2.0 kcal/mol and the Tripos standard field. This yielded a highly significant q 2 value of 0.840 (with S cv = 0.476, using four components), which indicates that it is a model with high statistical significance; a q 2 value of 0.6 is considered statistically significant in CoMFA and CoMSIA studies [39]. The conventional r 2 value of 0.988 and low standard error of estimate (SEE) value of 0.128 indicate the accuracy of the predictions of the model. High values of q 2 from the LOO analysis can be regarded as a necessary, but not a sufficient, condition for a model to possess significant predictive power [40]. In addition to LOO, the internal predictive ability of the model was further assessed by a group cross-validation performed with 30 groups for 30 times. The mean of 30 readings was given as r 2 cv ; the mean r 2 cv of 0.837 also suggests that the model has good internal predictability.
Similar to CoMFA, the CoMSIA model with the steric, electrostatic, hydrophobic, donor and acceptor fields result in satisfactory statistical data (Table 3, model A). However, the steric and donor fields (10% and 18% contribution, respectively) did not significantly contribute to affinity in the all-fields model. Moreover, presumably correlated fields increase the complexity of contour maps and complicate their interpretation. Therefore, CoMSIA approaches with smaller subsets of fields were calculated (Table 3, columns B-E). The hydrophobic field was always considered in accordance with the hydrophobic groups in the pharmacophore model. The internal predictivities of the field combinations B-E were only slightly reduced compared to the all-fields model A. Among the five models, the combination of the electrostatic, hydrophobic, and acceptor field (model E) has the highest external predictivity (r 2 pred = 0.671). The conventional non-cross validated r 2 of 0.975 and the SEE value of 0.180 indicate that this model is statistically highly significant. Analogous to the CoMFA, a group cross-validation was done to assess further the internal predictive ability of the model. The cross-validation was performed 30 times with 30 groups. The mean r 2 cv obtained was 0.864, which also indicated the robustness of this model.

Predictive Power of the Models
The cross-validated pK i values calculated by CoMFA and CoMSIA, and the residuals between the experimental and cross-validated pK i values of the compounds in the training set are listed in Table 4. The predictive powers of the CoMFA and CoMSIA models were further examined using a test set of 12 compounds not included in the training set. The predicted pK i values calculated by CoMFA and CoMSIA are also shown in Table 4.

Graphical Interpretation of the Fields
The CoMFA and CoMSIA contour maps of the PLS regression coefficients at each region grid point provide a graphical visualization of the various field contributions, which can explain the differences in the biological activities of each compound. These contour maps were generated using various field types of StDev*coefficients to show the favorable and unfavorable interactions between ligands and receptors in the active site.
In the CoMFA model, the fractions of steric and electrostatic fields are 46.0% and 54.0%, respectively. Favorable and unfavorable cutoff energies were set at the 80th and 20th percentiles for the steric contributions. The contour maps of the fields are shown in [ Figure 4(a)], with the higher affinity compound 20 as the reference structure. The surfaces indicate the regions where the increase (green region) or decrease (yellow region) in steric effect would be important for the improvement of binding affinity. The large green isopleths upon the thiochromene part reflect a sharp increase in affinity for all the anchor moieties transferred into this area. Compound 20, with its large bulky phenyl group, coincide with the green isopleths. When the thiochromene group in compound 20 was replaced by 8-methyl-8-azaspiro decane-7, 9-dione (such as compounds 24 and 32), and the yellow region was occupied by the large bulky groups, and the antagonistic activity of these compounds evidently decreased. This was probably due to the insufficient space to accommodate these bulky groups in the anchor moieties ofthe receptor-binding site, which caused collision among groups and the reduced affinity of the compounds. [ Figure 4(b)] shows the electrostatic contributions with compound 20 as a template ligand. The electrostatic contour map shows regions of red polyhedra (contribution level: 15%), where electron-rich substituents are beneficial for the binding affinity, whereas the blue colored regions (contribution level: 85%) show the areas where positively charged groups enhance the antagonistic activity. The large blue area near the thiochromane moiety indicates a region where negatively charged groups decrease antagonistic activity. Compound 32 has a carbonyl in the blue area, which may be not conducive to the improvement of activity. The electrostatic contour map shows a region of red contours neighbor to the oxygens connects with benzene, indicating that electron-rich substituents (such as bromine, cyano group) are beneficial for the binding affinity. In the CoMSIA model, the fractions of the electrostatic, hydrophobic, and hydrogen-bond donor and acceptor fields were 34.7%, 39.9% and 25.4%, respectively. The CoMSIA contour plots are also exemplified by some ligands of high affinity. Similarly, red and blue isopleths (contribution levels: 15% and 85%, respectively) of the electrostatic fields [ Figure 5(a)] enclose regions, where negative and positive charges have favorable effects on pK i , respectively. The contour maps of electrostatic fields are similar to that of the CoMFA. Compared with the CoMFA model, its electronic field areas are smaller. The blue contours rather indicate where negatively charged substructures are unfavorable, such as methyl substitution and amino substitution should be placed to enhanced binding affinity, while the red contours neighboring to oxygens indicate electronegative groups will do a favorable effect of antagonistic potency. For hydrophobic field contributions, orange and gray isopleths are drawn at contribution levels of 85% and 15%, respectively. These enclose regions favorable for hydrophobic and hydrophilic groups, respectively [ Figure 5(b)]. Orange areas are mainly distributed in the thiochromene side, which is consistent with the hydrophobic center of the pharmacophore hypotheses, as well as with the contour maps of the steric fields generated by the CoMFA and CoMSIA models. Hydrophobic interactions may be dominant for ligand binding. These areas suggest the importance of hydrophobic interaction near the phenylchroman group of the most active compound in the training set. Halogen (except fluorine) substitution in orange region will enhance molecular hydrophobic field and increase the affinity of the compound. In addition, the introduction of hydrophilic groups (such as amide and amine) on the other side of the benzene ring (transparent gray) may increase α 1A -AR antagonistic activity, whereas chlorine substitution in compound 32 may reduce the binding activity of α 1A -AR, which is a representative of weak α 1A -AR receptor binding.
The graphical interpretation of the field contributions of the H-bond properties is shown in [ Figure 5(c)] (the H-bond acceptor field). Magenta isopleths (80%) encompass regions wherein hydrogen bond acceptor groups lead to improved α 1A -AR antagonist activity; while areas encompassed by red isopleths (20%) should be avoided from hydrogen bond acceptor, which will result in impaired biological activity. Small red isopleths are inlaid in the purple region. The oxygen atom connected to the phenyl group may be an important factor to high activity of compound 20. The introduction of a hydrogen bond acceptor in the area may further improve the activity.
3D structures of human α 1A -AR homology models have been successfully developed by Li's group based on the crystal structure of β 2 -AR [41]. This optimized homology model was retrieved from Protein Model Data Base (PMDB entry: PM0075211). The pharmacophore features were placed into the binding pocket of this model, as shown in Figure 6. According to previous reports [41,42], hydrophobic interactions may occur between the hydrophobic group of a ligand and pockets constituted by Trp102, Cys110, Tyr111, Ser158, Phe193, Phe288, Phe289 and Phe312. Asp106 seems to play a role in antagonist binding through strong electrostatic interactions with the protonated nitrogen atom. Ser188 and Ser192 may have H-bond effects with acceptor and donor features. In contrast to the most reported 3D-QSAR models mentioned in "Introduction", GALAHAD was used in the present models to superimpose ligands based on pharmacophore molecular alignment. The flexible superimposition afforded by the current models more closely coincides with the receptor binding sites and reduces the likelihood of confusion arising from compounds belonging to different classes. The generated models suggested that electrostatic, hydrophobic, and hydrogen bonding interactions play important roles between ligands and receptors in the active site. The q 2 , r 2 pred and contours obtained from the present CoMFA and CoMSIA model show strong predictability and application and provide detailed information about the molecular features of the ligands, which will contribute to the antagonistic potency.

Conclusions
CoMFA and CoMSIA models were established based on a training set of 32 ligands using pharmacophore-based molecular alignment by GALAHAD. The statistical significance of the models was evaluated. The q 2 values of the CoMFA and CoMSIA models are 0.840 (4 components) and 0.840 (3 components), respectively. These models also predict the biological activities of 12 ligands of α 1A -AR antagonists in the test set.
In this study, GALAHAD is a useful pharmacophore alignment tool that can yield a good 3D-QSAR model for α 1A -AR antagonists. The present results of the 3D-QSAR model provided a valuable tool for predicting the activities of novel α 1A -AR antagonists and a basis upon which more active compounds can be derived by targeted structural modifications.