- Email: [email protected]
3D-QSAR study of the phenylsulfamic acid derivatives as HPTPβ inhibitors
Journal of Molecular Structure 1186 (2019) 11e22
Contents lists available at ScienceDirect
Journal of Molecular Structure journal homepage: http://www.elsevier.com/locate/molstruc
3D-QSAR study of the phenylsulfamic acid derivatives as HPTPb inhibitors Wenjuan Zhang a, b, Zhao Wei c, Chunyu Lin a, Zhonghua Wang d, Zhibing Zheng b, *, Song Li a, b, ** a
School of Pharmaceutical Engineering, Shenyang Pharmaceutical University, Shenyang, 110016, China National Engineering Research Center for the Emergency Drug, Beijing Institute of Pharmacology and Toxicology, Beijing 100850, China Department of Medicinal Chemistry, School of Pharmacy, Fourth Military Medical University, Xi'an, 300071, China d School of Chemical and Environmental Engineering, Shanghai Institute of Technology, Shanghai, 201203, China b c
a r t i c l e i n f o
a b s t r a c t
Article history: Received 16 January 2019 Received in revised form 25 February 2019 Accepted 25 February 2019 Available online 28 February 2019
Human protein tyrosine phosphatase beta (HPTPb) inhibitors have been used for the treatment of sepsis, cancer and inflammatory diseases. The phenylsulfamic acid derivatives are a group of powerful inhibitors of HPTPb. To explore the relationship between their structures and biological activities, the threedimensional quantitative structure-activity relationship (3D-QSAR) model is first constructed. Two highly predictive 3D-QSAR models are established. These models are the comparative molecular field 2 0.611, R2ncv 0.999) model and the comparative molecular similarity index analysis analysis (CoMFA: qcv 2 0.588, R2ncv 0.993) model. The results show a quite good external predictive power for the (CoMSIA: qcv test set, with R2pre values of 0.833 and 0.775, respectively. Furthermore, the contour maps of the 3D-QSAR models are analysed together with the results of molecular docking. The analysed results are conducive to discovering new binding sites and provide a reference for the construction of new potent HPTPb inhibitors compounds. © 2019 Elsevier B.V. All rights reserved.
Keywords: 3D-QSAR HPTPb Phenylsulfamic acid Molecular docking
1. Introduction Human protein tyrosine phosphatase beta (HPTPb) was discovered by Krueger in 1990 [1]. HPTPb is a classical receptor-like protein tyrosine phosphatase (PTP) and is specifically expressed in endothelial cells [2]. The mouse orthologue of HPTPb is called vascular endothelial protein tyrosine phosphatase (VE-PTP). HPTPb is involved in the angiopoietin (ANG)etyrosine kinase receptor (TIE) pathway. Angiopoietin 1 (ANG1) is a potent TIE2 activator whereas Angiopoietin 2 (ANG2) is a context-dependent agonist/ antagonist [3]. One of the physiological functions of activated TIE2 is to maintain the vascular endothelial integrity and reduce the leakage of blood vessels [4]. When HPTPb is inhibited, the rate of dephosphorylation of TIE2 is decreased, which is equivalent to either triggering TIE2 or to the result of the amplification of the ANG1 signal. Therefore, HPTPb inhibitors can be used for treating
* Corresponding author. ** Corresponding author. School of Pharmaceutical Engineering, Shenyang Pharmaceutical University, Shenyang, 110016, China. E-mail addresses: [email protected] (Z. Zheng), [email protected] (S. Li). https://doi.org/10.1016/j.molstruc.2019.02.107 0022-2860/© 2019 Elsevier B.V. All rights reserved.
sepsis, organ transplantation, atherosclerosis and vascular complications of diabetes by the activation of TIE2 [5]. Peters et al. identified many sulfamic acid derivatives such as compound 1 by using high-throughput screening to screen the P&GP corporate repository against multiple phosphatase enzymes. However, these compounds displayed weakly inhibitory effect on PTP1B, and the value of the best affect was 322.5 mM against PTP1B [6]. Upon structural optimization, the sulfamic acid moiety [7] was directly linked to an aryl group to mimic the PTP substratephosphotyrosine (pTyr). The inhibitory effect of phenylalanine analogue 2 on PTP1B (39.6 mM) was improved. In subsequent research, these phenylsulfamic acid compounds were also found to exhibit great inhibition activity for HPTPb. Subsequently, they regarded HPTPb as a target that previously attracted little research attention. Therefore, a series of compounds such as 3 with good activity against HPTPb were designed and synthesized, [8] and a potent, non-peptidic phenylsulfamic acid derivative was identified as a clinical candidate for treating diabetic macular edema (DME) [9] (see Fig. 1). 3D-QSAR is based on using the three-dimensional structure information and the biological activity of compounds to establish a reasonable mathematical model. The generated model is used to
12
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
Fig. 1. Development of phenylsulfamic acid derivatives.
determine the relationship between the structures and bioactivities. Furthermore, the model is also used to understand the mode of ligand interaction with the receptor. Comparative molecular field analysis (CoMFA) [10] and comparative molecular similarity indices analysis (CoMSIA) [11] are two common methods for building the models. In CoMFA and CoMSIA, the interaction between the receptor and the ligand is considered to be dependent on the differences in the molecular fields around the compounds. Therefore, the quantified molecular fields' parameters are regarded as independent variables, and the regression analysis of bioactivity is performed to represent the interaction between the ligand and the receptor. The CoMSIA is the further development of the CoMFA. The Lennard-Jones and Coulomb potential functions of the CoMFA are replaced with a distance-related Gaussian function. The function can effectively avoid the steepness of the potentials close to the molecular surface in the CoMFA [12,13]. The steric field, the electrostatic field, the hydrophobic field and the H-bond donor/ acceptor field are considered in the CoMSIA. Although a more stable 3D-QSAR model can be obtained by the CoMSIA, the results may not be better than the CoMFA. To accurately understand and use the 3D-QSAR model, two models were derived to complement and verify each other [14]. Here, 54 phenylsulfamic acid derivatives with a range of HPTPb inhibition activities were selected from the literature [15] as a set. The set was used to establish a reliable 3D-QSAR model by CoMFA and CoMSIA methods carried out using the SYBYL-X2.0 software (Tripos Associates, St. Louis, MO).
2. Methods and materials 2.1. The set and biological activity The selected 54 phenylsulfamic acids in the set had diverse structural features and showed a wide range of HPTPb inhibitory activities. Forty-one compounds were assigned to the training set for the model generation, and 13 compounds were selected as the test set for the model validation. The IC50 values of the HPTPb inhibitor were converted into the corresponding pIC50 (-logIC50) values. The structures and their actual or predicted pIC50 values of the compounds are listed in Table 1.
2.2.2. Acquiring the active configurations of compounds Obtaining accurate active conformations of the molecules corresponding to their bioactivities is a key step in the process of the 3D-QSAR model construction. Molecular docking is used to imitate and analyze the mode of the interaction between the ligand and the receptor. Using the docking results, the efficient active conformations of molecules can be determined. Therefore, the lowest energy conformations of the compounds obtained in the last step were collected and docked with the HPTPb protein one by one. The other parameters were set to default values in the SYBYL-X 2.0 program. The maximum number of poses per conformation was 20. Among the large number of docking poses, the pose with the highest C_Score value was chosen as the active conformation. When the highest C_Score values were the same, the pose with the highest Total_Score value was selected as the active conformation. All of the active conformations were collected and saved as the training set. The protein structure was obtained from the Protein Data Bank (PDB ID: 2I4H) [19]. The protein was treated by deleting all water molecules, extracting the ligand, adding hydrogen atoms to the protein and repairing the skeleton of the protein prior to its use. The ligand of the 2I4H complex underwent the same process to obtain the global lowest energy conformation. The acquired conformation of the ligand was re-docked into the HPTPb rigid protein by SurflexDock to validate the reliability of this docking method. As shown in Fig. S1 (supporting information), the retrieving pose of the ligand could overlay with the original confirmation and could interact properly with the receptor. 2.3. Compounds alignment Molecular alignment is crucial for the establishment of a 3DQSAR model [20]. Thus, two alignment methods were investigated in this study. The first alignment method was based on the receptor. The chosen docking poses were directly used for forming the 3DQSAR model. The result of this alignment is shown in Fig. 2A. The second alignment method was the template ligand-based alignment. The most potent compound was regarded as the template, and the 4-ethyl-phenylsulfamic acid was chosen as the common skeleton to be aligned by the “align database” function. The result is shown in Fig. 2B. 2.4. 3D-QSAR studies
2.2. Acquisition of molecular configurations 2.2.1. Global minimum energy configurations of compounds To explore the minimum energy conformations, all of the compounds were first optimized using the SYBYL-X 2.0 program by the Minimize with the Powell method, Tripos Force Field [16] and Gasteiger-Hückel charge [17] to obtain the local low energy conformation of the molecule. The calculation parameters are listed in Table S1 (supporting information). Then, these molecular conformations were continuously calculated by the simulated annealing [18] method with molecular dynamics to obtain the global minimum energy conformations. The calculation parameters are shown in Table S2.
The CoMFA and the CoMSIA are currently the two most widely used 3D-QSAR methods and differ mainly in their use of different fields to characterize the non-covalent interactions between the ligand and the receptor. To obtain the CoMFA and CoMSIA descriptor fields, the compounds were located in a 3D cubic lattice with the grid spacing of 2 Å in the x, y and z directions. A sp3 carbon atom probe with a van der Waals radius of 1.52 Å and a charge of þ1.0 wanders in this region. In the CoMFA method, the LennardJones and Coul potential energy functions were used to calculate the steric field energy and electrostatic field energy. In the CoMSIA method, a distance-dependent Gaussian function was introduced and used to calculate the distance between the probe atom and
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
13
Table 1 Compound structure and their actual and predicted pIC50 values.
IC50a(nM)
pIC50b
pIC50c
pIC50d
1*
6.000
2.222
2.526
3.205
2
1.000
3.000
2.997
3.233
3
0.2000
3.699
3.684
3.627
4
0.00005
7.301
7.250
7.144
5
0.1000
4.000
4.168
4.162
6
0.0800
4.097
4.076
4.222
7
1.000
3.000
2.996
3.044
8
0.2000
3.690
3.758
3.631
9
0.0020
5.699
5.518
5.235
No.
R1
R2
14
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
10
0.00005
7.301
7.365
7.423
11
1.000
3.000
2.945
3.276
12*
0.4000
3.398
3.300
3.005
13*
0.2000
3.699
4.415
5.230
14*
0.0800
4.097
4.209
5.984
15
0.0007
6.155
6.143
6.146
16
0.00005
7.301
7.342
7.207
17*
0.00005
7.301
6.399
7.881
18
1.000
3.000
2.963
2.949
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
15
19*
0.3000
3.523
4.032
3.772
20
20.00
1.699
1.731
1.699
21*
3.000
2.523
2.029
2.583
22
0.3000
3.523
3.532
3.543
23
0.3000
3.523
3.528
3.472
24
0.00005
7.301
7.331
7.463
25
49.00
1.310
1.281
1.349
26
85.00
1.071
1.068
0.894
27
266.0
0.575
0.533
0.719
28
584.0
0.234
0.233
0.274
29
113.0
0.947
0.944
0.703
30
110.0
0.959
0.957
0.900
31
138.0
0.860
0.862
0.737
16
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
32
98.00
1.009
1.091
1.440
33*
27.00
1.569
0.809
2.278
34
180.0
0.745
0.644
0.515
35
644.0
0.191
0.300
0.486
36
132.0
0.879
0.915
0.822
37
555.0
0.256
0.202
0.206
38
253.0
0.597
0.608
0.335
39
12.00
1.921
1.936
2.059
40
28.00
1.553
1.587
1.614
41
75.00
1.125
1.097
1.233
42*
56.00
1.252
1.752
1.777
43
40.00
1.398
1.394
1.188
44
14.00
1.854
1.821
1.956
45*
2.000
2.699
2.154
2.221
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
each molecule atom. Five similarity indices (steric, electrostatic, hydrophobic, hydrogen bond donor and acceptor fields) were calculated using a sp3 hybridized carbon atom with a þ1.0 charge, þhydrophobicity, þH-bond donor and þH-bond acceptor properties. Both the CoMFA and the CoMSIA use Gasteiger-Huckel charges. Since a large amount of molecular fields properties are regarded as independent variables, the partial least squares (PLS) [21,22]
17
method was adopted to treat the linear correlation between the fields and the bioactivities. Cross-validation analysis was performed with the leave-one-out (LOO) methodology to determine the optimal number of principal components (N) and the predictive 2 power of the model (qcv ). Furthermore, the non-cross-validation analysis was executed based on the N to evaluate the fitting ability of the 3D-QSAR model. Finally, the test set was imported, and its bioactivity was predicted by the 3D-QSAR model to assess the
18
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
Fig. 2. (A) Alignment based on the receptor. (B) Alignment based on the template molecule.
robustness of the derived models.
3. Results and discussion 3.1. 3D-QSAR statistical results The ligand-based alignment rule and the docking-based alignment rule were applied to build the 3D-QSAR models [23]. The basic requirements for the model are as follows: the cross valida2 tion coefficient (qcv ) must be greater than 0.5, the standard error of estimate (SEE) must be close to zero and R2ncv must be greater than 0.9. The results demonstrated that the latter alignment rule was obviously better than the former in the CoMFA and CoMSIA models. Therefore, the docking-based alignment rule was chosen to establish reliable 3D-QSAR models. In the CoMFA, the obtained independent variables of the training set were subjected to cross2 validated PLS analysis to identify the value of qcv of 0.611 for 9 components, and then a non-cross-validated PLS analysis was performed to obtain SEE ¼ 0.063, R2ncv ¼ 0.999 and F values ¼ 4913.304. For the CoMFA model, only steric and electrostatic field contributions were calculated, and their values were 0.409 and 0.591, respectively. The contribution of the electrostatic field was greater than that of the steric field for increasing the bioactivity. The possible CoMFA results are shown in Table 2. For the CoMSIA model, there are five different descriptor fields that are related to each other [24,25]. Therefore, all possible combinations of fields were calculated to determine the best predictive model. First, the effect of each single field on CoMSIA model was investigated. The result showed that the steric field had the highest 2 qcv value. Second, two fields were considered and one of the fields must be steric field. The predictive power of the CoMSIA model was increased when the steric and electrostatic fields were taken into account simultaneously. Based on the results of the second step,
Table 2 Possible CoMFA results. Descriptors
q2cva
Nb
SEEc
R2ncvd
Fe
Steric (S) Electrostatic (E) SþE
0.572 0.551 0.611
9 9 9
0.057 0.131 0.063
0.999 0.997 0.999
5961.632 1146.260 4913.304
a c b d e
Cross-validated correlation coefficient using the leave-one-out methods. Standard error of estimate. Optimum number of components. Non-cross-validated correlation coefficient. F-test value.
other three fields were coupled with the steric and electrostatic 2 fields to determine the best qcv values. Following a series of cal2 culations, a satisfactory result was obtained with qcv ¼ 0.588, R2ncv ¼ 0.993, SEE ¼ 0.198 and F ¼ 638.967. All of the calculated results are listed in Table 3. Four fields were included in the best CoMSIA model. The contributions of the steric, electrostatic, hydrophobic, hydrogen-bond acceptor fields were 0.130 0.398, 0.298 and 0.174, respectively. Electrostatic and hydrophobic descriptor fields played the dominant roles. The final results of the CoMFA model and the CoMSIA model are listed in Table 4. The bioactivities of the test and training sets were predicted by the CoMFA and CoMSIA models to evaluate the predictive power of the 3D-QSAR model. The correlation between the actual activity and the predicted values is displayed in Fig. 3. These results demonstrated that the biological activity forecast by the 3D-QSAR model is in good agreement with the experimental data. Therefore, the constructed CoMFA and CoMSIA models are robust and can be used to design new compound and predict the activities of new compounds in the future. 3.2. Analysis of the CoMFA and CoMSIA contour maps Analysis of the visual contour maps derived from the CoMFA and CoMSIA models is favourable for deeply understanding the correlation between the structure and bioactivity, and is helpful for identifying the special area where modifications will significantly affect the affinity of the ligand binding to the enzyme. Therefore,
Table 3 Possible CoMSIA results. Descriptors
q2a cv
Nb
SEEc
R2d ncv
Fe
Steric (S) Electrostatic (E) Hydrophobic (H) H-bond donor (D) H-bond acceptor (A) SþE SþH SþD SþA SþEþA SþEþH SþEþD SþEþAþH SþEþAþD SþEþHþD EþSþAþDþH
0.547 0.468 0.299 0.148 0.437 0.555 0.402 0.369 0.475 0.573 0.587 0.500 0.588 0.548 0.503 0.523
7 2 3 1 7 6 3 6 7 5 5 7 7 9 8 9
0.261 0.988 1.033 1.725 0.527 0.328 0.897 0.489 0.339 0.338 0.319 0.263 0.198 0.200 0.170 0.131
0.988 0.791 0.777 0.345 0.948 0.979 0.832 0.954 0.979 0.977 0.980 0.987 0.993 0.993 0.995 0.997
286.377 71.724 42.995 20.552 86.462 269.321 60.998 117.891 215.894 302.566 341.403 362.455 638.967 488.235 766.427 1142.701
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
19
deviation and coefficient (StDev*Coeff) is chosen to view the output of the CoMFA and CoMSIA analysis.
Table 4 Calculated data for the 3D-QSAR model. Data
CoMFA
CoMSIA
q2a cv R2d ncv SEEc e F R2pre Nb Field contribution Steric Electrostatic Hydrophobic H-bond donor H-bond acceptor
0.611 0.999 0.063 4913.304 0.833 9
0.588 0.993 0.198 638.967 0.775 8
0.409 0.591
0.130 0.398 0.298 0.174
the contour maps around compound 10 that has the highest activity value are presented in pictures. The field type of the standard
3.2.1. Electrostatic contour map The electrostatic contours of the CoMFA and CoMSIA models are shown in Fig. 4A and B respectively. The contour maps of the electrostatic field are shown in the red and blue-coloured region. The positively charged groups are present in the blue region and the negatively charged groups appear in the red region, contributing to improving bioactivity. The red region appears in the nitrogen and sulfur heteroaromatic ring of the R1 substituents, indicating that the negatively charged groups are favourable for improving the bioactivity. Red and blue contours are observed around the R2 substituents, so that the groups should be carefully selected at this position. The spatial orientation of compounds 4, 9, 10, 15, 16, 22, 24, 25, 29, 30, 32 and 46 are almost similar where the electronegative R1 substitutions are located in most red contours
Fig. 3. Scatter plot of the predicted values versus the actual pIC50 values based on the CoMFA model (A) and the CoMSIA model (B).
Fig. 4. (A) and (B) CoMFA and CoMSIA electrostatic contour maps around compound 10, respectively: favored (blue) and disfavoured disfavored (red); (C) and (D) Steric contour maps of CoMFA and CoMSIA respectively: favored (green) and disfavoured (yellow).
20
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
and the electropositive R2 substitutions are found on the blue contours. The activities of compounds 4, 9, 10, 15, 16, 22 and 24 are much higher than those of compounds 25, 30 and 32, because the (S)-2-methoxycarbonylamino-3-phenyl-1-propanoyl group with comprehensive properties is more beneficial for improving the bioactivity than the only electron-rich 2,2-dimethyl-1-propanoyl group at R2. Compound 24 bearing a strong electron-rich (4-(thiophen-2-yl)thiazol)-2-yl group at R1 shows higher activity than compound 22 with a weakly electron-rich (4-ethylthiazol)-2-yl group. Compounds 39, 41, 43, 44, 47 and 49 with medium bioactivities reveal a similar spatial orientation in which R1 and R2 substitutions occupy the contours opposite to the highly potent compounds 4, 9, 10, 15, 16, 22 and 24. Therefore, the activity of compound 46 is much better than that of 47. Compounds 26, 27, 28, 29, 30, 31, 32, 34, 35, 36, 37 and 38 with lower bioactivities show various spatial orientations. For the same spatial direction, compound 27 is more potent than compound 28, because the (4hydroxymethyl)thiazol)-2-yl substituent is more electronegative than the (4-ethoxycarbonyl)thiazol)-2-yl group. In addition, the activities of compounds 34 and 31 are superior to those of compounds 35 and 38, respectively, due to the different electronegativity of their corresponding R1 substitutions. According to the preceding analysis, the presence of electropositive and electronegative groups in R2 substituents is favourable for increased bioactivity. This result confirms the conclusion that the dramatic difference in the activities between compounds 15 and 53 is due to the presence of electron-rich only 5-phenyl-1,3,4thiadiazol-2-yl substituent at R2 in compounds 53. The electronrich diaryl groups are suitable for R1 substituents, as observed from compounds 10 and 18. The (4-(thiophen-2-yl)thiazol)-2-yl substituent is much more electronegative than the (4ethylthiazol)-2-yl group. This result is consistent with the docking studies. The presence of Tyr1733, Asn1735 and His1945 amino residues leads to the p-p stacking interaction and electrostatic interaction with the diaryl groups.
3.2.2. Steric contour map The steric contour maps of the CoMFA and CoMSIA models are shown in Fig. 4C and D in yellow and green colours. The green contour map denotes the regions where a bulky substituent would improve the activity, and the yellow contour map indicates region where steric bulk groups would reduce the activity. It is clear that the R1 and R2 substitutions are surrounded by most of the yellow areas. The phenomenon demonstrates that bulky groups are unfavourable for increasing the activity. However, small green contours are present at the terminal of the R1 and R2 substitutions. Since the yellow contours are found next to the green contours at the terminal of the R1 substitutions, the steric bulk groups must be selected judiciously. The steric bulk groups may be randomly chosen at the terminal of R2 substitutions, because the terminal of R2 substituent extends to the outside of the active pocket as observed from the molecular docking results. Because the 4-(3methoxyphenyl)thiazol-2-yl group is more bulky than the phenylthiazol-2-yl substituent at R1, compound 25 is more effective than compound 30. In the terminal of the R1 substitutions, a barrier is present and comprised by the Arg1809, Val 1810 and Lys1811 amino residues, and therefore, binding steric bulk groups are not conducive to improving the activity. The activity of compound 46 exceeds that of compound 40 due to the R2 substitution 3-(3chlorophenyl)propanoyl reaching to the green contours. Following this reasoning, compounds 18, 20 and 26 show the following binding affinity order: compound 18 ((S)-2-acetamido-3phenylpropaoyl) > compound 20 ((2-methoxycarbonylamino)-
acetonyl) > compound 26 (2,2-dimethyl-propanonyl). 3.2.3. Hydrophobic contour map Fig. 5A depicts the hydrophobic contour maps of the CoMSIA model. The yellow contours indicate that the hydrophobic group is related to increasing bioactivity, whereas, white contours indicate that binding hydrophobic groups would contribute to the opposite effect on the bioactivity. An examination of Fig. 5B reveals that some yellow contours overlap on the red region of the electrostatic contour map at the location of the R1 substituent. This phenomenon may account for the greater activity of compound 4 compared to that of compound 5. This finding is because the hydrophobicity of 4-cyclohexylthiazol-2-yl group is stronger than that of the 4ethyl-5-methylthiazol-2-yl group. Compound 5 is also more potent than compound 8 with the hydrophilic 4-(ethoxycarbonyl) thiazol-2-yl group. Some yellow contours are observed next to the red region of the electrostatic contour map near the methoxycarbonylamino of the R2 substituent, explaining the difference in the bioactivity between compounds 46 and 47 with the 3-(3chlorophenyl)propanoyl group and the 3-(3-methoxyphenyl) propanoyl substituent in R2, respectively. Other contours align with green contours of the steric contour map at the location of the end of R2 substituents. This discovery explains why the order of the R2 groups favourability is (S)-2-methoxycarbonylamino-3-phenyl-1propanoyl group > tert-butoxycarbonyl group > 2,2-dimethylpropanonyl group for the same R1 group. The white area of the hydrophilic contour map overlaps with the blue region of the electrostatic contour map. This result could be used to explain why the activity of compound 54 with the electropositive hydrophilic 5((methoxycarbonyl)methyl)-1,3,4thiadiazol-2-yl group is better than that of compound 53 with the electronegative hydrophobic 5phenyl-1,3,4-thiadiazol-2-yl group. 3.2.4. H-bond acceptor contour map For the H-bond acceptor contour maps (Fig. 5C), the magenta region indicates that a hydrogen bond acceptor is favored for improving the activity. The red region signifies that a hydrogen bond acceptor produces the opposite effect. Both of those regions are mainly distributed around the R2 substituents except for a small fraction of the magenta contours located in the R1 substituents. Based on the docking results, the presence of Tyr1733, Asn1734, Asn1735, Lys1811, His 1871, His1945 and Gln 1948 amino acids provides H-bond donor/acceptor for the compounds to generate hydrogen bonds to increase their biological activity. For example, the H-bond is formed between the trifluoromethyl of the R1 group of compound 6 and Lys1811, contributing to its bioactivity being greater than that of compound 3 with the t-butyl group at R1. Compounds 44 and 49 show different activities because of the different numbers of formed hydrogen bonds. In compound 44, the R1 and R2 substituents form three H-bonds with the Asn1735 and Arg1809 amino acids. However, the R1 and R2 substituents form five H-bonds with the Asn1734 and Asn1735 amino acids. Hence, it is necessary to introduce substituents containing H-bond donor or Hbond acceptor in these locations in order to increase the bioactivity. According to the above analysis of the CoMFA and CoMSIA contour maps, the deduced optimization strategy for improving the inhibitory activities of the compounds against HPTPb is summarized and presented in Fig. 6. 1) R1: Electronegative, H-bond acceptor, and hydrophobic groups are favourable, and the bulky group is unfavourable. 2) R2: Electronegative, electropositive, hydrophilic, hydrophobic and H-bond donor/acceptor groups are beneficial for increasing
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
21
Fig. 5. (A) Hydrophobic contour maps of CoMSIA: favored (yellow) and disfavoured (white). (B) Overlap of the hydrophobic and electrostatic contour maps of CoMSIA. (C) H-bond acceptor contour maps of CoMSIA: favored (magenta) and disfavoured (red).
Fig. 6. Structural modification tips derived from 3D-QSAR studies.
bioactivity, and the binding of the steric bulk group is suitable in the appropriate location. 4. Conclusion In this paper, a series of phenylsulfamic acid derivatives acting as HPTPb inhibitors are collected, optimized and calculated to construct the 3D-QSAR model. The generated CoMFA and CoMSIA models possess good internal and external predictive power (q2cv: 0.611, respectively, 0.588 and R2pred: 0.833, 0.775, respectively). The contour maps of the CoMFA and CoMSIA models provide insight into the correlation between the different fields and bioactivities, and provide a reference for the design of new compounds. The R1
substituents have the electronegativity, hydrophobicity and Hbond acceptor characteristics that are favourable for increasing the compound bioactivity. The R2 substituents should have comprehensive characteristics or maintain and simply modify the structure of the R2 substituent of the molecular template. Analysis of the 3D-QSAR model is coupled with the results of molecular docking calculations, providing better understanding of the binding mode of the inhibitors to the active site of HPTPb and determining the amino acid residues at the periphery of the active site. Therefore, these models not only can be applied for predicting the activity of a new HPTPb inhibitor but also can provide guidance for the further design of inhibitors.
22
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
Conflicts of interest The authors declare that this article content has no conflict of interest. Acknowledgements The study was supported by the Major New Drugs Innovation and Development (grant numbers 2015ZX09102003). Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.molstruc.2019.02.107. References [1] N.X. Krueger, M. Streuli, H. Saito, Structural diversity and evolution of human receptor-like protein tyrosine phosphatases, EMBO J. 9 (1990) 3241e3252. [2] G. Fachinger, U. Deutsch, W. Risau, Functional interaction of vascular endothelial-protein-tyrosine phosphatase with the Angiopoietin receptor Tie2, Oncogene 18 (1999) 5948e5953. [3] H.T. Yuan, E.V. Khankin, S.A. Karumanchi, S.M. Parikh, Angiopoietin 2 is a partial Agonist/Antagonist of Tie2 signaling in the endothelium, Mol. Cell Biol. 29 (2009) 2011e2022. [4] T. Mammoto, S.M. Parikh, A. Mammoto, D. Gallagher, B. Chan, G. Mostoslavsky, D.E. Ingber, V.P. Sukhatme, Angiopoietin-1 Requires p190 RhoGAP to protect against vascular leakage in vivo, J. Biol. Chem. 282 (2007) 23910e23918. [5] M. Frye, M. Dierkes, V. Küppers, M. Vockel, J. Tomm, D. Zeuschner, J. Rossaint, A. Zarbock, Y.K. Gou, K. Peters, A.F. Nottebaum, D. Vestweber, Interfering with VE-PTP stabilizes endothelial junctions in vivo via Tie-2 in the absence of VEcadherin, J. Exp. Med. 212 (2015) 2267e2287. [6] S.R. Klopfenstein, A.G. Evdokimov, A.O. Colson, N.T. Fairweather, J.J. Neuman, M.B. Maier, J.L. Gray, G.S. Gerwe, G.E. Stake, B.W. Howard, J.A. Farmer, M.E. Pokross, T.R. Downs, B. Kasibhatla, K.G. Peters, 1,2,3,4-Tetrahydroisoquinolinyl sulfamic acids as phosphatase PTP1B inhibitors, Bioorg. Med. Chem. Lett. 16 (2006) 1574e1578. [7] S. Sankar, R.K. Raheja, M.J. Newman, A.S. Bhat, D.H. Slee, K.J. Lee, Y.N. Lee, S.J. Mcconnell, E.A. Coats, T. Nguyen, Sulfamic Acids as Inhibitors of Human Cytoplasmic Protein Tyrosin Phosphatases, Ontogen Corporation, 2003, p. 2003. WO03082263A1. [8] K.K.D. Amarasinghe, A.G. Evidokimov, K. Xu, C.M. Clark, M.B. Maier, A. Srivastava, A.O. Colson, G.S. Gerwe, G.E. Stake, B.W. Howard, M.E. Pokross, Design and synthesis of potent, non-peptidic inhibitors of HPTPb, Bioorg. Med. Chem. Lett 16 (2006) 4252e4256.
[9] J. Shen, M. Frye, B.L. Lee, J.L. Reinardy, J.M. McClung, K. Ding, M. Kojima, H.M. Xia, C. Seidel, R.L. Silva, A. Dong, S.F. Hackett, J.X. Wang, B.W. Howard, D. Vestweber, C.D. Kontos, K.G. Peters, P.A. Campochiaro, Targeting VE-PTP activates TIE2 and stabilizes the ocular vasculature, J. Clin. Investig. 124 (2014) 4564e4576. [10] R.D. Crarner, J.D. Bunce, D.E. Patterson, Crossvalidation, bootstrapping, and partial least squares compared with multiple regression in conventional QSAR Studies, Quant. Struct. Act. Relat. 7 (1988) 18e25. [11] G. Klebe, U. Abraham, T. Mietzner, Molecular similarity indices in a comparative analysis (CoMSIA) of drug molecules to correlate and predict their biological activity, J. Med. Chem. 37 (1994) 4130e4146. [12] R.D. Cramer, D.E. Patterson, J.D. Bunce, Comparative molecular field analysis (CoMFA). Part 1. effect of shape on binding of steroids to carrier proteins, J. Am. Chem. Soc. 110 (1988) 5959e5967. [13] R.R. Mittal, R.A. Mckinnon, M.J. Sorich, Effect of steric molecular field settings on CoMFA predictivity, J. Mol. Model. 14 (2008) 59e67. [14] F. Wang, B. Zhou, Toward the identification of a reliable 3D-QSAR model for the protein tyrosine phosphatase 1B inhibitors, J. Mol. Struct. 1158 (2018) 75e87. [15] K.G. Peters, R. Shalwitz, Compositions and Methods for Treating Ocualr Edema, Neovascularization and Related Diseasea, Aerpio Therapeutics INC, 2012. WO2012047966. [16] M. Clark, R.D. Cramer, N. Van Opdenbosch, Validation of the general purpose Tripos 5.2 Force field, J. Comput. Chem. 10 (1989) 982e1012. [17] J. Gasteiger, M. Marsili, Iterative partial equalization of orbital electronegativity-a rapid access to atomic charges, Tetrahedron 36 (1980) 3219e3228. [18] S.R. Wilson, W. Cui, F. Guarnieri, Chapter 16 Conformational analysis of flexible molecules, Data Handling Sci. Technol. 15 (1995) 351e367. [19] A.G. Evdokimov, M. Pokross, R. Walter, M. Mekel, B. Cox, C. Li, R. Bechard, F. Genbauffe, R. Andrews, C. Diven, B. Howard, V. Rastogi, J. Gray, M. Maier, K.G. Peters, Engineering the catalytic domain of human protein tyrosine phosphatase for structure-based drug discovery, Acta. Crystallogr. D. Biol. Crystallogr. 62 (2006) 1435e1445. [20] C.L. Waller, T.I. Oprea, A. Giolitti, G.R. Marshall, Three-dimensional QSAR of human immunodeficiency virus (I) protease inhibitors. 1. A CoMFA study employing experimentally-determined alignment rules, J. Med. Chem. 36 (1993) 4152e4160. [21] L. Stahle, S. Wold, 6 Multivariate data analysis and experimental design in biomedical research, Prog. Med. Chem. 25 (1988) 291e338. [22] B.L. Bush, R.B. Nachbar, Sample-distance partial least squares: PLS optimized for many variables, with application to CoMFA, J. Comput. Aided Mol. Des. 7 (1993) 587e619. [23] K.C. Tsai a, Y.C. Chen, N.W. Hsiao, C.L. Wang, C.L. Lin, Y.C. Lee, M.Y. Li, B.H. Wang, A comparison of different electrostatic potentials on prediction accuracy in CoMFA and CoMSIA studies, Eur. J. Med. Chem. 45 (2010) 1544e1551. [24] D.M. Hawkins, The problem of overfitting, J. Chem. Inf. Comput. Sci. 44 (2004) 1e12. [25] A. Golbraikh, A. Tropsha, Beware of q2!, J. Mol. Graph. Model. 20 (2002) 269e276.
Contents lists available at ScienceDirect
Journal of Molecular Structure journal homepage: http://www.elsevier.com/locate/molstruc
3D-QSAR study of the phenylsulfamic acid derivatives as HPTPb inhibitors Wenjuan Zhang a, b, Zhao Wei c, Chunyu Lin a, Zhonghua Wang d, Zhibing Zheng b, *, Song Li a, b, ** a
School of Pharmaceutical Engineering, Shenyang Pharmaceutical University, Shenyang, 110016, China National Engineering Research Center for the Emergency Drug, Beijing Institute of Pharmacology and Toxicology, Beijing 100850, China Department of Medicinal Chemistry, School of Pharmacy, Fourth Military Medical University, Xi'an, 300071, China d School of Chemical and Environmental Engineering, Shanghai Institute of Technology, Shanghai, 201203, China b c
a r t i c l e i n f o
a b s t r a c t
Article history: Received 16 January 2019 Received in revised form 25 February 2019 Accepted 25 February 2019 Available online 28 February 2019
Human protein tyrosine phosphatase beta (HPTPb) inhibitors have been used for the treatment of sepsis, cancer and inflammatory diseases. The phenylsulfamic acid derivatives are a group of powerful inhibitors of HPTPb. To explore the relationship between their structures and biological activities, the threedimensional quantitative structure-activity relationship (3D-QSAR) model is first constructed. Two highly predictive 3D-QSAR models are established. These models are the comparative molecular field 2 0.611, R2ncv 0.999) model and the comparative molecular similarity index analysis analysis (CoMFA: qcv 2 0.588, R2ncv 0.993) model. The results show a quite good external predictive power for the (CoMSIA: qcv test set, with R2pre values of 0.833 and 0.775, respectively. Furthermore, the contour maps of the 3D-QSAR models are analysed together with the results of molecular docking. The analysed results are conducive to discovering new binding sites and provide a reference for the construction of new potent HPTPb inhibitors compounds. © 2019 Elsevier B.V. All rights reserved.
Keywords: 3D-QSAR HPTPb Phenylsulfamic acid Molecular docking
1. Introduction Human protein tyrosine phosphatase beta (HPTPb) was discovered by Krueger in 1990 [1]. HPTPb is a classical receptor-like protein tyrosine phosphatase (PTP) and is specifically expressed in endothelial cells [2]. The mouse orthologue of HPTPb is called vascular endothelial protein tyrosine phosphatase (VE-PTP). HPTPb is involved in the angiopoietin (ANG)etyrosine kinase receptor (TIE) pathway. Angiopoietin 1 (ANG1) is a potent TIE2 activator whereas Angiopoietin 2 (ANG2) is a context-dependent agonist/ antagonist [3]. One of the physiological functions of activated TIE2 is to maintain the vascular endothelial integrity and reduce the leakage of blood vessels [4]. When HPTPb is inhibited, the rate of dephosphorylation of TIE2 is decreased, which is equivalent to either triggering TIE2 or to the result of the amplification of the ANG1 signal. Therefore, HPTPb inhibitors can be used for treating
* Corresponding author. ** Corresponding author. School of Pharmaceutical Engineering, Shenyang Pharmaceutical University, Shenyang, 110016, China. E-mail addresses: [email protected] (Z. Zheng), [email protected] (S. Li). https://doi.org/10.1016/j.molstruc.2019.02.107 0022-2860/© 2019 Elsevier B.V. All rights reserved.
sepsis, organ transplantation, atherosclerosis and vascular complications of diabetes by the activation of TIE2 [5]. Peters et al. identified many sulfamic acid derivatives such as compound 1 by using high-throughput screening to screen the P&GP corporate repository against multiple phosphatase enzymes. However, these compounds displayed weakly inhibitory effect on PTP1B, and the value of the best affect was 322.5 mM against PTP1B [6]. Upon structural optimization, the sulfamic acid moiety [7] was directly linked to an aryl group to mimic the PTP substratephosphotyrosine (pTyr). The inhibitory effect of phenylalanine analogue 2 on PTP1B (39.6 mM) was improved. In subsequent research, these phenylsulfamic acid compounds were also found to exhibit great inhibition activity for HPTPb. Subsequently, they regarded HPTPb as a target that previously attracted little research attention. Therefore, a series of compounds such as 3 with good activity against HPTPb were designed and synthesized, [8] and a potent, non-peptidic phenylsulfamic acid derivative was identified as a clinical candidate for treating diabetic macular edema (DME) [9] (see Fig. 1). 3D-QSAR is based on using the three-dimensional structure information and the biological activity of compounds to establish a reasonable mathematical model. The generated model is used to
12
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
Fig. 1. Development of phenylsulfamic acid derivatives.
determine the relationship between the structures and bioactivities. Furthermore, the model is also used to understand the mode of ligand interaction with the receptor. Comparative molecular field analysis (CoMFA) [10] and comparative molecular similarity indices analysis (CoMSIA) [11] are two common methods for building the models. In CoMFA and CoMSIA, the interaction between the receptor and the ligand is considered to be dependent on the differences in the molecular fields around the compounds. Therefore, the quantified molecular fields' parameters are regarded as independent variables, and the regression analysis of bioactivity is performed to represent the interaction between the ligand and the receptor. The CoMSIA is the further development of the CoMFA. The Lennard-Jones and Coulomb potential functions of the CoMFA are replaced with a distance-related Gaussian function. The function can effectively avoid the steepness of the potentials close to the molecular surface in the CoMFA [12,13]. The steric field, the electrostatic field, the hydrophobic field and the H-bond donor/ acceptor field are considered in the CoMSIA. Although a more stable 3D-QSAR model can be obtained by the CoMSIA, the results may not be better than the CoMFA. To accurately understand and use the 3D-QSAR model, two models were derived to complement and verify each other [14]. Here, 54 phenylsulfamic acid derivatives with a range of HPTPb inhibition activities were selected from the literature [15] as a set. The set was used to establish a reliable 3D-QSAR model by CoMFA and CoMSIA methods carried out using the SYBYL-X2.0 software (Tripos Associates, St. Louis, MO).
2. Methods and materials 2.1. The set and biological activity The selected 54 phenylsulfamic acids in the set had diverse structural features and showed a wide range of HPTPb inhibitory activities. Forty-one compounds were assigned to the training set for the model generation, and 13 compounds were selected as the test set for the model validation. The IC50 values of the HPTPb inhibitor were converted into the corresponding pIC50 (-logIC50) values. The structures and their actual or predicted pIC50 values of the compounds are listed in Table 1.
2.2.2. Acquiring the active configurations of compounds Obtaining accurate active conformations of the molecules corresponding to their bioactivities is a key step in the process of the 3D-QSAR model construction. Molecular docking is used to imitate and analyze the mode of the interaction between the ligand and the receptor. Using the docking results, the efficient active conformations of molecules can be determined. Therefore, the lowest energy conformations of the compounds obtained in the last step were collected and docked with the HPTPb protein one by one. The other parameters were set to default values in the SYBYL-X 2.0 program. The maximum number of poses per conformation was 20. Among the large number of docking poses, the pose with the highest C_Score value was chosen as the active conformation. When the highest C_Score values were the same, the pose with the highest Total_Score value was selected as the active conformation. All of the active conformations were collected and saved as the training set. The protein structure was obtained from the Protein Data Bank (PDB ID: 2I4H) [19]. The protein was treated by deleting all water molecules, extracting the ligand, adding hydrogen atoms to the protein and repairing the skeleton of the protein prior to its use. The ligand of the 2I4H complex underwent the same process to obtain the global lowest energy conformation. The acquired conformation of the ligand was re-docked into the HPTPb rigid protein by SurflexDock to validate the reliability of this docking method. As shown in Fig. S1 (supporting information), the retrieving pose of the ligand could overlay with the original confirmation and could interact properly with the receptor. 2.3. Compounds alignment Molecular alignment is crucial for the establishment of a 3DQSAR model [20]. Thus, two alignment methods were investigated in this study. The first alignment method was based on the receptor. The chosen docking poses were directly used for forming the 3DQSAR model. The result of this alignment is shown in Fig. 2A. The second alignment method was the template ligand-based alignment. The most potent compound was regarded as the template, and the 4-ethyl-phenylsulfamic acid was chosen as the common skeleton to be aligned by the “align database” function. The result is shown in Fig. 2B. 2.4. 3D-QSAR studies
2.2. Acquisition of molecular configurations 2.2.1. Global minimum energy configurations of compounds To explore the minimum energy conformations, all of the compounds were first optimized using the SYBYL-X 2.0 program by the Minimize with the Powell method, Tripos Force Field [16] and Gasteiger-Hückel charge [17] to obtain the local low energy conformation of the molecule. The calculation parameters are listed in Table S1 (supporting information). Then, these molecular conformations were continuously calculated by the simulated annealing [18] method with molecular dynamics to obtain the global minimum energy conformations. The calculation parameters are shown in Table S2.
The CoMFA and the CoMSIA are currently the two most widely used 3D-QSAR methods and differ mainly in their use of different fields to characterize the non-covalent interactions between the ligand and the receptor. To obtain the CoMFA and CoMSIA descriptor fields, the compounds were located in a 3D cubic lattice with the grid spacing of 2 Å in the x, y and z directions. A sp3 carbon atom probe with a van der Waals radius of 1.52 Å and a charge of þ1.0 wanders in this region. In the CoMFA method, the LennardJones and Coul potential energy functions were used to calculate the steric field energy and electrostatic field energy. In the CoMSIA method, a distance-dependent Gaussian function was introduced and used to calculate the distance between the probe atom and
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
13
Table 1 Compound structure and their actual and predicted pIC50 values.
IC50a(nM)
pIC50b
pIC50c
pIC50d
1*
6.000
2.222
2.526
3.205
2
1.000
3.000
2.997
3.233
3
0.2000
3.699
3.684
3.627
4
0.00005
7.301
7.250
7.144
5
0.1000
4.000
4.168
4.162
6
0.0800
4.097
4.076
4.222
7
1.000
3.000
2.996
3.044
8
0.2000
3.690
3.758
3.631
9
0.0020
5.699
5.518
5.235
No.
R1
R2
14
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
10
0.00005
7.301
7.365
7.423
11
1.000
3.000
2.945
3.276
12*
0.4000
3.398
3.300
3.005
13*
0.2000
3.699
4.415
5.230
14*
0.0800
4.097
4.209
5.984
15
0.0007
6.155
6.143
6.146
16
0.00005
7.301
7.342
7.207
17*
0.00005
7.301
6.399
7.881
18
1.000
3.000
2.963
2.949
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
15
19*
0.3000
3.523
4.032
3.772
20
20.00
1.699
1.731
1.699
21*
3.000
2.523
2.029
2.583
22
0.3000
3.523
3.532
3.543
23
0.3000
3.523
3.528
3.472
24
0.00005
7.301
7.331
7.463
25
49.00
1.310
1.281
1.349
26
85.00
1.071
1.068
0.894
27
266.0
0.575
0.533
0.719
28
584.0
0.234
0.233
0.274
29
113.0
0.947
0.944
0.703
30
110.0
0.959
0.957
0.900
31
138.0
0.860
0.862
0.737
16
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
32
98.00
1.009
1.091
1.440
33*
27.00
1.569
0.809
2.278
34
180.0
0.745
0.644
0.515
35
644.0
0.191
0.300
0.486
36
132.0
0.879
0.915
0.822
37
555.0
0.256
0.202
0.206
38
253.0
0.597
0.608
0.335
39
12.00
1.921
1.936
2.059
40
28.00
1.553
1.587
1.614
41
75.00
1.125
1.097
1.233
42*
56.00
1.252
1.752
1.777
43
40.00
1.398
1.394
1.188
44
14.00
1.854
1.821
1.956
45*
2.000
2.699
2.154
2.221
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
each molecule atom. Five similarity indices (steric, electrostatic, hydrophobic, hydrogen bond donor and acceptor fields) were calculated using a sp3 hybridized carbon atom with a þ1.0 charge, þhydrophobicity, þH-bond donor and þH-bond acceptor properties. Both the CoMFA and the CoMSIA use Gasteiger-Huckel charges. Since a large amount of molecular fields properties are regarded as independent variables, the partial least squares (PLS) [21,22]
17
method was adopted to treat the linear correlation between the fields and the bioactivities. Cross-validation analysis was performed with the leave-one-out (LOO) methodology to determine the optimal number of principal components (N) and the predictive 2 power of the model (qcv ). Furthermore, the non-cross-validation analysis was executed based on the N to evaluate the fitting ability of the 3D-QSAR model. Finally, the test set was imported, and its bioactivity was predicted by the 3D-QSAR model to assess the
18
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
Fig. 2. (A) Alignment based on the receptor. (B) Alignment based on the template molecule.
robustness of the derived models.
3. Results and discussion 3.1. 3D-QSAR statistical results The ligand-based alignment rule and the docking-based alignment rule were applied to build the 3D-QSAR models [23]. The basic requirements for the model are as follows: the cross valida2 tion coefficient (qcv ) must be greater than 0.5, the standard error of estimate (SEE) must be close to zero and R2ncv must be greater than 0.9. The results demonstrated that the latter alignment rule was obviously better than the former in the CoMFA and CoMSIA models. Therefore, the docking-based alignment rule was chosen to establish reliable 3D-QSAR models. In the CoMFA, the obtained independent variables of the training set were subjected to cross2 validated PLS analysis to identify the value of qcv of 0.611 for 9 components, and then a non-cross-validated PLS analysis was performed to obtain SEE ¼ 0.063, R2ncv ¼ 0.999 and F values ¼ 4913.304. For the CoMFA model, only steric and electrostatic field contributions were calculated, and their values were 0.409 and 0.591, respectively. The contribution of the electrostatic field was greater than that of the steric field for increasing the bioactivity. The possible CoMFA results are shown in Table 2. For the CoMSIA model, there are five different descriptor fields that are related to each other [24,25]. Therefore, all possible combinations of fields were calculated to determine the best predictive model. First, the effect of each single field on CoMSIA model was investigated. The result showed that the steric field had the highest 2 qcv value. Second, two fields were considered and one of the fields must be steric field. The predictive power of the CoMSIA model was increased when the steric and electrostatic fields were taken into account simultaneously. Based on the results of the second step,
Table 2 Possible CoMFA results. Descriptors
q2cva
Nb
SEEc
R2ncvd
Fe
Steric (S) Electrostatic (E) SþE
0.572 0.551 0.611
9 9 9
0.057 0.131 0.063
0.999 0.997 0.999
5961.632 1146.260 4913.304
a c b d e
Cross-validated correlation coefficient using the leave-one-out methods. Standard error of estimate. Optimum number of components. Non-cross-validated correlation coefficient. F-test value.
other three fields were coupled with the steric and electrostatic 2 fields to determine the best qcv values. Following a series of cal2 culations, a satisfactory result was obtained with qcv ¼ 0.588, R2ncv ¼ 0.993, SEE ¼ 0.198 and F ¼ 638.967. All of the calculated results are listed in Table 3. Four fields were included in the best CoMSIA model. The contributions of the steric, electrostatic, hydrophobic, hydrogen-bond acceptor fields were 0.130 0.398, 0.298 and 0.174, respectively. Electrostatic and hydrophobic descriptor fields played the dominant roles. The final results of the CoMFA model and the CoMSIA model are listed in Table 4. The bioactivities of the test and training sets were predicted by the CoMFA and CoMSIA models to evaluate the predictive power of the 3D-QSAR model. The correlation between the actual activity and the predicted values is displayed in Fig. 3. These results demonstrated that the biological activity forecast by the 3D-QSAR model is in good agreement with the experimental data. Therefore, the constructed CoMFA and CoMSIA models are robust and can be used to design new compound and predict the activities of new compounds in the future. 3.2. Analysis of the CoMFA and CoMSIA contour maps Analysis of the visual contour maps derived from the CoMFA and CoMSIA models is favourable for deeply understanding the correlation between the structure and bioactivity, and is helpful for identifying the special area where modifications will significantly affect the affinity of the ligand binding to the enzyme. Therefore,
Table 3 Possible CoMSIA results. Descriptors
q2a cv
Nb
SEEc
R2d ncv
Fe
Steric (S) Electrostatic (E) Hydrophobic (H) H-bond donor (D) H-bond acceptor (A) SþE SþH SþD SþA SþEþA SþEþH SþEþD SþEþAþH SþEþAþD SþEþHþD EþSþAþDþH
0.547 0.468 0.299 0.148 0.437 0.555 0.402 0.369 0.475 0.573 0.587 0.500 0.588 0.548 0.503 0.523
7 2 3 1 7 6 3 6 7 5 5 7 7 9 8 9
0.261 0.988 1.033 1.725 0.527 0.328 0.897 0.489 0.339 0.338 0.319 0.263 0.198 0.200 0.170 0.131
0.988 0.791 0.777 0.345 0.948 0.979 0.832 0.954 0.979 0.977 0.980 0.987 0.993 0.993 0.995 0.997
286.377 71.724 42.995 20.552 86.462 269.321 60.998 117.891 215.894 302.566 341.403 362.455 638.967 488.235 766.427 1142.701
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
19
deviation and coefficient (StDev*Coeff) is chosen to view the output of the CoMFA and CoMSIA analysis.
Table 4 Calculated data for the 3D-QSAR model. Data
CoMFA
CoMSIA
q2a cv R2d ncv SEEc e F R2pre Nb Field contribution Steric Electrostatic Hydrophobic H-bond donor H-bond acceptor
0.611 0.999 0.063 4913.304 0.833 9
0.588 0.993 0.198 638.967 0.775 8
0.409 0.591
0.130 0.398 0.298 0.174
the contour maps around compound 10 that has the highest activity value are presented in pictures. The field type of the standard
3.2.1. Electrostatic contour map The electrostatic contours of the CoMFA and CoMSIA models are shown in Fig. 4A and B respectively. The contour maps of the electrostatic field are shown in the red and blue-coloured region. The positively charged groups are present in the blue region and the negatively charged groups appear in the red region, contributing to improving bioactivity. The red region appears in the nitrogen and sulfur heteroaromatic ring of the R1 substituents, indicating that the negatively charged groups are favourable for improving the bioactivity. Red and blue contours are observed around the R2 substituents, so that the groups should be carefully selected at this position. The spatial orientation of compounds 4, 9, 10, 15, 16, 22, 24, 25, 29, 30, 32 and 46 are almost similar where the electronegative R1 substitutions are located in most red contours
Fig. 3. Scatter plot of the predicted values versus the actual pIC50 values based on the CoMFA model (A) and the CoMSIA model (B).
Fig. 4. (A) and (B) CoMFA and CoMSIA electrostatic contour maps around compound 10, respectively: favored (blue) and disfavoured disfavored (red); (C) and (D) Steric contour maps of CoMFA and CoMSIA respectively: favored (green) and disfavoured (yellow).
20
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
and the electropositive R2 substitutions are found on the blue contours. The activities of compounds 4, 9, 10, 15, 16, 22 and 24 are much higher than those of compounds 25, 30 and 32, because the (S)-2-methoxycarbonylamino-3-phenyl-1-propanoyl group with comprehensive properties is more beneficial for improving the bioactivity than the only electron-rich 2,2-dimethyl-1-propanoyl group at R2. Compound 24 bearing a strong electron-rich (4-(thiophen-2-yl)thiazol)-2-yl group at R1 shows higher activity than compound 22 with a weakly electron-rich (4-ethylthiazol)-2-yl group. Compounds 39, 41, 43, 44, 47 and 49 with medium bioactivities reveal a similar spatial orientation in which R1 and R2 substitutions occupy the contours opposite to the highly potent compounds 4, 9, 10, 15, 16, 22 and 24. Therefore, the activity of compound 46 is much better than that of 47. Compounds 26, 27, 28, 29, 30, 31, 32, 34, 35, 36, 37 and 38 with lower bioactivities show various spatial orientations. For the same spatial direction, compound 27 is more potent than compound 28, because the (4hydroxymethyl)thiazol)-2-yl substituent is more electronegative than the (4-ethoxycarbonyl)thiazol)-2-yl group. In addition, the activities of compounds 34 and 31 are superior to those of compounds 35 and 38, respectively, due to the different electronegativity of their corresponding R1 substitutions. According to the preceding analysis, the presence of electropositive and electronegative groups in R2 substituents is favourable for increased bioactivity. This result confirms the conclusion that the dramatic difference in the activities between compounds 15 and 53 is due to the presence of electron-rich only 5-phenyl-1,3,4thiadiazol-2-yl substituent at R2 in compounds 53. The electronrich diaryl groups are suitable for R1 substituents, as observed from compounds 10 and 18. The (4-(thiophen-2-yl)thiazol)-2-yl substituent is much more electronegative than the (4ethylthiazol)-2-yl group. This result is consistent with the docking studies. The presence of Tyr1733, Asn1735 and His1945 amino residues leads to the p-p stacking interaction and electrostatic interaction with the diaryl groups.
3.2.2. Steric contour map The steric contour maps of the CoMFA and CoMSIA models are shown in Fig. 4C and D in yellow and green colours. The green contour map denotes the regions where a bulky substituent would improve the activity, and the yellow contour map indicates region where steric bulk groups would reduce the activity. It is clear that the R1 and R2 substitutions are surrounded by most of the yellow areas. The phenomenon demonstrates that bulky groups are unfavourable for increasing the activity. However, small green contours are present at the terminal of the R1 and R2 substitutions. Since the yellow contours are found next to the green contours at the terminal of the R1 substitutions, the steric bulk groups must be selected judiciously. The steric bulk groups may be randomly chosen at the terminal of R2 substitutions, because the terminal of R2 substituent extends to the outside of the active pocket as observed from the molecular docking results. Because the 4-(3methoxyphenyl)thiazol-2-yl group is more bulky than the phenylthiazol-2-yl substituent at R1, compound 25 is more effective than compound 30. In the terminal of the R1 substitutions, a barrier is present and comprised by the Arg1809, Val 1810 and Lys1811 amino residues, and therefore, binding steric bulk groups are not conducive to improving the activity. The activity of compound 46 exceeds that of compound 40 due to the R2 substitution 3-(3chlorophenyl)propanoyl reaching to the green contours. Following this reasoning, compounds 18, 20 and 26 show the following binding affinity order: compound 18 ((S)-2-acetamido-3phenylpropaoyl) > compound 20 ((2-methoxycarbonylamino)-
acetonyl) > compound 26 (2,2-dimethyl-propanonyl). 3.2.3. Hydrophobic contour map Fig. 5A depicts the hydrophobic contour maps of the CoMSIA model. The yellow contours indicate that the hydrophobic group is related to increasing bioactivity, whereas, white contours indicate that binding hydrophobic groups would contribute to the opposite effect on the bioactivity. An examination of Fig. 5B reveals that some yellow contours overlap on the red region of the electrostatic contour map at the location of the R1 substituent. This phenomenon may account for the greater activity of compound 4 compared to that of compound 5. This finding is because the hydrophobicity of 4-cyclohexylthiazol-2-yl group is stronger than that of the 4ethyl-5-methylthiazol-2-yl group. Compound 5 is also more potent than compound 8 with the hydrophilic 4-(ethoxycarbonyl) thiazol-2-yl group. Some yellow contours are observed next to the red region of the electrostatic contour map near the methoxycarbonylamino of the R2 substituent, explaining the difference in the bioactivity between compounds 46 and 47 with the 3-(3chlorophenyl)propanoyl group and the 3-(3-methoxyphenyl) propanoyl substituent in R2, respectively. Other contours align with green contours of the steric contour map at the location of the end of R2 substituents. This discovery explains why the order of the R2 groups favourability is (S)-2-methoxycarbonylamino-3-phenyl-1propanoyl group > tert-butoxycarbonyl group > 2,2-dimethylpropanonyl group for the same R1 group. The white area of the hydrophilic contour map overlaps with the blue region of the electrostatic contour map. This result could be used to explain why the activity of compound 54 with the electropositive hydrophilic 5((methoxycarbonyl)methyl)-1,3,4thiadiazol-2-yl group is better than that of compound 53 with the electronegative hydrophobic 5phenyl-1,3,4-thiadiazol-2-yl group. 3.2.4. H-bond acceptor contour map For the H-bond acceptor contour maps (Fig. 5C), the magenta region indicates that a hydrogen bond acceptor is favored for improving the activity. The red region signifies that a hydrogen bond acceptor produces the opposite effect. Both of those regions are mainly distributed around the R2 substituents except for a small fraction of the magenta contours located in the R1 substituents. Based on the docking results, the presence of Tyr1733, Asn1734, Asn1735, Lys1811, His 1871, His1945 and Gln 1948 amino acids provides H-bond donor/acceptor for the compounds to generate hydrogen bonds to increase their biological activity. For example, the H-bond is formed between the trifluoromethyl of the R1 group of compound 6 and Lys1811, contributing to its bioactivity being greater than that of compound 3 with the t-butyl group at R1. Compounds 44 and 49 show different activities because of the different numbers of formed hydrogen bonds. In compound 44, the R1 and R2 substituents form three H-bonds with the Asn1735 and Arg1809 amino acids. However, the R1 and R2 substituents form five H-bonds with the Asn1734 and Asn1735 amino acids. Hence, it is necessary to introduce substituents containing H-bond donor or Hbond acceptor in these locations in order to increase the bioactivity. According to the above analysis of the CoMFA and CoMSIA contour maps, the deduced optimization strategy for improving the inhibitory activities of the compounds against HPTPb is summarized and presented in Fig. 6. 1) R1: Electronegative, H-bond acceptor, and hydrophobic groups are favourable, and the bulky group is unfavourable. 2) R2: Electronegative, electropositive, hydrophilic, hydrophobic and H-bond donor/acceptor groups are beneficial for increasing
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
21
Fig. 5. (A) Hydrophobic contour maps of CoMSIA: favored (yellow) and disfavoured (white). (B) Overlap of the hydrophobic and electrostatic contour maps of CoMSIA. (C) H-bond acceptor contour maps of CoMSIA: favored (magenta) and disfavoured (red).
Fig. 6. Structural modification tips derived from 3D-QSAR studies.
bioactivity, and the binding of the steric bulk group is suitable in the appropriate location. 4. Conclusion In this paper, a series of phenylsulfamic acid derivatives acting as HPTPb inhibitors are collected, optimized and calculated to construct the 3D-QSAR model. The generated CoMFA and CoMSIA models possess good internal and external predictive power (q2cv: 0.611, respectively, 0.588 and R2pred: 0.833, 0.775, respectively). The contour maps of the CoMFA and CoMSIA models provide insight into the correlation between the different fields and bioactivities, and provide a reference for the design of new compounds. The R1
substituents have the electronegativity, hydrophobicity and Hbond acceptor characteristics that are favourable for increasing the compound bioactivity. The R2 substituents should have comprehensive characteristics or maintain and simply modify the structure of the R2 substituent of the molecular template. Analysis of the 3D-QSAR model is coupled with the results of molecular docking calculations, providing better understanding of the binding mode of the inhibitors to the active site of HPTPb and determining the amino acid residues at the periphery of the active site. Therefore, these models not only can be applied for predicting the activity of a new HPTPb inhibitor but also can provide guidance for the further design of inhibitors.
22
W. Zhang et al. / Journal of Molecular Structure 1186 (2019) 11e22
Conflicts of interest The authors declare that this article content has no conflict of interest. Acknowledgements The study was supported by the Major New Drugs Innovation and Development (grant numbers 2015ZX09102003). Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.molstruc.2019.02.107. References [1] N.X. Krueger, M. Streuli, H. Saito, Structural diversity and evolution of human receptor-like protein tyrosine phosphatases, EMBO J. 9 (1990) 3241e3252. [2] G. Fachinger, U. Deutsch, W. Risau, Functional interaction of vascular endothelial-protein-tyrosine phosphatase with the Angiopoietin receptor Tie2, Oncogene 18 (1999) 5948e5953. [3] H.T. Yuan, E.V. Khankin, S.A. Karumanchi, S.M. Parikh, Angiopoietin 2 is a partial Agonist/Antagonist of Tie2 signaling in the endothelium, Mol. Cell Biol. 29 (2009) 2011e2022. [4] T. Mammoto, S.M. Parikh, A. Mammoto, D. Gallagher, B. Chan, G. Mostoslavsky, D.E. Ingber, V.P. Sukhatme, Angiopoietin-1 Requires p190 RhoGAP to protect against vascular leakage in vivo, J. Biol. Chem. 282 (2007) 23910e23918. [5] M. Frye, M. Dierkes, V. Küppers, M. Vockel, J. Tomm, D. Zeuschner, J. Rossaint, A. Zarbock, Y.K. Gou, K. Peters, A.F. Nottebaum, D. Vestweber, Interfering with VE-PTP stabilizes endothelial junctions in vivo via Tie-2 in the absence of VEcadherin, J. Exp. Med. 212 (2015) 2267e2287. [6] S.R. Klopfenstein, A.G. Evdokimov, A.O. Colson, N.T. Fairweather, J.J. Neuman, M.B. Maier, J.L. Gray, G.S. Gerwe, G.E. Stake, B.W. Howard, J.A. Farmer, M.E. Pokross, T.R. Downs, B. Kasibhatla, K.G. Peters, 1,2,3,4-Tetrahydroisoquinolinyl sulfamic acids as phosphatase PTP1B inhibitors, Bioorg. Med. Chem. Lett. 16 (2006) 1574e1578. [7] S. Sankar, R.K. Raheja, M.J. Newman, A.S. Bhat, D.H. Slee, K.J. Lee, Y.N. Lee, S.J. Mcconnell, E.A. Coats, T. Nguyen, Sulfamic Acids as Inhibitors of Human Cytoplasmic Protein Tyrosin Phosphatases, Ontogen Corporation, 2003, p. 2003. WO03082263A1. [8] K.K.D. Amarasinghe, A.G. Evidokimov, K. Xu, C.M. Clark, M.B. Maier, A. Srivastava, A.O. Colson, G.S. Gerwe, G.E. Stake, B.W. Howard, M.E. Pokross, Design and synthesis of potent, non-peptidic inhibitors of HPTPb, Bioorg. Med. Chem. Lett 16 (2006) 4252e4256.
[9] J. Shen, M. Frye, B.L. Lee, J.L. Reinardy, J.M. McClung, K. Ding, M. Kojima, H.M. Xia, C. Seidel, R.L. Silva, A. Dong, S.F. Hackett, J.X. Wang, B.W. Howard, D. Vestweber, C.D. Kontos, K.G. Peters, P.A. Campochiaro, Targeting VE-PTP activates TIE2 and stabilizes the ocular vasculature, J. Clin. Investig. 124 (2014) 4564e4576. [10] R.D. Crarner, J.D. Bunce, D.E. Patterson, Crossvalidation, bootstrapping, and partial least squares compared with multiple regression in conventional QSAR Studies, Quant. Struct. Act. Relat. 7 (1988) 18e25. [11] G. Klebe, U. Abraham, T. Mietzner, Molecular similarity indices in a comparative analysis (CoMSIA) of drug molecules to correlate and predict their biological activity, J. Med. Chem. 37 (1994) 4130e4146. [12] R.D. Cramer, D.E. Patterson, J.D. Bunce, Comparative molecular field analysis (CoMFA). Part 1. effect of shape on binding of steroids to carrier proteins, J. Am. Chem. Soc. 110 (1988) 5959e5967. [13] R.R. Mittal, R.A. Mckinnon, M.J. Sorich, Effect of steric molecular field settings on CoMFA predictivity, J. Mol. Model. 14 (2008) 59e67. [14] F. Wang, B. Zhou, Toward the identification of a reliable 3D-QSAR model for the protein tyrosine phosphatase 1B inhibitors, J. Mol. Struct. 1158 (2018) 75e87. [15] K.G. Peters, R. Shalwitz, Compositions and Methods for Treating Ocualr Edema, Neovascularization and Related Diseasea, Aerpio Therapeutics INC, 2012. WO2012047966. [16] M. Clark, R.D. Cramer, N. Van Opdenbosch, Validation of the general purpose Tripos 5.2 Force field, J. Comput. Chem. 10 (1989) 982e1012. [17] J. Gasteiger, M. Marsili, Iterative partial equalization of orbital electronegativity-a rapid access to atomic charges, Tetrahedron 36 (1980) 3219e3228. [18] S.R. Wilson, W. Cui, F. Guarnieri, Chapter 16 Conformational analysis of flexible molecules, Data Handling Sci. Technol. 15 (1995) 351e367. [19] A.G. Evdokimov, M. Pokross, R. Walter, M. Mekel, B. Cox, C. Li, R. Bechard, F. Genbauffe, R. Andrews, C. Diven, B. Howard, V. Rastogi, J. Gray, M. Maier, K.G. Peters, Engineering the catalytic domain of human protein tyrosine phosphatase for structure-based drug discovery, Acta. Crystallogr. D. Biol. Crystallogr. 62 (2006) 1435e1445. [20] C.L. Waller, T.I. Oprea, A. Giolitti, G.R. Marshall, Three-dimensional QSAR of human immunodeficiency virus (I) protease inhibitors. 1. A CoMFA study employing experimentally-determined alignment rules, J. Med. Chem. 36 (1993) 4152e4160. [21] L. Stahle, S. Wold, 6 Multivariate data analysis and experimental design in biomedical research, Prog. Med. Chem. 25 (1988) 291e338. [22] B.L. Bush, R.B. Nachbar, Sample-distance partial least squares: PLS optimized for many variables, with application to CoMFA, J. Comput. Aided Mol. Des. 7 (1993) 587e619. [23] K.C. Tsai a, Y.C. Chen, N.W. Hsiao, C.L. Wang, C.L. Lin, Y.C. Lee, M.Y. Li, B.H. Wang, A comparison of different electrostatic potentials on prediction accuracy in CoMFA and CoMSIA studies, Eur. J. Med. Chem. 45 (2010) 1544e1551. [24] D.M. Hawkins, The problem of overfitting, J. Chem. Inf. Comput. Sci. 44 (2004) 1e12. [25] A. Golbraikh, A. Tropsha, Beware of q2!, J. Mol. Graph. Model. 20 (2002) 269e276.