Structural and molecular insights into tacrine-benzofuran hybrid induced inhibition of amyloid-β peptide aggregation and BACE1 activity

Abstract Amyloid-β (Aβ) aggregation and β-amyloid precursor protein cleaving enzyme 1 (BACE1) are the potential therapeutic drug targets for Alzheimer’s disease (AD). A recent study highlighted that tacrine-benzofuran hybrid C1 displayed anti-aggregation activity against Aβ42 peptide and inhibit BACE1 activity. However, the inhibition mechanism of C1 against Aβ42 aggregation and BACE1 activity remains unclear. Thus, molecular dynamics (MD) simulations of Aβ42 monomer and BACE1 with and without C1 were performed to inspect the inhibitory mechanism of C1 against Aβ42 aggregation and BACE1 activity. In addition, a ligand-based virtual screening followed by MD simulations was employed to explore potent new small-molecule dual inhibitors of Aβ42 aggregation and BACE1 activity. MD simulations highlighted that C1 promotes the non aggregating helical conformation in Aβ42 and destabilizes D23–K28 salt bridge that plays a vital role in the self-aggregation of Aβ42. C1 displays a favourable binding free energy (–50.7 ± 7.3 kcal/mol) with Aβ42 monomer and preferentially binds to the central hydrophobic core (CHC) residues. MD simulations highlighted that C1 strongly interacted with the BACE1 active site (Asp32 and Asp228) and active pockets. The scrutiny of interatomic distances among key residues of BACE1 highlighted the close flap (non-active) position in BACE1 on the incorporation of C1. The MD simulations explain the observed high inhibitory activity of C1 against Aβ aggregation and BACE1 in the in vitro studies. The ligand-based virtual screening followed by MD simulations identified CHEMBL2019027 (C2) as a promising dual inhibitor of Aβ42 aggregation and BACE1 activity. Communicated by Ramaswamy H. Sarma


Introduction
Alzheimer's disease (AD) is a neurodegenerative disease that is mainly characterized by extracellular senile plaques comprising deposits of amyloid-b (Ab) peptide and intracellular neurofibrillary tangles (NFTs) consisting of hyperphosphorylated tau protein aggregates in the AD affected brain (LaFerla et al., 2007;Stancu et al., 2014).The World Health Organization (WHO) reported that approximately 55 million people worldwide have dementia and the number will increase to 78 million by 2030 (Gauthier et al., 2021).Although AD is a multifactorial disease, however, blocking Ab aggregation was highlighted as the most prominent therapeutic approach for delaying the onset of AD (Ahn et al., 2019;Di Santo et al., 2012;Novick et al., 2012).The sequential endoproteolytic splitting of amyloid precursor protein (APP) by BACE1 and c-secretases yield Ab, a peptide consisting of 37-43 amino acids (Kang et al., 1987;Mattson, 1997).Among various Ab fragments, Ab 42 is the predominant component in the extracellular senile plaques (Findeis, 2007).
The splitting of APP by BACE1 results in the formation of a soluble N-terminal ectodomain along with a 99-residue Cterminal fragment (CTF99), which is cleaved by c-secretase to produce Ab (Yan & Vassar, 2014).BACE1 has been recognized as the first rate-limiting enzyme in Ab biogenesis (Xu et al., 2012).BACE1 comprises N-terminal ectodomain, C-terminal and an aspartic protease domain.The active site consisting of Asp32 and Asp228 is located at the interface of N-and Cterminal lobes.A long, flexible b-hairpin loop consisting of residues Val67-Glu77 known as a flap in N-terminal lobe is present over the aspartic dyad.Apart from the aspartic dyad and flap region of BACE1, various loops (insert À A, À B, À C, À D, À E and À F) of BACE1 play a critical role in the recognition of the substrate and its subsequent binding (Barman et al., 2011).Lastly, a b-hairpin loop (Lys9 À Tyr14) known as a 10s loop can assume 'up' or 'down' conformation that is affected by its interaction with Thr232 (McGaughey et al., 2007).
Zha et al. identified a MTDL, C1 (tacrine-benzofuran hybrid), as a potent inhibitor of Ab aggregation (61.3%),BACE1 activity (IC 50 ¼ 1.35 lM) and human acetylcholinesterase (hAChE) (IC 50 ¼ 0.86 nM) (Zha et al., 2016).Among a library of 26 hybrids, C1 displayed good inhibitory activity against Ab aggregation, and BACE1 activity, displayed no significant hepatotoxicity and was noted to be significantly safer than tacrine.However, the underlying interaction mechanism of C1 with Ab and BACE1 remains unclear.In this work, MD simulations have been used to explore the key interactions and inhibition mechanism of C1 against Ab 42 aggregation and BACE1 activity.Additionally, lead compounds were identified from the library of bioactive compounds and approved drugs as dual inhibitors of Ab 42 aggregation and BACE1 activity using ligand-based virtual screening with C1 as a reference followed by MD simulations.

System preparation
The chemical structures of compounds were drawn using ChemDraw Ultra (Figure 1a; Mills, 2006).The optimization of compounds was done at the HF level using a 6-31(G) basis set by Gaussian09 (Frisch et al., 2009).The GROMOS96 parameters for compounds were generated using Automated Topology Builder (ATB) (Malde et al., 2011), which has been extensively used for the generation of topologies and parameters of various ligands compatible with GROMOS (Plazinski et al., 2016;Qu et al., 2016;Shuaib et al., 2017).
In this work, MD simulations of Ab 42 monomer and BACE1 alone as well as in the presence of C1 and lead compounds (C2 and C3) have been performed (Table 1).The PDB ID: 1Z0Q (Tomaselli et al., 2006;Figure 1b) was selected for MD simulations of Ab 42 monomer, which has been extensively employed to explore the conformational dynamics of Ab 42 peptide in the presence of various inhibitors, metal ions, etc. (Kalhor & Jabbari, 2017;Khatua et al., 2019;Kumari et al., 2020;Zou et al., 2019).The PDB ID: 1FKN (Hong et al., 2002), chain A (Figure 1c) has been chosen for MD simulations of BACE1, which has been employed in several studies to explore the interactions of various inhibitors with BACE1 (Ellis & Shen, 2015;Kumar et al., 2017;2019;Manoharan & Ghoshal, 2018).

Ligand-based virtual screening
The ligand-based virtual screening was performed using the SwissSimilarity tool with C1 as a reference compound (Zoete et al., 2016).The SwissSimilarity is a web-based tool used to screen small-molecule databases like ChEMBL, ChEBI, DrugBank, Ligand Expo, ZINC, etc. (Degtyarenko et al., 2008;Feng et al., 2004;Gaulton et al., 2012;Irwin & Shoichet, 2005;Wishart et al., 2006).The various small-molecule databases were screened by applying electroshape and spectrophores descriptors in the SwissSimilarity tool.The electroshape provides shape and electrostatic information, whereas spectrophores is a one-dimensional descriptor based on the property fields.A library of 707 compounds was generated based on the similarity score (only those compounds with a score higher than 0.6 or a similarity index of 60% were included in the library).The PyRx-virtual screening software was used to generate the docking files (PDBQT).The energy minimization of compounds was performed using Open Babel in PyRx before the evaluation of binding energy using AutoDock Vina (Dallakyan & Olson, 2015;Trott & Olson, 2010).Two metrics, the enrichment factor (EF) and receiver operating characteristic (ROC) curve, were employed to evaluate the enrichment success of the virtual screening strategy (Bender & Glen, 2005;Feher et al. 2003;Jain, 2004).The EF was evaluated using the following equation: Enrichment factor ðEFÞ ¼ ½TP � D�=½Ht � A� where true positive (TP)¼ the number of inhibitors correctly identified as actives, D¼ the total number of compounds in the database, Ht ¼ the total number of compounds retrieved by the model, and A ¼ the total number of active compounds in the dataset.The area under the curve (AUC) for ROC plots was computed by the trapezoidal integration method as implemented in the NCSS software (NCCS, 2022).

MD simulations and analysis
The conformational ensembles of all systems were generated using GROMACS (version 5.0.5;Abraham et al., 2015) with GROMOS96 54a7 force field (Lin & van Gunsteren, 2013;Van Der Spoel et al., 2005).Several studies have employed the GROMOS96 force field to explore the structural changes and conformational dynamics of peptides (de Oliveira et al., 2021;   Pal et al., 2020;Zein et al., 2022).The protonation state of Ab 42 monomer was chosen at physiological pH.The N-terminus of Ab 42 monomer was taken as NH 3 þ and the C-terminus was taken as COO -.All systems were soaked with simple point-charge (SPC) (Berendsen et al., 1981) water molecules (11,824 for Ab 42 monomer, 11,625 for Ab 42 monomer-C1, 11,632 for Ab 42 monomer-C2 and 11,636 for Ab 42 monomer-C3) in a cubic box of size 7.17 nm � 7.17 nm � 7.17 nm.The Na þ ions were added to the simulation box to keep the overall system neutral.
The BACE1 was protonated to match the experimental pH 4.5 using pKa values of the titratable residues.The PROPKA web tool was used to find the pKa values (Li et al., 2005;Søndergaard et al., 2011).The protonation state of Asp, Glu, Lys, His and Arg was assigned based on pKa values.The N-and C-terminal was taken as NH 3 þ and COO -, respectively.All systems were soaked with SPC water molecules (23,609 for apo-BACE1, 23,493 for BACE1-C1, 23,594 for BACE1-C2 and 23,597 for BACE1-C3) in a cubic box of dimension 9.21 nm � 9.21 nm � 9.21 nm.
For calculating long-range electrostatic interactions, a method known as particle mesh Ewald (PME) was used (Darden et al., 1993;Essmann et al., 1995).The LINCS algorithm was used to constrain the bond lengths (Hess et al., 1997).The equilibration of systems was performed first under NVT conditions for 100 ps at 300 K and then under NPT conditions for 100 ps.The modified Berendsen thermostat was used to keep up the temperature at 300 K (Bussi et al., 2007), and pressure was kept at 1 bar by the Parrinello-Rahman barostat (Parrinello & Rahman, 1981).
The GROMACS tools, dictionary of secondary structure of proteins (DSSP) program (Kabsch & Sander, 1983) PyMOL (DeLano, 2002), and visual molecular dynamics (VMD) (Humphrey et al., 1996) were employed to analyze the MD simulation trajectories.The chemical shifts (d sim ) for the Ca and Cb atoms of the representative conformation of the most-populated conformational cluster of Ab 42 monomer and BACE1 have been evaluated by SHIFTX2 (Han et al., 2011).The 3 J NHÀ Ha coupling constants were evaluated from dihedral angles / and w using parameters reported by Vuister and Bax (1993) in the Karplus equation (Karplus, 1959).The microstates were clustered by Daura et al. algorithm (Daura et al., 1999).The GROMACS tool gmx distance was used to calculate the distance between D23-K28 residues involved in the salt bridge.The distance between residue pairs was evaluated using gmx mdmat.The interaction energy of C1 as well as lead compounds (C2, C3 and C4) with Ab 42 monomer and BACE1 was determined by molecular mechanics Poisson Boltzmann surface area (MM-PBSA) method using g_mmpbsa (Kumari et al., 2014) in accordance with the previous study (Kaur et al., 2021).The GROMACS tool gmx sham was employed for the free energy landscape (FEL) maps.

Stabilization of Ab 42 monomer conformational ensemble in the presence of C1: insights from MD simulations
The molecular docking results indicated favourable binding (À 7.4 kcal/mol) of C1 with Ab 42 monomer (Table S1).C1 forms a hydrogen bond with Asp23 of Ab 42 monomer (Table S1, Figure S1a).C1 displayed hydrophobic contacts with Gln11, His14, Gln15, Val18, Phe19, Phe20, Glu22 and Asn27 of Ab 42 monomer (Table S1, Figure S1c).Further, MD simulations were carried out to explore the internal dynamics and conformational changes in Ab 42 peptide on the incorporation of C1.

Validation of Ab 42 monomer conformational ensemble with NMR data
The consistency of the MD simulation data was evaluated by comparing the chemical shifts and 3 J NHÀ Ha coupling constants of simulated and experimental data (Hou et al., 2004).A good correlation (R ¼ 0.96 for Ca and R ¼ 0.99 for Cb atoms) was observed between the calculated chemical shifts (d sim ) and experimental chemical shifts (d exp ) of Ab 42 monomer (Hou et al., 2004;Figure S2a-b).Additionally, Ab 42 monomer conformational ensemble was characterized using J-coupling ( 3 J NHÀ Ha ) constants (Vuister & Bax, 1993).For Ab 42 monomer, 3 J NHÀ Ha coupling constants displayed an average value of 5.0 Hz (Figure S3), which closely matched the average experimental value of 5.4 Hz (Tomaselli et al., 2006).A good correlation between d sim and d exp , 3 J NHÀ Ha coupling constants indicates that the conformational ensemble of Ab 42 monomer generated by MD matches with the experimental data.

C1 increases the structural stability of Ab 42 monomer
A noteworthy decrease in the total number of microstates from 95 to 30 in Ab 42 monomer was noted on the incorporation of C1 (Table S2), which indicates higher thermodynamic stability of Ab 42 peptide on the incorporation of C1.The population of the top three microstates in Ab 42 monomer was noted to be 58.3,18.2 and 8.4% (Figure 2), whereas it increased to 68.9, 21.1 and 9.0%, respectively, in Ab 42 monomer-C1 complex, which highlights enhanced stability of Ab 42 structure when C1 was added to it (Figure 2).The representative member of the top three microstates depicts higher helical content in Ab 42 peptide in the presence of C1, which suggests lower aggregation tendency of Ab 42 (Figure 2).The RMSD fluctuates with a value of �1.14 nm in Ab 42 monomer, whereas a lower value (�1.05 nm) was observed in Ab 42 monomer-C1 complex (Figure 3a).The RMSD values match with the RMSD of Ab 42 monomer with and without bleomycin (Kumari et al., 2020).The R g remains stable at a value of �1.15 nm for Ab 42 monomer and a marginally lower value (�1.09 nm) was noted for Ab 42 monomer on the incorporation of C1 (Figure 3b).The repeat simulations of Ab 42 monomer with and without C1 depicts almost identical RMSD, which indicates the reliability of MD simulations (Figure S4a, b).
The Ab 42 residues Asp1, Ala2, Arg5 À Phe20 and Ile31 À Met35 displayed lower fluctuations in the presence of C1 (Figure 3c), which highlight higher conformational stability of Ab 42 peptide in the presence of C1.The RMSF analysis is consistent with Ghorbani et al., which reported reduced fluctuations of Ab 42 monomer in the presence of phenolic compounds (Ghorbani et al., 2020).

C1 promotes the sampling of helical conformation in the conformational ensemble of Ab 42 monomer
The content of various secondary structures was analyzed in Ab 42 monomer with and without C1.The sampling of a-helix (34%) conformation is dominant followed by the random coil (21%) in Ab 42 monomer (Table 2, Figure S5).The population of turn and bend conformations was 15 and 11%, respectively (Table 2).On adding C1, the a-helix and coil conformations increased from 34 to 50% and 21 to 22%, whereas turn and bend conformations decreased from 15 to 11% and 11 to 8%, respectively.Notably, a remarkable increase in the a-helix content of Ab 42 monomer-C1 complex highlights the conservation of a-helix in Ab 42 peptide.
In Ab 42 monomer, the N-terminal residues Asp1-Tyr10  undergo a conformational transition into the coil and turn conformations (Figure S5, upper panel), whereas Asp1-Tyr10 sampled non aggregation-prone a-helix conformation in presence of C1 (Figure S5, lower panel).

C1 destabilizes D23-K28 salt bridge interaction and disrupts side chains contacts in Ab 42 monomer
The loop region that connects two b-sheets in Ab fibrillar models is stabilized by D23-K28 salt bridge (Berhanu & Hansmann, 2012;Tarus et al., 2006).In the case of Ab 42 alone, a distance peak at �0.32 nm indicated D23-K28 salt bridge interaction (Figure 4a), which is well below the reported distance of �0.46 nm (Truong et al., 2014).For Ab 42 monomer-C1 complex, a peak was observed at a higher value (�1.04 nm) and, notably, no peak was seen at �0.32 nm (Figure 4a).The disruption of D23-K28 salt bridge in Ab 42 peptide on the incorporation of C1 depicts a lower aggregation probability of Ab 42 peptide.The results match with Liu et al. that reported D23-K28 salt bridge disruption in Ab 42 monomer on the incorporation of Edaravone (Liu et al., 2020).
The intramolecular contacts between CHC, mid domain and C-terminal (Ans27 À Gly33) regions of Ab peptide lead to the formation of b-sheet conformation (Coskuner & Wise-Scira, 2013).A loss of intrapeptide side-chain contacts in Ab 42 monomer was noted on the incorporation of C1 (Figure S6).The contacts between the N-terminal, CHC, mid domain and C-terminal regions of Ab 42 peptide significantly decreased in presence of C1.The reduced intrapeptide contacts in Ab 42 monomer on the incorporation of C1 depicted reduced sampling of b-sheet conformation.

Impact of C1 on the hairpin formation in Ab 42 monomer and SASA analysis
The stabilization of hairpin conformation within the 14-36 segment of Ab 42 drives the self-assembly of Ab 42 peptide (Maity et al., 2017;Ma & Nussinov, 2002).To assess the effect of C1 on the hairpin conformation in the Ab 42 monomer, the angle (h) between Ca atoms of K28, A30 and I32 residues of Ab 42 monomer was evaluated.For Ab 42 monomer, the average value of angle (h) was noted to be 58.81� , whereas a significantly higher value (89.26 � ) of angle (h) was observed in Ab 42 monomer-C1 complex (Figure S7a), which, in turn, highlights disruption of the hairpin conformation of Ab 42 monomer in the presence of C1.Further, distances between D23 À L34 and A21 À V36 residue pairs in Ab 42 monomer were evaluated to assess the effect of C1 on the hairpin conformation in the Ab 42 monomer.The average value of the distance between D23 À L34 and A21 À V36 in the Ab 42 monomer was noted to be �1.32 and �1.44 nm, respectively (Figure S7b-c and Table 3).In contrast, a larger value of the distances between D23 À L34 (�1.56 nm) and A21 À V36 (�1.76 nm) was observed in Ab 42 monomer-C1 complex (Figure S7b-c and Table 3), which highlights destabilization of the hairpin conformation of Ab 42 peptide on the incorporation of C1.These are consistent with Ghorbani et al., which reported that phenolic compounds (Morin and Luteolin) destabilized the hairpin or U-shaped conformation of Ab 42 monomer by increasing the angle between Ca atoms of K28, A30, I32, and D23 À L34, A21 À V36 distances (Ghorbani et al., 2020).
To assess the effect of C1 on the structural and conformational preferences of Ab 42 monomer, solvent accessible surface area (SASA) was evaluated.The SASA for Ab 42 monomer-C1 complex (44.22 nm 2 ) was noted to be higher as compared to Ab 42 monomer (40.20 nm 2 ) (Figure S8a), which highlights that interactions of C1 with Ab 42 monomer weakened the intrapeptide hydrophobic contacts, which, in turn, decreased the aggregation tendency of Ab 42 monomer.Liu et al. (2020) also reported a higher SASA value of 40 nm 2 for Ab 42 monomer-Edaravone complex as compared to 35 nm 2 in Ab 42 monomer.A lower value of SASA for CHC residues of Ab 42 peptide was noted in the presence of C1 (Figure S8b).Thus, CHC residues remained buried and protected from solvent (Senguen et al., 2011), which is   consistent with the observation that C1 remains bound to the CHC region of Ab 42 peptide during simulation.
The per-residue binding free energy analysis highlighted the significant contribution of His14, CHC residues (Leu17, Val18, Phe20), Ala21 and Leu34 towards binding to C1 (Figure 5).The interaction energy analysis depicts the contribution of CHC residues of Ab 42 peptide in binding to C1.

Impact of C1 on the hydrogen bonds and FEL of Ab 42 monomer
An increase in the hydrogen bonds from �29.80 in Ab 42 monomer to �34.98 in Ab 42 monomer À C1 complex was noted (Figure S9a), which highlights stabilization of Ab 42 monomer structure on the incorporation of C1.Further, an average of �1.21 hydrogen bonds between Ab 42 monomer and C1 were observed in Ab 42 monomer À C1 complex (Figure S9b).In Ab 42 monomer, the lowest energy conformations i and ii are populated with 51 and 44% a-helix content (Figure 6a, Table 5), whereas significantly higher a-helix content (64%) was sampled in the lowest energy conformation of Ab 42 monomer À C1 complex (Figure 6b, Table 5).The FEL results highlight that C1 stabilized non aggregating helical conformation of Ab 42 peptide.

Identification of new inhibitors of Ab 42 aggregation
In this work, a ligand-based virtual screening approach has been employed to identify new inhibitors of Ab 42 aggregation from small-molecule databases using C1 as a reference (Figure 1a).A library consisting of 707 compounds was generated by ligand-based virtual screening of bioactive compounds and drugs from various small-molecule databases using SwissSimilarity tool (Zoete et al., 2016).Only those compounds that display a similarity score higher than 0.6 or a similarity index of 60% were included in the library.The EF and ROC curves were used to validate the virtual screening of the compound library.The enrichment factor (EF) is a tool to evaluate the quality of the docking and scoring as compared to the random selection in the virtual screening strategies.Enrichment factors may range from 1, which implies the random sorting of the compounds with no 'enrichment', to >100, in which only a slight percentage of the library needs to be screened to identify a higher number of active compounds.A higher EF value indicates the greater ability of a pharmacophore in identifying the active compounds from the database.The EF value was noted to be 74.9, which highlights the ability of the employed virtual screening methodology to correctly identify the active compounds.Yan et al. reported an EF value of 96 for a pharmacophore model, which indicates that the model is highly efficient in predicting the active molecules from the database (Yan et al., 2019).Further, the ROC curve depicted the AUC value of 0.67 and 0.58 for the Ab 42 monomer and BACE1, respectively (Figure S10), which indicated the reliability of the virtual screening strategy in the present study.
Next, the most promising compounds were selected from the library based on binding energy as well as key interactions in the aggregation-prone regions of Ab 42 monomer (Table S3, Figure S11).The molecular docking results indicated that C2 (À 8.2 kcal/mol), C3 (À 7.7 kcal/mol) and C4 (À 7.6 kcal/mol) bind to Ab 42 monomer with more affinity as compared to C1 (À 7.4 kcal/mol) and display hydrogen bonds, hydrophobic contacts with the key residues of Ab 42 monomer involved in the aggregation (Table S3, Figure S11).
The secondary structure content was analyzed in Ab 42 monomer on the incorporation of C2 and C3 (Table S5, Figure S14).In the presence of C2, the a-helix and coil conformations increased from 34 to 53% and 21 to 24%, whereas turn and bend conformations decreased from 15 to 13% and 11 to 9%, respectively (Table 2 and S5).A significant increase in the a-helix content of Ab 42 monomer in the presence of C2 highlights the conservation of a-helix in Ab 42 peptide.In the case of Ab 42 monomer-C3, the a-helix content was decreased from 34 to 29% with a concomitant increase in the p-helix content from 15 to 27% (Table 2 and  S5).The coil conformations increased from 21 to 31%, whereas turn and bend conformations decreased from 15 to 5% and 11 to 8%, respectively (Table 2 and S5).Notably, a remarkable increase in the p-helix content of Ab 42 monomer on the incorporation of C3 highlights the conservation of p-helix in Ab 42 peptide.
The loop region that connects two b-sheets in Ab fibrillar models is stabilized by D23-K28 salt bridge (Berhanu & Hansmann, 2012;Tarus et al., 2006).A distance peak at �0.32 nm indicated D23-K28 salt bridge interaction in Ab 42 monomer (Figure 4a).For Ab 42 monomer-C2 and Ab 42 monomer-C3, the peaks were observed at higher values (�1.12 and �0.50 nm) and, notably, no peak was seen at �0.32 nm (Figure S15).The disruption of D23-K28 salt bridge in Ab 42 peptide on the incorporation of C2 and C3 depicts a lower aggregation probability of Ab 42 peptide.The side chain-side chain contacts among CHC, mid domain and Cterminal (Ans27 À Gly33) residues of Ab peptide results in the formation of b-sheet conformation (Coskuner & Wise-Scira, 2013).The contacts between the N-terminal, CHC, mid domain and C-terminal regions of Ab 42 peptide significantly decreased in presence of C2 and C3 which depicts a reduced sampling of b-sheet conformation (Figure S6 and S16).The integrated computational approach employing ligand-based virtual screening followed by MD simulations identified C2 and C3 as promising lead compounds against Ab 42 aggregation.Further, C2 can be experimentally tested for its anti-  aggregation activity based on its high binding affinity, preservation of the a-helix content and strong tendency to destabilize the D23-K28 salt bridge in Ab 42 monomer.

C1 affects the structural stability of BACE1
The average RMSD for apo-BACE1 was noted to be 0.21 nm (Figure 7a), which matches with RMSD (0.21 nm) reported by Kumar et al. (2017).In the case of BACE1-C1 complex, the average RMSD fluctuates at a lower value (0.18 nm) during simulation that depicts higher structural stability of BACE1 in the presence of C1.The R g of BACE1 was observed to be marginally lower on the incorporation of C1 (Figure 7b).Further, RMSD and R g analysis for the repeat simulations highlighted almost similar values that indicate the reliability of MD simulations (Figure S18a-d).
The flap region and loops (insert À A, À B, À C, À D, À E and À F) situated close to the active site of apo-BACE1 possess conformational flexibility and exhibit substantial movements during substrate binding (Barman et al., 2011;Kumar et al., 2016;Manoharan & Ghoshal, 2018;Spronk & Carlson, 2011;Xiong et al., 2004).The RMSF analysis highlighted lower fluctu- Ser328 À Thr331, Glu339 and Phe341 À Phe374 residues of BACE1 in the presence of C1 (Figure 7c).The binding of C1 reduced the conformational fluctuations in flap residues and loops (inserts-A, À B, À C, À D, À E and À F) of BACE1, which compacted the active site and hindered the access of the substrate to the active site.

C1 influences the flap dynamics in BACE1
In BACE1, the flap movement plays the main role in the entry of the ligand to the aspartic dyad active site and the tight binding of the ligand in the active site (Hong & Tang, Table 6.The most promising compounds against Ab 42 aggregation ranked based on binding free energy (kcal/mol).The residues of the CHC region of Ab 42 monomer displaying significant contribution in binding with the lead compounds are depicted in bold.

Most promising compounds
2004; Shimizu et al., 2008).Thus, it is important to monitor the flap movement on the incorporation of ligand to evaluate its inhibitory potential on BACE1 activity.Xu et al. reported a value of �1.0 nm for the distance between the flap and Asp228 in BACE1, and mentioned that the flap is in open conformation during simulation (Xu et al., 2012).For BACE1, the distance between the flap and Asp228 was �1.03 nm (Figure 8), which highlights the open flap position in BACE1.In contrast, the distance was significantly reduced to �0.71 nm in BACE1 À C1 complex (Figure 8), which indicates the closed flap position of BACE1 on the incorporation of C1.
Several studies reported monitoring of interatomic distance between Ca Thr72 and Cb Asp32 atoms as an indicator of flap dynamics in BACE1 (Barman et al., 2011;Kumar et al., 2016;Manoharan & Ghoshal, 2018;Shimizu et al., 2008).An increase in the Ca Thr72-Cb Asp32 distance from 1.0 to 1.4 nm indicates movement from a close flap position to an open flap position (Dhanabalan et al., 2017;Manoharan & Ghoshal, 2018).In apo-BACE1, the interatomic distance between Ca Thr72 and Cb Asp32 atoms was noted to be 1.30 nm, which highlights the semi-open flap position (Figure 9a).For BACE1-C1 complex, the distance significantly decreased to �1.00 nm during simulation that depicts close flap (non-active) position (Figure 9a).The distance between Ca Tyr71 and Cb Asp32 atoms was calculated as Tyr71 was reported at the tip of the flap region in various studies (Barman et al., 2011;Hong et al., 2000;Spronk & Carlson, 2011).The distance between Ca Tyr71 and Cb Asp32 atoms of apo-BACE1 was noted to be �1.33 nm, whereas a gradual decrease in the interatomic distance to �0.94 nm was observed in the BACE1-C1 complex during simulation that indicates close flap position in BACE1 À C1 complex (Figure 9b).
Further, the distance between Ca Thr72 and Ca Thr329 atoms is considered an important parameter to monitor flap movement as it determines the gap between the C-terminal region and the tip of the flap of BACE1.In apo-BACE1, the Ca Thr72-Ca Thr329 distance was �1.02 nm (Figure 9c), which highlights the open flap position.In comparison, a significantly lower interatomic distance with an average value of �0.80 nm was observed in BACE1-C1 complex, which indicates a close flap position.The BACE1 residues Thr72 and Arg235 are known to interact with inhibitors (Hong et al., 2000), thus the distance between Oc1 Thr72 and NH1 Arg235 atoms has been evaluated in BACE1 with and without C1.A decrease in the Oc1 Thr72-NH1 Arg235 distance from 0.65 nm in apo-BACE1 to 0.61 nm in BACE1-C1 complex was observed (Figure 9d).The distance analyses between  important residues of BACE1 depict close flap position in BACE1 on the incorporation C1 that can also be visualized in the conformational snapshots of BACE1-C1 complex (Figure S19).
The phenyl ring of C1 displayed p À p stacking interactions with Tyr71 located in the flap region of BACE1.The average COM distance was observed to be 0.38 nm (Figure 10a), which highlights p À p stacking interactions as the distance is less than 0.44 nm (Alvarez, 2013;Janiak, 2000).Additionally, p À p stacking interactions were also observed between Phe108 of BACE1 and the phenyl ring of C1, and the average COM distance was observed to be 0.38 nm (Figure 10b).C1 is involved in hydrophobic contacts as well as p À p stacking interactions with flap, active site residues of BACE1 that restrict the flap dynamics and aid in the close (non-active) flap conformation (Figure 10c).Thus, C1 bind with active pocket and flap residues of BACE1, which, in turn, inhibits the BACE1 activity.

C1 alters the FEL of BACE1
The FEL of apo-BACE1 displayed three regions of minimum energy, whereas two minimum energy regions were noted in the BACE1 À C1 complex (Figure 11a).The characteristic conformations from minimum energy regions of BACE1 À C1 complex displayed key interactions of C1 with Asp32, Tyr71, Thr72, Phe108, Asp228 and Gly230 (Figure 11b), which matches with docking and binding free energy analyses.
The Ca Thr72-Cb Asp32 interatomic distance has been recognized as an important parameter to assess the flap movement in BACE1 (Dhanabalan et al., 2017;Kumar et al., 2017;Shimizu et al., 2008).For apo-BACE1, the distance between Ca Thr72 À Cb Asp32 atoms was noted to be 1.20, 1.19 and 1.17 nm for the minimum energy conformations i, ii and iii, respectively, which indicates semi-open flap position (Figure 11a).For BACE1 À C1 complex, the interatomic distance was 1.07 and 1.05 nm for i and ii, respectively, that depict close flap position (Figure 11b).Thus, flap adopt a semi-open position in apo-BACE1 and close position (nonactive) in the presence of C1.The FEL analysis confirmed strong binding of C1 in the BACE1 active site and thus responsible for the inhibition of BACE1 enzymatic activity.

Ligand-based virtual screening to identify new inhibitors of BACE1
The ligand-based virtual screening using C1 as a reference afforded a library of 707 compounds and this library was further employed to identify new compounds as inhibitors of BACE1.Among 707 compounds, the most promising compounds were selected based on binding energy (strong binding affinity) as well as key interactions with the residues of the flap region and various loops (inserts-A, À B, À C, À D, À E and À F) of BACE1.The compounds CHEMBL2019027 (C2), CHEMBL199670 (C3) and CHEMBL199861 (C4) bind to BACE1 with binding energies of À 10.6, À 10.4 and À 9.3 kcal/mol, respectively, as compared to À 8.5 kcal/mol noted in the case of C1 and display interactions with the residues of the flap region and various loops of BACE1 (Table S3).The binding free energies of C2, C3 and C4 with BACE1 are noted to be À 69.4 ± 6.2, À 66.5 ± 9.2 and À 39.6 ± 3.3 kcal/mol (Table 7), respectively, which highlight that C2 and C3 bind to BACE1 with more affinity as compared to C1 (-64.5 ± 18.0 kcal/mol).
The binding free energy analyses revealed Tyr71 of BACE1 as the most important residue that displays a significant contribution to the binding with C1, C2 and C3 (Table S6 and Figure S22).This is consistent with the molecular docking analysis that highlighted interactions of Tyr71 of BACE1 with C1, C2 and C3.Based on the binding free energy, C2 and C3 are selected for the MD simulation analysis to decipher interactions with BACE1.
In BACE1, the flap dynamics provide access to the ligand into the aspartic dyad active site (Hong & Tang, 2004;Shimizu et al., 2008).For BACE1, the distance between the flap and Asp228 was �1.03 nm (Figure 8), which highlights the open flap position in BACE1.In contrast, the distance was significantly reduced to �0.49 and �0.85 nm in BACE1 À C2 and BACE1 À C3 (Figure S24), respectively, which indicates the closed flap position of BACE1 on the incorporation of C2 and C3.In apo-BACE1, the interatomic distance between Ca Thr72 and Cb Asp32 atoms was noted to be 1.30 nm, which highlights the semi-open flap position (Figure 9a).The distance significantly decreased to �1.00 and �1.10 nm in BACE1-C2 and BACE1-C3, respectively, which indicate close flap (non-active) position (Figure S25a).The residue Tyr71 lies at the tip of the flap region and its distance with Asp32 was evaluated to monitor flap dynamics (Barman et al., 2011;Hong et al., 2000;Spronk & Carlson, 2011).The distance between Ca Tyr71 and Cb Asp32 atoms of apo-BACE1 was noted to be �1.33 nm (Figure 9b), whereas a gradual decrease in the interatomic distance to �0.97 and �1.00 nm was observed in the BACE1-C2 and BACE1-C3 during simulation that indicates close flap position in BACE1 À C2 and BACE1-C3 (Figure S25b).
To monitor flap movement, the distance between Ca Thr72 and Ca Thr329 atoms was evaluated.In apo-BACE1, the Ca Thr72-Ca Thr329 distance was �1.02 nm (Figure 9c), which highlights the open flap position.In comparison, a significantly lower interatomic distance with an average value of �0.78 and �0.96 nm was observed in BACE1-C2 and BACE1-C3, respectively, which indicates the close flap position (Figure S25c).Further, a decrease in the Oc1 Thr72-NH1 Arg235 distance from 0.65 nm in apo-BACE1 to 0.26 and 0.36 nm in BACE1-C2 and BACE1-C3, respectively, was observed (Figure 9d and Figure S25d).The distance analyses between important residues of BACE1 depict close flap position in BACE1 on the incorporation of C2 and C3.Thus, MD simulations highlighted C2 and C3 as promising lead compounds for the inhibition of BACE1 activity.

Conclusions
In this work, the inhibition mechanism of C1 against Ab 42 aggregation and BACE1 activity was explored using MD simulations.Further, new dual inhibitors of Ab 42 aggregation Table 7.The different components of binding free energies (kcal/mol) of C2, C3, and C4 with BACE1 evaluated by molecular mechanics Poisson Boltzmann surface area (MM À PBSA) method.
The molecular docking results indicated that C1 binds strongly (À 8.5 kcal/mol) with BACE1.The interatomic distance analysis between the flap and Asp228 implies a closed flap position (non-active) of BACE1 in BACE1 À C1 complex and an open flap position (active) in BACE1.The MM À PBSA analysis confirmed that C1 binds strongly with BACE1 (DG binding ¼ À 64.5 ± 18.0 kcal/mol) with a significant contribution of active pockets, aspartic dyad and flap residues.The FEL analysis depicted a semi-open flap position in apo-BACE1 and a close flap position in BACE1 À C1 complex.The MD simulations have provided key insights into the dynamics information as well as the interaction mechanism of C1 with Ab 42 monomer and BACE1.Further, ligand-based virtual screening followed by MD simulations identified a bioactive compound CHEMBL2019027 (C2) that binds to Ab 42 monomer with a higher binding affinity (DG binding ¼ À 65.1 ± 5.6 kcal/mol) than C1, preserve of the a-helix content of Ab 42 monomer and display strong tendency to destabilize the D23-K28 salt bridge in Ab 42 monomer.Notably, C2 also binds to BACE1 with high affinity (DG binding ¼ À 69.4 ± 6.2 kcal/mol) and the flap of BACE1 adopts a close position (non-active) in the presence of C2.MD simulations revealed that Val18 of Ab 42 monomer and Tyr71 of BACE1 are the most important residues that displayed significant contribution to the binding with C1 and C2.Thus, C2 can be further explored as a dual inhibitor of Ab 42 and BACE1 activity using in vitro and in vivo studies.

Figure 2 .
Figure 2. The central members of the top three microstates of Ab 42 monomer (panel a) and Ab 42 monomer À C1 complex (panel b).

Figure 3 .
Figure 3.The RMSD (panel a) and R g (panel b) of Ab 42 monomer with and without C1.The RMSF of each residue in Ab 42 monomer with and without C1 (panel c).
b-content is the sum of b-sheet and b-bridge.

Figure 4 .
Figure 4.The salt bridge probability in Ab 42 peptide in the absence and presence of C1 (panel a).The conformational snapshot depicting the distance between D23 and K28 in Ab 42 monomer À C1 complex (panel b).

Figure 5 .
Figure 5.The residue-wise binding energy of Ab peptide in Ab 42 monomer À C1 complex.

Figure 6 .
Figure 6.The FEL of Ab 42 monomer (panel a) and Ab 42 monomer À C1 complex (panel b).The lowest energy conformations are shown in the cartoon in cyan underneath FEL.

Figure 7 .
Figure 7.The time evolution of RMSD (panel a) and R g (panel b) for BACE1 with and without C1.The RMSF of each residue of BACE1 in the absence and presence of C1 (panel c).

Figure 8 .
Figure 8.The distance between Ca atoms of flap residues (Val67 À Glu77) and carboxyl O atoms of Asp228 is shown during simulation.

Figure 9 .
Figure 9.The interatomic distances among key residues of BACE1 in the absence (black) and presence (red) of C1.

Figure 10 .
Figure 10.The COM distance between Tyr71-C1 (panel a) and Phe108-C1 (panel b) during simulation.The conformational snapshot from the minimum energy basin of BACE1 À C1 complex displays p À p stacking interactions between BACE1 residues (Tyr71 and Phe108) and C1 (panel c).

Figure 11 .
Figure 11.The FEL for apo À BACE1 (panel a) and BACE1 À C1 complex (panel b).The high energy conformations are shown in orange and blue correspond to minimum energy conformations.The interatomic distance (nm) between Ca Thr72 and Cb Asp32 atoms of BACE1 in the representative conformations extracted from FEL is shown.

Table 1 .
Details of the simulated systems in this work.

Table 2 .
The secondary structure occupancies in Ab 42 peptide with and without C1.

Table 3 .
The angle (h) between Ca atoms of K28, A30 and I32 residues of Ab 42 monomer, and the D23 À L34 and A21 À V36 distances in Ab 42 monomer with and without C1.
a DE MM

Table 5 .
The secondary structure content of lowest energy conformations of Ab 42 monomer with and without C1.