Deep Learning and Random Forest Approach for Finding the Optimal

Subscriber access provided by TULANE UNIVERSITY

Computational Biochemistry

Deep Learning and Random Forest Approach for Finding Optimal TCM Formula for Treatment of Alzheimer’s Disease Hsin-Yi Chen, JianQiang Chen, Jun-Yan Li, Hung-Jin Huang, Xi Chen, Hao-Ying Zhang, and Calvin Yu-Chian Chen J. Chem. Inf. Model., Just Accepted Manuscript • DOI: 10.1021/acs.jcim.9b00041 • Publication Date (Web): 19 Mar 2019 Downloaded from http://pubs.acs.org on March 19, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

The Random Forest model and Deep Learning model were utilized to build AD predicted models to calculate the possibility of our TCM formula. The result showed Methyl 3-O-feruloylquinate contained in Phellodendron amurense and Cynanogenin A contained in Cynanchum atratum are capable of forming stable interactions with GSK3β.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Deep Learning and Random Forest Approach for

2

Finding Optimal TCM Formula for Treatment of

3

Alzheimer’s Disease

4 5

Hsin-Yi Chen†#, Jian-Qiang Chen†#, Jun-Yan Li†#, Hung-Jin Huang†#, Xi Chen†,

6

Hao-Ying Zhang†, Calvin Yu-Chian Chen†,‡,§*

7 8

†

School of Intelligent Systems Engineering, Sun Yat-sen University, Shenzhen, 510275, China

9

‡

Department of Medical Research, China Medical University Hospital, Taichung 40447,

10

Taiwan

11

§Department

12

Taiwan

of Bioinformatics and Medical Engineering, Asia University, Taichung 41354,

13 14

#

15

* Corresponding Authors

16

Calvin Yu-Chian Chen, Ph.D.

17

School of Intelligent Systems Engineering, Sun Yat-sen University, Guangzhou 510275, China.

18

TEL: 15626413023

19

E-mail: [email protected]

Equal contribution

20 1 ACS Paragon Plus Environment

Page 2 of 64

Page 3 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2

Abstract:

3

It has demonstrated that Glycogen synthase kinase-3β (GSK3β) is related with Alzheimer’s

4

disease (AD). Based on the world largest traditional Chinese medicine (TCM) database, network

5

pharmacology-based approach were utilized to investigate the TCM candidates that can dock well

6

in multiple targets. Support Vector Machine(SVM), multiple linear regression (MLR) methods

7

were utilized to obtain predicted models, Specially, The Deep Learning and The Random Forest

8

algorithm were adopted, We achieved R2 of 0.927 on training set and 0.862 on the test set in Deep

9

Learning, The latter’s R2 on the training set is 0.869 and on the test set is 0.890.Besides,

10

comparative molecular similarity indices analysis (CoMSIA) was performed to get a predicted

11

model. All the training models achieved good results on test set.100ns MD simulation evaluated

12

the stability of GSK3β protein-ligand complexes. Methyl 3-O-feruloylquinate and Cynanogenin

13

A both induced more compactness to GSK3β complex and stable condition among all simulation

14

times, the GSK3β complex also had no substantial fluctuation after simulation time of 5 ns. For

15

TCM molecules, we used trained models to calculate predicted bioactivity values.Thus, the

16

optimum TCM candidates were obtained by ranking the predicted values.The result showed

17

Methyl 3-O-feruloylquinate contained in Phellodendron amurense and Cynanogenin A contained

18

in Cynanchum atratum are capable of forming stable interactions with GSK3β.

19 20

Keywords: Glycogen synthase kinase-3β (GSK3β); TCM formula; Network pharmacology;

21

Deep Learning; Random Forest; MD simulation

22

2 ACS Paragon Plus Environment


Page 4 of 64

1

1. Introduction

2

Alzheimer’s disease (AD) is one of the most common neurodegenerative diseases of the central

3

nervous system (CNS) that affects people over the age of 65. 1 The feature of Alzheimer’s disease

4

(AD) is the abnormal accumulation and processing of mutant or damaged intra and extracellular

5

proteins.2While the detailed molecular mechanism underlying the onset and progression of AD

6

remains elusive, previous studies have correlated it with the accumulation of amyloid beta (Aβ),

7

phosphorylated tau, 3 mitochondrial dysfunction, and neurofibrillary tangles in brains. 4 Targeting

8

of β-amyloid (Aβ) has been the focus in the development of improved treatments for AD, but

9

these have not yet shown any clinical benefit.5

10

Regarding the Aβ pathway, Aβ peptides accumulation formed by sequential cleavage of β-

11

amyloid precursor protein (APP) induces dementia syndromes in patients. 6 The cleaving enzyme,

12

β-site amyloid precursor protein cleaving enzyme 1 (BACE1), cleaves APP at two β-sites to

13

generate β peptide. 7 Inhibition of GSK-3 may help in the treatment of various diseases, such as

14

Alzheimer's disease. 8 A known contributor to this pathway is the glycogen synthase kinase 3 beta

15

(GSK3β), a proline-directed serine/threonine protein kinase,9 mainly known for its function in

16

regulating glycogen metabolism by inactivating glycogen synthase through phosphorylation.

17

Aberrant high phosphorylation activity of GSK3β has been linked to elevated BACE1 expression

18

and Aβ production.10 Therefore, inhibiting the phosphorylation activity of GSK3β can be a viable

19

and promising strategy to develop novel drugs and therapeutics for AD. 11Using multi-layer neural

20

network

21

(RF) model can successfully identify natural proteins from bait proteins. 14

22

Over the past few years,

23

2011, and we have developed computer-aided drug design15 webservers such as iScreen16 and

24

iSMART17, which help the related institutions and researchers improve their efficiency.

25

In this effort, we sought to develop novel formula with natural active ingredients capable of

12

,Deep Learning model shows a satisfactory predictive ability. 13 The Random Forest

we established the TCM database (TCM Database@Taiwan) since


Page 5 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

inhibiting GSK3β and other key components/proteins in the pathogenic pathway of AD by

2

exploiting the TCM database.

3

drug strategy may not affect entire biological system19-21, a specific target still not fully to discover

4

effectively drugs for treating diseases in organisms. Recent researches of drug design have been

5

based on multi-target concept as strategy in the developments of effective drugs22, 23. Network

6

analysis of signaling pathway were wildly used to identify potential factors for many diseases. 24-

7

26

8

AD was explored by utilizing the STRING database and the associated analyzing tools and a

9

collection of promising target proteins along with GSK3β was pinned down. TCM database

10

virtual screening, bioactivity prediction (QSAR modeling) and drug-target analysis were applied

11

to identify lead compounds in TCM database. Molecular docking was ensued to find out the

12

possible interaction modes between the active ingredients (small molecules) and the protein

13

targets. Then the stability of protein-ligand complex was verified by molecular dynamics (MD)

14

simulations and the detailed interactions between the ligands and the target proteins were

15

delineated. In the end, we presented these verified compounds from the TCM database that

16

possesses potential bioactivities against AD. (Figure 1)

17

2. Material and methods

18

2.1 Pathway enrichment analysis and network construction

19

The biological pathways was constructed using STRING database (STRING v10.5, http://string-

20

db.org/).27 It provides functional association between two proteins, the protein–protein interaction

21

(PPI) information were retrieved from physical interaction databases and databases of curated

22

biological pathway knowledge, such as IntAct, MINT, BioGRID, DIP, KEGG, Reactome, GO,

23

etc. It can be queried from inside the Cytoscape software framework. The PPI network of AD was

24

generated and visualized using the Cytoscape STRING app with data associated with Homo

18

The trend of drug development has reveal that one target one

Specifically, the protein-protein interaction (PPI) network associated with the pathogenesis of



Page 6 of 64

1

sapiens. Top 10 proteins with confidence score of 0.8 was set to construct the PPI network and

2

the network was expanded by adding new nodes. Functional enrichment analysis was conducted

3

using the information from KEGG database and the results were visualized using Cytoscape in

4

order to infer the target-disease relationship.

5

2.2 Screening and molecular docking

6

The crystal structure of potential targets for Alzheimer’s disease such as GSK3β, Nicastrin and

7

APBB1 were obtained from PDB database (PDB accession code: 4ACC, 5FN5 and 3D8D

8

respectively).28, 29 The Prepare Protein module of Accelrys Discovery Studio 2.5.5.9350 (DS 2.5)

9

was employed to insert missing atoms, model missing loops, and remove crystal waters.

10

Protonation state of each protein was assigned at pH 7.4. Total number of 61,000 TCM

11

compounds obtained from the website of TCM Database@Taiwan

12

criteria and Lipinski's rule of five. A number of 18776 of the resulting compounds were then

13

prepared for docking by energy minimization using the smart minimizer algorithm before the

14

subsequent docking studies, in which Monte-Carlo sampling with the CHARMm force field

15

(Chemistry at HARvard Molecular Mechanics) was chosen for ligand conformation generation,

16

using the LigandFit module of DS 2.5. 30 The active site of GSK3β ATP binding site was defined

17

based on the interaction of 7YG binding pose (the compound 23 in Berg's study,

18

(PDB code: 4ACC), the 7YG is selective for GSK-3β ATP binding site. During the docking

19

process, DREIDING force field, containing Gasteiger charges was chosen to calculate interaction

20

energy for each ligand.

21

2.3 Construction of QSAR model


18

were filtered by ADMET

28

to GSK3β

Page 7 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

2.3.1 2D-QSAR model

2

30 GSK3β inhibitors including 24 training molecules and 6 test molecules with bioactivity (Ki

3

value) from Berg's study 28 were utilized to build QSAR models (Table 1). All Ki values (nM)

4

were converted to pKi values using equation (1):

5

pKi = 9 - log(Ki)

(1)

6

The 2D structure of each molecule was created using ChemBioOffice 2010 and the molecular

7

properties of GSK3β inhibitors, the independent property, were calculated with the Calculate

8

Molecular Properties module of DS 2.5. The pKi of each ligand was treated as a dependent

9

property, as opposed to the calculated properties from chemical scaffold as independent ones in

10

the regression step. The Genetic Function Approximation (GFA) algorithm was employed to

11

select suitable properties related to activity value. Then the 2D-QSAR models were used to

12

predict the bioactivity of TCM candidates with the selected properties and bioactivity of each

13

GSK3β inhibitors.

14

2.3.1.1 Deep Learning model

15

.Aiming to enhance the accuracy of predicted activity value, the Deep Learning (DL) method

16

was applied to build a new model. We construct a simple 4-layers full connected neural network

17

using ReLu (Rectified Linear Units) function as activation function. Dropout technique was also

18

applied in the second and third layer (with rate 0.4 for the second and 0.6 for the third) to reduce

19

over-fitting. The loss function is mean squared error (MSE). We use the Adam optimizer

20

(Kingma & Ba, 2014) with learning rate 0.0001 and other parameters used in the original paper.

21

2.3.1.2 Random Forest model

22

Similarly, the Random Forest was also applied in this paper. About the 204 properties calculated

23

by the Genetic Function Approximation (GFA) algorithm, we used Pearson correlation

24

coefficient matrix to analyze the relationship between these properties. The result was presented 6 ACS Paragon Plus Environment


Page 8 of 64

1

as a picture with the python package Yellowbrick. About the 30 data set, we used principal

2

component analysis including 2D and 3D to analyze the relationship between properties. The

3

preprocessing of the data set involved the following steps. Firstly, the numerical value of

4

properties was be zoomed to 0 to 1, then the 54 properties with variance greater than 0.05 were

5

be chosen. About these 54 properties, we did normalization to get the data with mean is 0 and

6

variance is 1. Finally, we got 9 properties using Lasso feature selection. The Pearson correlation

7

coefficient matrix heatmap of 9 properties shows the correlation coefficient of 9 properties

8

selected is small enough.

9

2.3.1.3 SVM and MLR model

10

The regression of models were generated by support vector machine (SVM) and multiple linear

11

regression (MLR), which were generated by LibSVM and MATLAB, respectively. The equation

12

(2) expresses the squared correlation coefficient (r 2) of SVM: 31

13

14 15

2

r =

̅ ̅ ̅ (l̅ ∑li=1 f(xi )yi −∑li=1 f(xi ) ∑li=1 yi )

2

̅ ̅ ̅ ̅ (l̅ ∑li=1 f(xi )2 −(∑li=1 f(xi ))2)(l̅ ∑li=1 yi 2 −(∑li=1 yi )2)

(2)

The expression of MLR model is described by equation (3): Y = b0 + b1×X1 + b2×X2 + ...+ b7×X7

(3)

16

In particular, lead drugs of GSK3β ATP binding site were also utilized to generate predict model

17

for candidate validation, a number of 117 chemical compounds (Table 1, Table S1) were used to

18

build property-DockScore relationship models.

19

were calculated with the activity site of GSK3β ATP binding site by using docking study.

20

2.3.2 3D-QSAR model

21

26 GSK3β inhibitors including 21 training molecules and 5 test molecules with bioactivity (Ki

22

value) from Berg's study 28 were utilized to build CoMSIA model with Sybyl-X 1.1. (Table 1).

23

The CoMSIA model was constructed according to their activity and structural characteristics

32, 33

The DockScore of all chemical compounds


Page 9 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

including stereoscopic fields, electric fields, hydrophobic fields, and hydrogen bonds for

2

acceptors and donors. The models were evaluate through cross-validation (CV). The residual

3

between the measured and predicted values were calculated. The square sum of explained

4

(SSE), F-test, coefficient of determination (R2), and q2 for cross-validation were referred as

5

evaluation indexes of CoMSIA model. Various field combinations had been studied separately.

6

2.4 MD simulation

7

To prepare the protein-ligand system for MD simulations, the topology files and parameters for

8

the ligands were obtained from SwissParam web server. 34 The periodic cubic box for simulation

9

had the margin of 1.2 nm and was solvated using TIP3 water molecules. To mimic the

10

physiological condition, we included 0.145 M Na+ and Cl-ions to neutralize the system charge.

11

The systems were minimized for 5,000 steps using steepest descent followed by equilibration in

12

canonical ensemble (NVT, T = 310K) for 100 ps with positional constraints applied and stages of

13

equilibration in isothermal-isobaric ensemble (NPT) to release the constraints. We extended the

14

NPT production runs to 100,000 ps. The van der Waals interactions were evaluated by Lennard-

15

Jones potential 35. The bond lengths between each pair of atoms were constrained by the linear

16

constraint solver (LINCS) algorithm. We sampled every 20 ps over the simulations to analyze the

17

variation of our protein-ligand systems using GROMACS 4.5.5. All simulations were performed

18

using the charmm27 force field. 36, 37

19

3. Results and discussion

20

3.1 Potential target proteins for Alzheimer’ disease form PPI network analysis

21

Protein-protein interaction (PPI) is known to be a critical component in cell signaling transduction,

22

thereby, plays important roles in the pathogenesis of many diseases. Over the years, researchers

23

have acquired a substantial amount of data on PPI in various forms, which contributed greatly to

24

the quest for understanding the molecular mechanism for diseases and the ultimate remedies. To 8 ACS Paragon Plus Environment


Page 10 of 64

1

identify target proteins for Alzheimer’s disease, we utilized website of STRING v10.5 to construct

2

a PPI network for target protein selection. (Figure 2) We then performed functional enrichments

3

analysis with this PPI network to sort out the protein interactions associated with Alzheimer’s

4

disease. The resulting proteins were ranked/categorized by interaction score, amongst them the

5

GSK3β, NCSTN, MAPT and PSEN1 were related with Alzheimer’s disease from Kyoto

6

Encyclopedia of Genes and Genomes (KEGG) database (Table S2), with false discovery rate

7

(FDR) values all below 0.05. We then interrogated the relationships between the four proteins and

8

their related pathway, and constructed an association network using Cytoscape with the data of

9

the functional enrichment from KEGG database (Figure 3). It then became clear that the Notch

10

signaling pathway, Neurotrophin signaling pathway and Wnt signaling pathway (colored in red)

11

converged to AD though the four proteins. Many research have been reported that these three

12

signaling pathway play important roles in pathogenesis of AD. 38-41 Therefore, GSK3β, NCSTN,

13

MAPT and PSEN1 are potential effective molecular targets for disrupting the pathogenic pathway

14

of Alzheimer’s disease, and we believe it is advisable to focus our drug selection on them.

15

3.2 Results of predicted model generation

16

3.2.1 Deep Learning model discussion

17

In this paper, we applied deep learning method on our research and achieve R2 of 0.927 on

18

training set and 0.862 on the test set. We used a simple 4-layers full connected neural network

19

for our research (Figure 4). The activation function is ReLu (Rectified Linear Units). Our

20

training samples are fewer, when training neural networks, the trained models are prone to over-

21

fitting. Dropout can be used as a trick to train deep neural networks. In each training batch,

22

over-fitting can be significantly reduced by ignoring half of the feature detectors (making half of

23

the hidden layer nodes have a value of 0).Dropout technique is applied in the second and third

24

layer (with rate 0.4 for the second and 0.6 for the third) to reduce over-fitting. The loss function 9 ACS Paragon Plus Environment

Page 11 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

is mean squared error (MSE). The Adam optimizer is suitable for unstable objective functions

2

and its hyper parameters are well interpreted and usually need little or no adjustment. By setting

3

Adam optimizer’s learning rate as 0.0001. We have done 120 times different attempts and

4

achieve credible results. (Figure 5).

5

3.2.2 Random Forest model discussion

6

As for Random Forest model, the Pearson Ranking of 204 Features illustrates that high correlation

7

between some features, which are bigger than 0.4(Figure 6). In order to reduce the number of

8

related features, replace the original variable with fewer variables. Principal component analysis

9

has no effect if the original variables are orthogonal to each other, ie, there is no correlation. The

10

Two-dimensional principal component analysis and Three-dimensional principal component

11

analysis demonstrate dimensionality reduction conduce to hunt out the most representative

12

features (Figure 7). We scaled the number of 204 properties and filtered from them with variance

13

greater than 0.05 to find out the most representative 54 properties. Then, we made the Lasso

14

feature selection to get 9 features with small correlation coefficient and good orthogonality

15

(Figure 8). Up to 20 trees are used for training, and samples that have not been extracted are used

16

as verification sets. The mean square error (MSE) on the training set is 0.122, on the test set is

17

0.049. The R2 on the training set is 0.869, on the test set is 0.890(Figure 9).

18

3.2.3 MLR and SVM models discussion

19

Among the property-activity relationship models we built by using Genetic Function

20

Approximation (GFA) algorithm, the most reasonable one (Table 2) contains eight molecular

21

properties: logD, Jurs_PPSA_3, Jurs_RPSA, Jurs_TPSA, Jurs_WPSA_3, Minimized_Energy,

22

PMI_X, and Shadow_XYfrac. The multiply linear function of this model is expressed in equation

23

(8)



Page 12 of 64

1

pKi = -28.251 + 0.79395 × logD + 1.292 × Jurs_PPSA_3 − 85.727 × Jurs_RPSA + 0.15375 ×

2

Jurs_TPSA − 1.2504 × Jurs_WPSA_3 + 0.14322 × Minimized_Energy + 0.0055567 × PMI_X +

3

25.714 ×Shadow_XYfrac

4

In equation (8), logD is the octanol-water partition coefficient in pH 7.4, accounting for the

5

ionization states of each ligand. Jurs_PPSA_3 is the positive surface area weighted by all

6

positively charged atoms of each ligand. Jurs_RPSA stands for the polar surface area divided by

7

all molecular of solvent-accessible surface area. Jurs_TPSA is the total polar surface area

8

(summation of solvent-accessible surface areas of atoms). Jurs_WPSA_3 represents the surface-

9

weighted charged partial surface areas. Minimized_Energy is the energy after a fast minimization

10

procedure by using clean force field. PMI_X exhibits the orientation and conformational rigidity

11

of each ligand. Shadow_XYfrac indicates the area of the molecular shadow of each ligand in the

12

xz plane. The properties obtained from GFA algorithm showed that properties of each ligand with

13

ionization states, solvent-accessible surface areas, and orientation of scaffold were significantly

14

related to experimental pKi values.

15

In order to construct a bioactivity predictive model for each molecule in TCM database, we

16

computed a score by using these eight properties for all the training and test molecules. We then

17

validated the prediction accuracy of these two models by using SVM and MLR models displayed

18

accurate prediction ability indicated by the high values of correlation coefficient (r 2) observed in

19

the SVM model (r2 = 0.7925) and MLR model (r2 = 0.8942) (Figure 10). Hence, we utilized the

20

SVM and MLR models to predict the bioactivities of compounds from TCM database.

21

We also utilize GFA algorithm to select reasonable property for property-DockScore relationship

22

models generation (Table 3), the best model contains fifteen molecular properties are related to

23

DockScote (R2 = 0.8494), and the multiply linear function of this model is expressed in equation

24

(9):

25

DockScore = -19.655 + 3.5212 *ALogP − 1.2877 × ES_Count_aaCH + 2.0698 × ES_Count_aaN

(8)


Page 13 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

+ 4.6548 × ES_Count_sssN + 0.31246 × ES_Sum_tN + 2.4761 × HBA_Count + 2.6922

×

2

Num_Rings + 2.1511 × Num_RotatableBonds + 0.15699 × Molecular_PolarSASA − 0.35 ×

3

Molecular_PolarSurfaceArea − 4.2065 × CHI_3_C − 11.063 × IAC_Mean − 1.4235 × Kappa_1

4

+ 0.97326 × Dipole_mag + 29.415 × Jurs_FPSA_2 + 392.97 × Jurs_FPSA_3 + 0.26191 ×

5

Jurs_PNSA_1 + 0.047399 × Jurs_WNSA_2 − 0.040872 × Strain_Energy + 18.979 ×

6

Shadow_XYfrac

7

In equation (9), ALogP is log of the octanol-water partition coefficient. The ES_Count property

8

is calculate the E-state count for atom type in the molecule. The ES_Sum property is calculate the

9

E-state sum for atom type in the molecule. HBA_Count is the number of hydrogen bond accepting

10

groups. Num_Rings is defined as the number of rings in the smallest set of smallest rings.

11

Num_RotatableBonds is the number of Rotatable bonds in the molecule. Molecular_PolarSASA

12

is the value of the polar solvent accessible surface area for each molecule. CHI_3_C is a class of

13

connectivity indices of topological descriptors. IAC_Mean is Graph-Theoretical InfoContent

14

descriptors. Kappa_1 means Kappa Shape Indices. Dipole_mag is 3D electronic descriptors in an

15

electrostatic field. Jurs_FPSA descriptor is fractional charged partial surface areas. Jurs_PNSA_1

16

is sum of the solvent-accessible surface areas of all negatively charged atoms. Jurs_WNSA_2

17

indicates surface-weighted charged partial surface areas. Strain_Energy is the difference between

18

Energy and Minimized_Energy. Shadow_XYfrac indicates the area of the molecular shadow of

19

each ligand in the xz plane.

20

The 117 inhibitors of GSK3β ATP binding site were used to calculate the value of the fifteen

21

properties. The SVM and MLR algorithm were employed to construct DockScore predicted

22

models. We further validate the prediction accuracy of these two models, the high value of

23

correlation coefficient (r2) displayed in the SVM model (r2 = 0.971) and MLR models (r2 = 0.822).

24

CoMSIA model of 3D-QSAR were constructed (Figure 11). The value of CoMSIA were

25

calculated and further to obtain the predicted activity. Partial least square (PLS) analysis and

(9)



Page 14 of 64

1

validation of CoMSIA model were provided the value of correlation coefficient (R2), standard

2

error of estimate (SEE) and F test (F ratio) for evaluating, as well as the cross validation (q2cv)

3

for assessing (Table 4). Large enough R2 and relatively small SEE could explain the rationality

4

of the model. And the q2cv which greater than 0.5 were worth being considered. CoMSIA model

5

suggested the hydrophobic group (blue area), hydrogen bonds donor (Magenta area) for drug

6

design. Areas that were inappropriate for introducing hydrogen bond donor would be warned

7

correspondingly (green area) (Figure 11B). External validations of the models were provided to

8

identify its accuracy.

9

3.3 TCM candidates selection via database screening and bioactivity prediction

10

The compounds from TCM database were used in the high throughput virtual screening against

11

GSK3β, Nicastrin and APBB1. Many small molecule GSK3β inhibitors have been reported in

12

previous studies.

13

prediction and filter out reasonable TCM candidates (Table 5). The ADMET description of each

14

candidate were calculated (Table 6). As Nicastrin and APBB1 have no sufficient experimental

15

data to determine the bio-activity of chemical compounds and there is no available control set to

16

compare with TCM candidates, the QSAR model and MD simulation cannot be used to analysis

17

chemical compounds. Therefore, the only way to investigate potent ligand is calculate binding

18

score through docking study (Table S3). The drug-target network showed Methyl 3-O-

19

feruloylquinate and Morusimic acid B have potential to target multi-protein (Figure 12). Methyl

20

3-O-feruloylquinate, Morusimic acid B and Cynanogenin A had high vote score among the top

21

eight candidates in the results of docking study and bioactivity prediction. Their chemical

22

scaffolds are displayed in Figure 13 and he docking pose for the chosen three candidates and

23

GSK3β complex are presented in Figure 14. These three TCM candidates have a common

24

property, a carboxylate group in their scaffolds. This functional group might increase the binding

25

affinity for R141 of GSK3β. Analyzing the docking poses revealed that all small molecules form

42

So they could be used to construct a QSAR model to performed the activity


Page 15 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

H-bond interactions with residues on GSK3β. Methyl 3-O-feruloylquinate interacts with Q185 of

2

GSK3β through H-bond interactions. Furthermore, residue I62 and R141, interact with docked

3

ligand, are close to Methyl 3-O-feruloylquinate. Morusimic acid B has H-bond interactions with

4

N186 and R141. Cynanogenin A forms two H-bonds with Q185 and R141. 7YG engages with

5

K85 and D133 through 2 H-bonds. Interestingly, there is no H-bond formed between 7YG and

6

residues of Q185 and R141. From the Table 5, the multiple QSAR validations (SVM, MLR,RF,DL

7

and CoMSIA), we integrate all the five methodology and select out the potent compounds(Table

8

7).

9

3.4 Stability of the protein-ligand complexes

10

To efficiently search over the large conformational space for the possible protein-ligand

11

interaction mode, most docking studies treat the receptor as a rigid body, which often results in

12

docking poses don’t necessarily represent the protein-ligand interactions under dynamic condition.

13

Therefore, we performed MD simulations with the protein-ligand complexes recovered from

14

docking to validate the stability of them by monitoring the structural variation over the

15

simulations. We plotted of root mean square deviation (RMSD) of each complex with reference

16

to the corresponding first frame in the simulation trajectory and shown them in Figure 15. The

17

fluctuation in the plot of RMSD tended to stabilize after 5 ns (Figure 15A), and the value of

18

RMSD is in average of 0.23 nm. The RMSD evolution profile of the GSK3β complex with Methyl

19

3-O-feruloylquinate, Morusimic acid B, or Cynanogenin A, resembles that of 7YG (Figure 16A),

20

suggesting an intrinsic stability of the protein structure in the complexes.

21

We also computed the radius of complex gyration for all the complexes over the simulations

22

(Figure 16B). The value of complex gyration is inversely proportional to the compactness of the

23

structure. The compactness of the GSK3β complexes with the three TCM candidates revealed less

24

variation than that of the GSK3β-7YG complex. The gyration values of GSK3β complexes with

25

Morusimic acid B, and Cynanogenin A did not change to a great extend during the simulation. 14 ACS Paragon Plus Environment


Page 16 of 64

1

The lowest value of gyration appeared in the GSK3β-Methyl 3-O-feruloylquinate complex after

2

80 ns of simulation, which indicated the protein structure of GSK3β become more compact at

3

final stage of MD simulation.

4

To analyze the energy of all simulation system with solvent and protein complexes, we calculated

5

total energy of simulation system during 100 ns (Figure 17), the value of total energy is sum of

6

potential energy and kinetic energy. We found that the total energy of all simulated systems are

7

within the average of -1.012 × 106 over all simulation time, this results denoted that situation of

8

simulation system is constant, this phenomenon also illustrated the solvent are not affect the

9

protein structure of GSK3β during MD process. Because all of the simulation systems of GSK3β

10

complexes are stable during MD simulation, we further analyzed the structural variation of each

11

ligand by RMSD calculation (Figure 18). We found that Cynanogenin A has large increased values

12

after 5ns. For 7YG, Methyl 3-O-feruloylquinate, and Morusimic acid B, the ligand RMSD value

13

were tend to stable during all simulation time. Therefore, we further measured distances between

14

interacted atoms in next studies to observe the variation of Cynanogenin A in GSK3β binding site.

15

3.5 H-bond distance analysis form dynamics protein-ligand complexes

16

To analyze the variation of Cynanogenin A, we found that two residues, R141 (ARG141) and

17

R144 (ARG144), can generate interactions. These residues are important for Methyl 3-O-

18

feruloylquinate, Morusimic acid B and Cynanogenin A to bind to GSK3β. Consequently, we

19

calculate the distance between Cynanogenin A and residues R141, R144 and D200 (ASP200) of

20

GSK3β throughout the simulation (Figure 19). About the H-bond distance analysis of

21

Cynanogenin A, we also found that the H-bond distances of R141 and R144 were increased after

22

30 ns and 60 ns, respectively. On the other hand, the distance between Cynanogenin A and D200

23

was decreased during 70 ns to 90 ns. The data illustrated that the Cynanogenin A was close to

24

D200 in GSK3β binding site; hence, we further analyzed the binding conformation of

25

Cynanogenin A in the binding site of GSK3β at initial stage (0 ns) and final stage (100 ns) of the 15 ACS Paragon Plus Environment

Page 17 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

MD simulation. The snapshots of Cynanogenin A were shown in Figure 14. As for snapshot of

2

initial stage, the Cynanogenin A interacted with D141 and D144 in GSK3β binding site, which

3

forms two H-bond interactions on both of the residues. At final stage of MD simulation time, the

4

position of Cynanogenin A was close to D200. Because the D200 located in GSK3β binding site,

5

the data of H-bond distance revealed that Cynanogenin A move into the binding pocket at final

6

stage of MD simulation time.

7

3.6 Tunnels prediction in protein structure and migration of docked ligands

8

Here, we analyzed the tunnels and channels in the structure of protein-ligand complexes for all

9

MD conformations (Figure 20), the possible channels were produced from the place of each 43

10

docked ligand by using CAVER 3.0 software.

11

individual clusters. The width of produced channels was very small in either of GSK3β-7YG

12

(Figure 20A) or Methyl 3-O-feruloylquinate (Figure 20B), which indicated that these two ligands

13

might not easily escape from the binding site. On the contrary, the produced channels from

14

GSK3β-Morusimic acid B (Figure 20C) performed more distinct predicted tunnels, suggesting

15

that Morusimic acid B has a higher chance to escape from the docking site. Interestingly, we found

16

the GSK3β-Cynanogenin A complex featured one channel and displayed from inside of protein

17

structure, which was colored in light blue channel in Figure 20D. This tunnel resides inside of

18

protein structure, suggesting that Cynanogenin A move from outside into GSK3β binding pockets.

19

We then utilized the DSSP analysis tool of GROMACS to delineate the evolution of the secondary

20

structure features of the protein-ligand complexes in our simulations (Figure 21). All our GSK3β-

21

ligand complexes displayed stable secondary structure during over the entire course of our

22

simulations (Figure 21A-C). It is worth noting that Cynanogenin A induced changes to the

23

secondary structure of residue 65 to 70 in the window of from 0 ns to 15 ns (Figure 21D),

24

suggesting that the binding of Cynanogenin A with the active site of GSK3β affected surrounding

25

residues, and altered the stability of secondary structure.

The calculated tunnels were group into five



Page 18 of 64

1

In the end, we measured the distance between protein and ligand to illustrate the migration of

2

each ligand over the simulations by using Ligplot plus (Figure 22).

3

GSK3β and Morusimic acid B exibits a large fluctuation 30 ns into the free MD simulation, which

4

suggested that Morusimic acid B cannot bind to GSK3β stably, the distance between protein and

5

ligand is in the average of 3.0 nm at final stage of MD simulation time. However, the distance

6

between GSK3β and Cynanogenin A reduced from 2.0 nm to 1.5 nm during 100 ns simulation

7

time. The phenomenon of Cynanogenin A is related with H-bond distance analysis, tunnels

8

prediction and DSSP analysis. During the process of Cynanogenin A migrated into the binding

9

site of GSK3β at final stage of MD simulation, the reduced distance of H-bond appeared between

10

D200 and Cynanogenin A, one possible channel was observed from the results of tunnels

11

prediction, and the results of DSSP analysis showed that secondary structure of binding region

12

was changed from residue 65 to 70. We also analyzed the MD snapshot with 2D diagram to

13

observer the variation of protein-ligand interactions at initial state (0 ns) and final state (100 ns)

14

in Figure 23. The key residues, R141 and R144, could interacted with each ligand at 0 ns, but

15

Morusimic acid B change interacted residues at 100 ns (Figure 23F), which was related to the

16

result of distance analysis from protein and each ligand. 7YG and Methyl 3-O-feruloylquinate did

17

not change the binding residues among all simulation times. The binding residues still had

18

hydrophobic and H-bond interactions with 7YG and Methyl 3-O-feruloylquinate. Since the

19

binding position of Cynanogenin A changed to binding site of GSK3β during dynamic condition,

20

the numbers of binding residues are increased at final stage of MD simulation (Figure 23H). These

21

results suggested that TCM candidates, Methyl 3-O-feruloylquinate and Cynanogenin A, might

22

be regarded as potential inhibitors for GSK3β. According to the information of herbal resource in

23

TCM database, Methyl 3-O-feruloylquinate could be isolated from Phellodendron amurense and

24

Stemona japonica,45, 46 Cynanogenin A from Cynanchum atratum.

25

Phellodendron amurense is Amur cork tree, which has been used to tread hematochezia, diabetes, 17 ACS Paragon Plus Environment

44

47

The distance between

The common name of

Page 19 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


48

1

mouth sores, and drench tur-bidity in medical therapy.

2

pathogenic microorganisms.

3

been used to treat skin inflammatory in folk medicine. 49

4

4. Conclusion

5

Identifying drug targets used to be an utterly empirical and experimental procedure, in which

6

researchers delved into literatures to fish out potential candidates and set up experiments to prove

7

or disprove, which were often laborious. Nowadays, as large-scale databases of protein-protein

8

interaction are established and ever increasing in the amount of information curated in them, the

9

quest for potential drug targets has become unprecedentedly efficient by integrating them into the

10

process. With the purpose of tackling Alzheimer’s disease, we performed the pathway network

11

analysis on PPI network to find out potential drug targets for it. Our PPI network studies revealed

12

the importance of GSK3β and Nicastrin in the pathogenesis of Alzheimer’s disease, the result of

13

network analysis highlighted the substantial involvement of these two proteins in Alzheimer’s

14

disease. With these the protein targets in hand, we opted to unearth natural compounds of potential

15

inhibitory effect to them from known herbs, mainly due to the readiness and safety associated

16

herb medicine application. Our group has made a major contribution in unifying traditional

17

Chinese medicine with modern medicinal chemistry by developing the TCM, which contains all

18

known herbs and the corresponding active compounds in it. Among the compounds from TCM

19

database, we determined that Methyl 3-O-feruloylquinate and Cynanogenin A were capable of

20

binding to GSK3β by predicting the bioactivities and performing molecular docking. We further

21

evaluated our results by MD simulations. The MD simulations shown that Cynanogenin A had

22

unstable Ligand RMSD after MD simulation, but the MD analyses showing that Cynanogenin A

23

migrated into the binding pocket after 5 ns. Methyl 3-O-feruloylquinate and Cynanogenin A both

24

induced more compactness to GSK3β complex and stable condition among all simulation times,

25

the GSK3β complex also had no substantial fluctuation after simulation time of 5 ns. In addition,

46

Stemona japonica has potential to anti-

Cynanchum atratum belong to Apocynaceae family, which has



Page 20 of 64

1

the MD simulation analysis displayed that these two GSK3β complexes can generate stable

2

conformations in the end of MD simulation time. We then identified that Methyl 3-O-

3

feruloylquinate could be isolated from Phellodendron amurense and Stemona japonica 45, 46 and

4

Cynanogenin A is available in Cynanchum atratum.

5

formula comprises Phellodendron amurense and Cynanchum atratum, to be applied in

6

Alzheimer’s disease treatment.

7

Acknowledgments

8

This work was supported by Guangzhou science and technology fund (Grant No 201803010072)

9

and China Medical University Hospital (DMR-106-151, DMR-106-071). We also acknowledge

47

Therefore, we propose this novel herb

10

the start-up funding from SYSU “Hundred Talent Program”.

11

Disclosure

12

The author reports no conflicts of interest in this work.

13 14 15 16 17 18 19 20 21 22 23 24 25


Page 21 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2

References 1. Smits, L. L.; Flapper, M.; Sistermans, N.; Pijnenburg, Y. A.; Scheltens, P.; van der

3

Flier, W. M., Apraxia in Mild Cognitive Impairment and Alzheimer's Disease: Validity and

4

Reliability of the Van Heugten Test for Apraxia. Dementia and Geriatric Cognitive Disorders

5

2014, 38 (1-2), 55-64.

6

2.

Tan, S. H.; Karri, V.; Tay, N. W. R.; Chang, K. H.;

7

S.;

Keh, H. W.; Candasamy, M., Emerging pathways to neurodegeneration: Dissecting the

8

critical molecular mechanisms in Alzheimer’s disease, Parkinson’s disease. Biomedicine &

9

Pharmacotherapy 2019, 111, 765-777.

Ah, H. Y.; Ng, P. Q.; Ho, H.

10

3.

Cavallini, A.; Brewerton, S.; Bell, A.; Sargent, S.; Glover, S.; Hardy, C.; Moore,

11

R.; Calley, J.; Ramachandran, D.; Poidinger, M.;

12

Szekeres, P.; Bose, S., An unbiased approach to identifying tau kinases that phosphorylate tau at

13

sites associated with Alzheimer disease. Journal of Biological Chemistry 2013, 288 (32), 23331-

14

47.

15

4.

16

beta-induced mitochondrial dysfunction involves the inhibition of GSK-3beta. J Alzheimers Dis

17

2012, 32 (4), 981-96.

18

5.

19

Alzheimer's disease. Bioorg Med Chem Lett 2019, 29 (2), 125-133.

20

6.

21

genetic perspective. Cell 2005, 120 (4), 545-55.

22

7.

23

generation level by BACE1 enzymatic activity and transcription. FASEB Journal 2006, 20 (2),

24

285-92.

25

8.

26

relationship (SAR) studies of synthetic glycogen synthase kinase-3β inhibitors: A critical review.

Karran, E.; Davies, P.; Hutton, M.;

Huang, H. C.; Xu, K.; Jiang, Z. F., Curcumin-mediated neuroprotection against amyloid-

Fish, P. V.; Steadman, D.; Bayle, E. D.; Whiting, P., New approaches for the treatment of

Tanzi, R. E.; Bertram, L., Twenty years of the Alzheimer's disease amyloid hypothesis: a

Li, Y.; Zhou, W.; Tong, Y.; He, G.; Song, W., Control of APP processing and Abeta

Xu, M.; Wang, S. L.; Zhu, L.;

Wu, P. Y.; Dai, W. B.; Rakesh, K. P., Structure-activity



Page 22 of 64

1

European Journal of Medicinal Chemistry 2019, 164, 448-470.

2

9.

3

Drug Discov 2004, 3 (6), 479-87.

4

10. Eldar-Finkelman, H., Glycogen synthase kinase 3: an emerging therapeutic target. Trends

5

Mol Med 2002, 8 (3), 126-32.

6

11. Reddy, P. H., Amyloid beta-induced glycogen synthase kinase 3beta phosphorylated VDAC1

7

in Alzheimer's disease: implications for synaptic dysfunction and neuronal damage. Biochimica

8

et Biophysica Acta 2013, 1832 (12), 1913-21.

9

12. Tsygvintsev, A., On the overfly algorithm in deep learning of neural networks. Applied

Cohen, P.; Goedert, M., GSK3 inhibitors: development and therapeutic potential. Nat Rev

10

Mathematics and Computation 2019, 349, 348-358.

11

13. Cai, C.;

12

Deep Learning-Based Prediction of Drug-Induced Cardiotoxicity. Journal of Chemical

13

Information and Modeling 2019.

14

14. Pei, J.; Zheng, Z.; Merz, K. M., Random Forest Refinement of the KECSA2 Knowledge-

15

Based Scoring Function for Protein Decoy Detection. Journal of Chemical Information and

16

Modeling 2019.

17

15. Chen, Y.-C., Beware of docking! 2014; Vol. 36.

18

16. tsung-Ying, T.; Chang, K.-W.; Yu-Chian Chen, C., IScreen: World's first cloud-computing

19

web server for virtual screening and de novo drug design based on TCM database@Taiwan. 2011;

20

Vol. 25, p 525-31.

21

17. Chang, K.-W.;

22

Sun, M.-F.;

23

Computing Web Server for Traditional Chinese Medicine for Online Virtual Screening, de novo

24

Evolution and Drug Design. Journal of Biomolecular Structure and Dynamics 2011, 29 (1), 243-

25

250.

Guo, P.;

Zhou, Y.;

Tsai, T.-Y.;

Chen, H.-Y.;

Zhou, J.;

Chen, K.-C.;

Wang, Q.;

Zhang, F.;

Yang, S.-C.;

Fang, J.; Cheng, F.,

Huang, H.-J.;

Chang, T.-T.;

Tsai, F.-J.; Chen, C. Y.-C., iSMART: An Integrated Cloud


Page 23 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

18. Chen, C. Y., TCM Database@Taiwan: the world's largest traditional Chinese medicine

2

database for drug screening in silico. PLoS One 2011, 6 (1), e15939.

3

19. Boran, A. D.; Iyengar, R., Systems approaches to polypharmacology and drug discovery.

4

Current opinion in drug discovery & development 2010, 13 (3), 297-309.

5

20. Zhang, B.;

6

Novel Macropinocytosis Inducer by a Network Target Approach. Frontiers in pharmacology

7

2018, 9, 10.

8

21. Asad, Y.; Ahmad, S.; Rungrotmongkol, T.; Ranaghan, K. E.; Azam, S. S., Immuno-

9

informatics driven proteome-wide investigation revealed novel peptide-based vaccine targets

10

against emerging multiple drug resistant Providencia stuartii. Journal of molecular graphics &

11

modelling 2018, 80, 238-250.

12

22. Talevi, A., Multi-target pharmacology: possibilities and limitations of the "skeleton key

13

approach" from a medicinal chemist perspective. Frontiers in pharmacology 2015, 6, 205.

14

23. Liu, Y. L.; Liu, J. H.;

15

Huang, C. M.; Lin, S. Y.; Chang, C. T.; Huang, C. C., Association of inflammatory cytokines

16

with mortality in peritoneal dialysis patients. BioMedicine 2017, 7 (1), 1.

17

24. Zhang, G. X.; Zhang, Y. Y.; Zhang, X. X.; Wang, P. Q.; Liu, J.;

18

Different network pharmacology mechanisms of Danshen-based Fangjis in the treatment of stable

19

angina. Acta pharmacologica Sinica 2018.

20

25. Li, S.;

21

of a Functional Formula of Three Chinese Medicinal Herbs: Experimental Evidence and Network

22

Pharmacology-Based Identification of Mechanism of Action and Potential Bioactive Components.

23

Molecules 2018, 23 (2).

24

26. Zhang, X.; Pi, Z.; Zheng, Z.; Liu, Z.; Song, F., Comprehensive investigation of in-vivo

25

ingredients and action mechanism of iridoid extract from Gardeniae Fructus by liquid

Wang, X.; Li, Y.;

Wu, M.;

Wang, S. Y.; Li, S., Matrine Is Identified as a

Wang, I. K.; Ju, S. W.;

Yu, T. M.; Chen, I. R.; Liu, Y. C.;

Liu, Q.; Wang, Z.,

Wang, N.; Hong, M.; Tan, H. Y.; Pan, G.; Feng, Y., Hepatoprotective Effects



Page 24 of 64

1

chromatography combined with mass spectrometry, microdialysis sampling and network

2

pharmacology. Journal of chromatography. B, Analytical technologies in the biomedical and life

3

sciences 2018, 1076, 70-76.

4

27. Szklarczyk, D.;

5

Santos, A.;

6

STRING database in 2017: quality-controlled protein-protein association networks, made broadly

7

accessible. Nucleic acids research 2017, 45 (D1), D362-D368.

8

28. Berg, S.;

Bergh, M.;

Hellberg, S.;

9

von Berg, S.;

Weigelt, T.;

Ormo, M.;

Morris, J. H.;

Doncheva, N. T.;

Cook, H.;

Roth, A.;

Kuhn, M.;

Bork, P.;

Simonovic, M.;

Jensen, L. J.; von Mering, C., The

Hogdin, K.; Xue, Y.;

Wyder, S.;

Lo-Alfredsson, Y.;

Tucker, J.;

Soderman, P.;

Neelissen, J.;

Jerning, E.;

10

Nilsson, Y.; Bhat, R., Discovery of novel potent and highly selective glycogen synthase kinase-

11

3beta (GSK3beta) inhibitors for Alzheimer's disease: design, synthesis, and characterization of

12

pyrazines. J Med Chem 2012, 55 (21), 9107-19.

13

29. Bai, X. C.;

14

conformational space of the catalytic subunit of human gamma-secretase. eLife 2015, 4.

15

30. Brooks, B. R.; Brooks, C. L., 3rd; Mackerell, A. D., Jr.; Nilsson, L.; Petrella, R. J.;

16

Roux, B.;

17

Cui, Q.;

18

Kuczera, K.; Lazaridis, T.; Ma, J.; Ovchinnikov, V.; Paci, E.; Pastor, R. W.; Post, C.

19

B.;

20

Yang, W.;

21

Comput Chem 2009, 30 (10), 1545-614.

22

31. Chang, C.-C.; Lin, C.-J., LIBSVM: A library for support vector machines. ACM Trans. Intell.

23

Syst. Technol. 2011, 2 (3), 1-27.

24

32. Saraswati, A. P.;

25

Glycogen synthase kinase-3 and its inhibitors: Potential target for various therapeutic conditions.

Rajendra, E.;

Won, Y.;

Archontis, G.;

Dinner, A. R.;

Pu, J. Z.;

Yang, G.;

Feig, M.;

Schaefer, M.;

Shi, Y.; Scheres, S. H., Sampling the

Bartels, C.; Fischer, S.;

Tidor, B.;

Boresch, S.; Gao, J.;

Venable, R. M.;

Caflisch, A.; Hodoscek, M.;

Woodcock, H. L.;

Caves, L.; Im, W.;

Wu, X.;

York, D. M.; Karplus, M., CHARMM: the biomolecular simulation program. J

Ali Hussaini, S. M.;

Krishna, N. H.;


Babu, B. N.; Kamal, A.,

Page 25 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

European journal of medicinal chemistry 2018, 144, 843-858.

2

33. Meijer, L.;

3

kinase 3. Trends in pharmacological sciences 2004, 25 (9), 471-80.

4

34. Zoete, V.; Cuendet, M. A.; Grosdidier, A.; Michielin, O., SwissParam: a fast force field

5

generation tool for small organic molecules. J Comput Chem 2011, 32 (11), 2359-68.

6

35. Jones, J. E., On the Determination of Molecular Fields. II. From the Equation of State of a

7

Gas. Proceedings of the Royal Society of London. Series A 1924, 106 (738), 463-477.

8

36. MacKerell, A. D.;

9

Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; Joseph-McCarthy, D.;

Flajolet, M.; Greengard, P., Pharmacological inhibitors of glycogen synthase

Bashford, D.;

Bellott, M.;

Dunbrack, R. L.;

Evanseck, J. D.; Kuchnir, L.;

10

Kuczera, K.;

Lau, F. T. K.; Mattos, C.; Michnick, S.; Ngo, T.; Nguyen, D. T.; Prodhom,

11

B.;

12

Watanabe, M.;

13

for Molecular Modeling and Dynamics Studies of Proteins. The Journal of Physical Chemistry B

14

1998, 102 (18), 3586-3616.

15

37. Foloppe, N.; MacKerell, J. A. D., All-atom empirical force field for nucleic acids: I.

16

Parameter optimization based on small molecule and condensed phase macromolecular target data.

17

J Comput Chem 2000, 21 (2), 86-104.

18

38. Hardy, J.; Israel, A., Alzheimer's disease: In search of [gamma]-secretase. Nature 1999, 398

19

(6727), 466-467.

20

39. De Ferrari, G. V.; Inestrosa, N. C., Wnt signaling function in Alzheimer's disease. Brain

21

research. Brain research reviews 2000, 33 (1), 1-12.

22

40. Inestrosa, N. C.; Varela-Nallar, L., Wnt signaling in the nervous system and in Alzheimer's

23

disease. Journal of molecular cell biology 2014, 6 (1), 64-74.

24

41. Hock, C.;

25

neurotrophin imbalances in Alzheimer disease: decreased levels of brain-derived neurotrophic

Reiher, W. E.;

Roux, B.;

Schlenkrich, M.;

Wiórkiewicz-Kuczera, J.;

Heese, K.;

Smith, J. C.;

Stote, R.;

Straub, J.;

Yin, D.; Karplus, M., All-Atom Empirical Potential

Hulette, C.;

Rosenberg, C.; Otten, U., Region-specific



Page 26 of 64

1

factor and increased levels of nerve growth factor in hippocampus and cortical areas. Archives of

2

neurology 2000, 57 (6), 846-51.

3

42. Valasani, K. R.;

4

synthesis, pharmacophore modeling, virtual screening, and molecular docking studies for

5

identification of novel cyclophilin D inhibitors. Journal of chemical information and modeling

6

2014, 54 (3), 902-12.

7

43. Chovancova, E.; Pavelka, A.;

8

Gora, A.; Sustr, V.; Klvana, M.;

9

J., CAVER 3.0: a tool for the analysis of transport pathways in dynamic protein structures. PLoS

Vangavaragu, J. R.;

Day, V. W.; Yan, S. S., Structure based design,

Benes, P.; Strnad, O.; Brezovsky, J.; Medek, P.;

Kozlikova, B.;

Biedermannova, L.; Sochor, J.; Damborsky,

10

Comput Biol 2012, 8 (10), e1002708.

11

44. Laskowski, R. A.; Swindells, M. B., LigPlot+: multiple ligand-protein interaction diagrams

12

for drug discovery. Journal of chemical information and modeling 2011, 51 (10), 2778-86.

13

45. Ida, Y.;

14

Phellodendron Amurense Bark. Phytochemistry 1994, 35 (1), 209-215.

15

46. Ge, F.;

16

X.; Zuo, J.; Ye, Y., Isolation of chlorogenic acids and their derivatives fromStemona japonica

17

by preparative HPLC and evaluation of their anti-AIV (H5N1) activityin vitro. Phytochemical

18

Analysis 2007, 18 (3), 213-218.

19

47. Bai, H.; Li, W.;

20

novel pregnane glycosides from Cynanchum atratum. Tetrahedron 2005, 61 (24), 5797-5811.

21

48. Akihisa, T.;

22

Koike, K.; Abe, M., Melanogenesis-Inhibitory and Cytotoxic Activities of Limonoids, Alkaloids,

23

and Phenolic Compounds from Phellodendron amurense Bark. Chemistry & biodiversity 2017,

24

14 (7).

25

49. Choi, Y. Y.; Kim, M. H.; Lee, H.;

Satoh, Y.;

Ohtsuka, M.;

Nagasao, M.; Shoji, J., Phenolic Constituents of

Ke, C.; Tang, W.; Yang, X.; Tang, C.; Qin, G.; Xu, R.;

Li, T.; Chen,

Koike, K.; Satou, T.; Chen, Y. J.; Nikaido, T., Cynanosides A-J, ten

Yokokawa, S.;

Ogihara, E.;

Matsumoto, M.;

Zhang, J.;

Kikuchi, T.;

Ahn, K. S.; Um, J. Y.; Lee, S. G.; Kim, J.; Yang, 25


Page 27 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

W. M., Cynanchum atratum inhibits the development of atopic dermatitis in 2,4-

2

dinitrochlorobenzene-induced mice. Biomedicine & pharmacotherapy = Biomedecine &

3

pharmacotherapie 2017, 90, 321-327.

4



Page 28 of 64

Table 1. Chemical scaffold of GSK3β inhibitors of Berg’s study for QSAR model generation. Core

Comp. R1

R2

R3

pKi

0

--

--

7.39

H

--

7.13

--

6.52

--

1 2*

H

3*

F

H

8.89

4

CH3

H

9.34

5

CF3

H

8.96

6

CH3

CH3

8.82

7

H

F

8.05

8

H

CH3

8.20

9

H

H

8.70

10

H

H

7.74

11

H

H

8.31

12*

H

H

8.36

13*

H

H

9.40

14

H

H

15*

H

H 27


SO2NH2

7.92 9.17

Page 29 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


16

H

H

7.70

17*

H

H

8.08

18

H

H

8.51

19

H

H

7.80

20

C

C

NH2

6.92

21*

N

C

NH2

7.19

22

N

N

H

6.43

23

--

9.66

24

--

9.32

25#

--

7.92

26

--

6.16

27#

--

9.00

28#

--

7.13

29#

--

7.05

30#

--

7.66

Compound 1-30 without * for SVM, MLR model, compound 0-21 for CoMSIA model, * test compound for SVM, MLR model, # test compound for CoMSIA.


Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46

Page 30 of 64

Table 2. The top 10 models to predicted related properties of chemicals of Berg’s study by GFA algorithm R2

Equation of Model Model_1=

-28.251 + 0.79395 * logD + 1.292 * Jurs_PPSA_3 − 85.727 * Jurs_RPSA + 0.15375 * Jurs_TPSA − 1.2504 * Jurs_WPSA_3 + 0.14322 * Minimized_Energy + 0.0055567 * PMI_X + 25.714 * Shadow_XYfrac

Model_2 = -113.98 + 0.79395 * logD + 1.292 * Jurs_PPSA_3 + 85.727 * Jurs_RASA + 0.15375 * Jurs_TPSA − 1.2504 * Jurs_WPSA_3 + 0.14322 * Minimized_Energy + 0.0055567 * PMI_X + 25.714 * Shadow_XYfrac Model_3 = -113.71 + 0.67595 * logD + 3.5206 * Jurs_PPSA_3 + 12.758 * Jurs_RPSA + 0.13481 * Jurs_SASA − 4.6307 * Jurs_WPSA_3 + 0.13264 * Minimized_Energy + 0.0053795 * PMI_X + 21.569 * Shadow_XYfrac Model_4 = -100.95 + 0.67595 * logD + 3.5206 * Jurs_PPSA_3 − 12.758 * Jurs_RASA + 0.13481 * Jurs_SASA − 4.6307 * Jurs_WPSA_3 + 0.13264 * Minimized_Energy + 0.0053795 * PMI_X + 21.569 * Shadow_XYfrac Model_5 = -22.063 + 0.57547 * logD − 0.15598 * Num_Chains + 259.93 * Jurs_FPSA_3 − 68.494 * Jurs_RPSA + 0.12083 * Jurs_TPSA + 0.15364 * Minimized_Energy + 0.0061851 * PMI_X + 27.547 * Shadow_XYfrac Model_6 = -90.557 + 0.57547 * logD − 0.15598 * Num_Chains + 259.93 * Jurs_FPSA_3 + 68.494 * Jurs_RASA + 0.12083 * Jurs_TPSA + 0.15364 * Minimized_Energy + 0.0061851 * PMI_X + 27.547 * Shadow_XYfrac Model_7 = -26.332 + 0.72646 * logD − 1070.6 * Jurs_FPSA_3 + 3.7603 * Jurs_PPSA_3 + 0.024697 * Jurs_TPSA − 2.6268 * Jurs_WPSA_3 + 0.14699 * Minimized_Energy + 0.0061307 * PMI_X + 25.083 * Shadow_XYfrac Model_8 = -33.207 + 0.68468 * logD − 1390.9 * Jurs_FPSA_3 − 0.050945 * Jurs_PNSA_3 + 4.9397 * Jurs_PPSA_3 − 3.4294 * Jurs_WPSA_3 + 0.1476 * Minimized_Energy + 0.0062177 * PMI_X + 29.63 * Shadow_XYfrac


0.863

0.863

0.861

0.861

0.861

0.861

0.860

0.859

Page 31 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46


Model_9 = -24.483 − 0.94715 * ES_Count_sNH2 + 0.67404 * logD − 0.12282 * Num_Chains + 0.41936 * Jurs_PPSA_3 + 0.023046 * Jurs_TPSA + 0.1639 * Minimized_Energy + 0.0073735 * PMI_X + 28.03 * Shadow_XYfrac Model_10 = -31.01 − 1.1799 * ES_Count_sNH2 + 0.86682 * logD + 0.78994 * Jurs_PPSA_3 + 0.029885 * Jurs_TPSA − 0.4277 * Jurs_WPSA_3 + 0.15829 * Minimized_Energy + 0.0073822 * PMI_X + 27.477 * Shadow_XYfrac


0.856

0.856


Page 32 of 64

Table 3. The top 10 models to predicted related properties of chemicals of Berg’s study by GFA algorithm R2

Equation of Model Model_1=

-28.251 + 0.79395 * logD + 1.292 * Jurs_PPSA_3 − 85.727 * Jurs_RPSA + 0.15375 * Jurs_TPSA − 1.2504 * 0.863 Jurs_WPSA_3 + 0.14322 * Minimized_Energy + 0.0055567 * PMI_X + 25.714 * Shadow_XYfrac

Model_2 = -113.98 + 0.79395 * logD + 1.292 * Jurs_PPSA_3 + 85.727 * Jurs_RASA + 0.15375 * Jurs_TPSA − 1.2504 * 0.863 Jurs_WPSA_3 + 0.14322 * Minimized_Energy + 0.0055567 * PMI_X + 25.714 * Shadow_XYfrac Model_3 = -113.71 + 0.67595 * logD + 3.5206 * Jurs_PPSA_3 + 12.758 * Jurs_RPSA + 0.13481 * Jurs_SASA − 4.6307 * 0.861 Jurs_WPSA_3 + 0.13264 * Minimized_Energy + 0.0053795 * PMI_X + 21.569 * Shadow_XYfrac Model_4 = -100.95 + 0.67595 * logD + 3.5206 * Jurs_PPSA_3 − 12.758 * Jurs_RASA + 0.13481 * Jurs_SASA − 4.6307 * 0.861 Jurs_WPSA_3 + 0.13264 * Minimized_Energy + 0.0053795 * PMI_X + 21.569 * Shadow_XYfrac Model_5 = -22.063 + 0.57547 * logD − 0.15598 * Num_Chains + 259.93 * Jurs_FPSA_3 − 68.494 * Jurs_RPSA + 0.12083 * 0.861 Jurs_TPSA + 0.15364 * Minimized_Energy + 0.0061851 * PMI_X + 27.547 * Shadow_XYfrac Model_6 = -90.557 + 0.57547 * logD − 0.15598 * Num_Chains + 259.93 * Jurs_FPSA_3 + 68.494 * Jurs_RASA + 0.12083 * 0.861 Jurs_TPSA + 0.15364 * Minimized_Energy + 0.0061851 * PMI_X + 27.547 * Shadow_XYfrac


Page 33 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46


Model_7 = -26.332 + 0.72646 * logD − 1070.6 * Jurs_FPSA_3 + 3.7603 * Jurs_PPSA_3 + 0.024697 * Jurs_TPSA − 2.6268 * 0.860 Jurs_WPSA_3 + 0.14699 * Minimized_Energy + 0.0061307 * PMI_X + 25.083 * Shadow_XYfrac Model_8 = -33.207 + 0.68468 * logD − 1390.9 * Jurs_FPSA_3 − 0.050945 * Jurs_PNSA_3 + 4.9397 * Jurs_PPSA_3 − 3.4294 * 0.859 Jurs_WPSA_3 + 0.1476 * Minimized_Energy + 0.0062177 * PMI_X + 29.63 * Shadow_XYfrac Model_9 = -24.483 − 0.94715 * ES_Count_sNH2 + 0.67404 * logD − 0.12282 * Num_Chains + 0.41936 * Jurs_PPSA_3 + 0.856 0.023046 * Jurs_TPSA + 0.1639 * Minimized_Energy + 0.0073735 * PMI_X + 28.03 * Shadow_XYfrac Model_10 = -31.01 − 1.1799 * ES_Count_sNH2 + 0.86682 * logD + 0.78994 * Jurs_PPSA_3 + 0.029885 * Jurs_TPSA − 0.4277 * 0.856 Jurs_WPSA_3 + 0.15829 * Minimized_Energy + 0.0073822 * PMI_X + 27.477 * Shadow_XYfrac



Page 34 of 64

Table 4 PLS Analysis and Validation of CoMSIA Models. q2cv

ONC

r2

SEE

F ratio

S

0.363

5

0.942

0.214

48.840

H D A

0.523 -0.120 0.325

5 5 5

0.973 0.566 0.835

0.147 0.587 0.362

106.686 3.919 15.211

S+H S+D S+A

0.526 0.307 0.375

5 5 5

0.971 0.940 0.920

0.142 0.218 0.252

99.824 47.212 34.643

H+D H+A D+A

0.545 0.477 0.336

5 5 5

0.977 0.950 0.889

0.134 0.199 0.297

129.898 57.303 23.977

S+H+D S+H+A S+D+A

0.500 0.466 0.394

5 5 5

0.973 0.952 0.935

0.146 0.195 0.227

108.700 59.928 43.147

H+D+A S+H+D+A

0.462 0.468

5 5

0.952 0.955

0.196 0.189

59.928 63.420

Parameter

q2cv: Correlation Coefficient (cross validation) r2: Correlation Coefficient (non-cross validation) ONC: Optimal number of components SEE: Standard error of estimate F ratio: F-test value


Page 35 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46


S: Steric H: Hydrophobic D: Hydrogen bond donor A: Hydrogen bond acceptor H+D Filed proportion: H: 90.0%; D:10.0%



Page 36 of 64

Table 5 Docking results of top ten candidates and 7YG, activity value of each chemical molecule was predicted by QSAR models. Predicted value (pKi)

d

SVM

MLR

RF

DL

Nodifloridin A

8.325

9.128

8.064

8.073

137.700

0

4

0

1

4

2,2'-[Benzene-1,4diylbis(methanediylo xybenzene-4,1-

8.366

8.806

7.451

8.501

126.304

0

4

1

1

4

7.561

8.933

8.591

6.063

105.216

0

3

0

0

4

4

methoxyphenyl) propionic acid] Prostaglandin E1

7.689

8.245

8.077

10.670

98.746

0

3

0

1

4

5 6 7

Tianshic acid Methyl 3-O-feruloylquinate Morusimic acid B

7.908 9.664 8.148

9.560 7.584 9.233

7.833 8.584 8.820

7.655 8.291 9.098

95.808 94.757 82.394

0 1 0

3 4 3

0 0 0

0 1 0

4 4 4

8

Cynanogenin_A 7-Methoxy-aristolochiac acid

8.694 8.082

9.443 8.327

8.296 8.487

6.523 10.304

81.451 79.165

0 0

3 4

0 1

0 1

4 3

1

2

3

9

Name

diyl]bis(oxoacetic acid) Dihydroferulic acid[3-(4-hydroxy-3-

a

Absorption bBBB


c

CYP2D6

Hepatoto xicity

e

Dock Score

Ranking

Solubility Level

Page 37 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46


10 Control set

3,15,19-trihydroxy8(17),13-ent-labdad ien-16-oic acid

8.367

7.381

8.623

5.935

71.920

0

3

0

1

4

7YG

7.499

7.341

7.852

7.304

67.118

1

4

0

1

3

Notes: a

Predicted by Structure-binding affinity predictive model

Abbreviations: SVM, support vector machine; MLR, multiple linear regressions; CoMSIA, comparative molecular similarity indices analysis.



Page 38 of 64

Table 6 The predicted ADMET descriptors of top eight candidates and 7YG.

Name

a

Absorption

b

c

BBB

CYP2D6

d

Hepatotoxicity

e

Solubility Level

Nodifloridin A

0

4

0

1

4

2,2'-[Benzene-1,4diylbis(methanediylo

0

4

1

1

4

xybenzene-4,1-diyl]bis(oxoacetic acid) Dihydroferulic

0

3

0

0

4

acid[3-(4-hydroxy-3-methoxyphenyl) propionic acid] Prostaglandin E1 Tianshic acid

0 0

3 3

0 0

1 0

4 4

Methyl 3-O-feruloylquinate

1

4

0

1

4

Morusimic acid B

0

3

0

0

4

Cynanogenin_A 7-Methoxy-aristolochiac acid

0 0

3 4

0 1

0 1

4 3

3,15,19-trihydroxy-8(17),13-ent-

0

3

0

1

4

labdad 37 ACS Paragon Plus Environment

Page 39 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46


ien-16-oic acid 7YG

1

4

0

1

3

a

Absorption: Good absorption = 0; Moderate absorption = 1; Low absorption = 2;

b

BBB (Blood Brain Barrier): Very high penetration = 0; High penetration = 1; Medium penetration = 2; Low penetration = 3; Undefined

penetration = 4 c

CYP2D6: Non-inhibitor = 0, Inhibitor = 1

d

Hepatotoxicity: Non-inhibitor = 0, Inhibitor = 1

e

Level 3: Drug-likeness with good solubility; Level 4: Drug-likeness with optimal solubility



Page 40 of 64

Table 7 Vote score of top eight candidates

Vote score

pKi

Dock Score

SVM

MLR

RF

DL

1

1

0

0

0

0

0

0

2

1

0

0

1

0

0

0

0

2

0

0

1

0

0

1

0

0

2

Prostaglandin E1 Tianshic acid Methyl 3-O-feruloylquinate

0 0 1

0 1 0

0 0 1

1 0 1

0 0 1

0 0 0

1 0 0

0 0 1

2 1 5

Morusimic acid B Cynanogenin_A

0 1

1 1

1 1

1 0

1 0

0 1

1 0

1 0

6 4

Name Nodifloridin A 2,2'-[Benzene-1,4diylbis(methanediyloxybenzene4,1-diyl]bis(oxoacetic acid) Dihydroferulic acid[3-(4-hydroxy3-methoxyphenyl) propionic acid]

CoMSIA

SVM

MLR

Multi- Totaltarget score

Vote score: For all activity values predicted by one algorithm, the top 50% were voted 1 point, and others were voted as 0 point.


Page 41 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure legends

Figure 1. Flowchart of experiment process. Abbreviations: TCM, traditional Chinese medicine; AD, Alzheimer's disease; GFA, genetic function approximation; SVM, support vector machine; MLR, multiple linear regressions; PLS, partial least square; CoMSIA, comparative molecular similarity indices analysis. MD, molecular dynamics;



Figure 2. PPI networks obtained from the Cytoscape STRING app. Each network nodes represent proteins. Colorful network nodes are first shell of interactors for AD.


Page 42 of 64

Page 43 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3. Association network generated form functional enrichments data under Cytoscape tool. The investigated matching proteins were ellipse, and correlative pathways for the matching proteins were rectangle. The width of each edge in this network was depended on the value of FDR score, a thick edge has low FDR score than a thin edge.



Figure 4. Deep learning network diagram including four layers, the activation is ReLu, and the optimizer is Adam Optimizer with learning rate is 0.0001.


Page 44 of 64

Page 45 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 5. Scatter plots to present the results of 120 experiments.



Figure 6. The relation of 204 features ranked by Pearson Correlation Coefficient.


Page 46 of 64

Page 47 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


A

B

Figure 7. The Principal component analysis (PCA) visualization by Yellowbrick . (A) 2D, (B) 3D. 46 ACS Paragon Plus Environment


Figure 8. The Pearson correlation coefficient matrix heatmap of 9 selected features.


Page 48 of 64

Page 49 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 9. Scatter plots to identify relationships between observed activity and predicted activity for Random Forest model, all unit of activity value were represented by pKi. The values of MSE of the model on train and on test are 0.122 and 0.049, respectively. The values of R2 of the model on train and on test are 0.869 and 0.890, respectively.



Figure 10. Scatter plots to identify relationships between observed activity and predicted activity for (A) SVM and (B) MLR models, all unit of activity value were represented by pKi. The values of R2 of the two models are 0.793 and 0.894, respectively. 49 ACS Paragon Plus Environment

Page 50 of 64

Page 51 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 11. CoMSIA model constructed by GSK3β inhibitors. (A) Scatter plots to validate CoMSIA model, (B) CoMSIA model.



Figure 12. The results of drug-target network analyze from database screening, the protein targets were denoted by red rectangle.


Page 52 of 64

Page 53 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 13. The chemical scaffold of (A) Methyl 3-O-feruloylquinate, (B) Morusimic acid B, (C) Cynanogenin A and (D) 7YG



Figure 14. The docking poses of (A) Methyl 3-O-feruloylquinate (green), (B) Morusimic acid B (orange), (C) Cynanogenin A (purple) and (D) 7YG (yellow) in GSK3β binding site.


Page 54 of 64

Page 55 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 15. The snapshots of Cynanogenin A in GSK3β binding site at (A) 0ns and (B) 100ns, the docked ligand are colored in purple.



Figure 16.

Figure 16. The MD trajectory analysis of (A) Protein RMSD and (B) Radius of Gyration for GSK3β and TCM ligands: 7YG, Methyl 3-O-feruloylquinate, Morusimic acid B, and Cynanogenin A during 100ns simulation time.


Page 56 of 64

Page 57 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 17. The total energy calculations of GSK3β and TCM ligands: 7YG, Methyl 3-O-feruloylquinate, Morusimic acid B, and Cynanogenin A during 100ns simulation time.



Figure 18. The MD trajectory analysis of Ligand RMSD for 7YG, Methyl 3-Oferuloylquinate, Morusimic acid B, and Cynanogenin A in GSK3β structure among 100ns simulation time.


Page 58 of 64

Page 59 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 19. The H-bond distance between Cynanogenin A and residues of ARG141, ARG144 and ASP200 during 100ns simulation times of GSK3β.



Figure 20. The possible tunnels identified from all protein-ligand complexes with: (A) 7YG (B) Methyl 3-O-feruloylquinate, (C) Morusimic acid B, and (D) Cynanogenin A, by 100ns MD simulation times. All of tunnels grouped into three clusters, which are colored by red, blue and green, respectively.


Page 60 of 64

Page 61 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 21. DSSP analysis of five snapshots among 20 ns simulation times: proteinligand complexes of (A) 7YG, (B) Methyl 3-O-feruloylquinate, (C) Morusimic acid B, (D) Cynanogenin A. Residues from 35 to 85 on GSK3β are colored in blue, which were used to observe the variation of secondary structures during all MD simulations.



Figure 22. Distance between GSK3β and docked ligands during 100 ns simulation times.


Page 62 of 64

Page 63 of 64 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 23. Snapshot analysis of GSK3β complexes for comparison between initial (0 ns) and final states (100 ns) from MD simulation: (A) Initial state and (B) final state of 7YG. (C) Initial state and (D) final state of Methyl 3-O-feruloylquinate. (E) Initial state and (F) final state of Morusimic acid B. (G) initial state and (H) final state of Cynanogenin A. Blue chemical scaffold were used to display each ligand in 2D 62 ACS Paragon Plus Environment


diagram. Hydrophobic residues were represented by eyelash. H-bond interactions were indicated by green dash line. All difference of each residue from comparison between initial and final states was marked by red circle.


Page 64 of 64

Deep Learning and Random Forest Approach for Finding the Optimal

Recommend Documents