Experimental and computational evaluation of dipeptidyl peptidase III inhibitors based on quinazolinone-Schiff’s bases

Abstract Dipeptidyl peptidase III (DPP III) is a zinc-dependent enzyme that sequentially hydrolyzes biologically active peptides by cleaving dipeptides from their N-termini. Although its fundamental role is not been fully elucidated, human DPP III (hDPP III) has been recognized in several pathophysiological processes of interest for drug development. In this article 27 quinazolinone-Schiff’s bases were studied for their inhibitory activity against hDPP III combining an in vitro experiment with a computational approach. The biochemical assay showed that most compounds exhibited inhibitory activity at the 100 μM concentration. The best QSAR model included descriptors from the following 2D descriptor groups: information content indices, 2D autocorrelations, and edge adjacency indices. Five compounds were found to be the most potent inhibitors with IC50 values below 10 µM, while molecular docking predicted that these compounds bind to the central enzyme cleft and interact with residues of the substrate binding subsites. Molecular dynamics simulations of the most potent inhibitor (IC50=0.96 µM) provided valuable information explaining the role of PHE109, ARG319, GLU327, GLU329, and ILE386 in the mechanism of the inhibitor binding and stabilization. This is the first study that gives insight into quinazolinone-Schiff’s bases binding to this metalloenzyme. Communicated by Ramaswamy H. Sarma


Introduction
Dipeptidyl peptidase III (DPP III; EC 3.4.14.4) is a zinc enzyme belonging to the metallopeptidase family M49, which sequentially hydrolyzes a number of biologically active peptides by cleaving dipeptides from their N-termini (Bar� sun et al., 2007;Chen & Barrett, 2004).Family M49 has been recognized as a distinct group based on the unique structural motifs HEXXGH and EEXR (K)AE(D) which are necessary for the binding of zinc cation (Zn 2þ ) to the enzyme active site and catalytic activity.Zinc cation is responsible for the activation of a water molecule that in the hydrolysis act as a nucleophile and attack the carbonyl carbon atom of the scissile peptide bond (Chen & Barrett, 2004).Insight into the crystal structure of human DPP III (hDPP III) without ligand showed that the molecule has an elongated shape (50 � 151 x 54 Å) with two domains separated by a wide cleft (Figure 1 left).The upper domain was mostly made of a-helix and the lower domain with a mixed a-helices and b-sheets.The active site with the zinc cation is located in the lower part of the upper domain and is directed towards the inter-domain cleft (Bezerra et al., 2012).According to previous experimental and computational approaches, it has been shown that the enzyme can be transformed from an extended (,,open") form (Figure 1 left) into a compact (,,closed") form (Figure 1 right), while the conformation of both protein domains has been preserved (Bezerra et al., 2012;Kumar et al., 2016;Tomi� c et al., 2012;2015).In addition, MD simulations of the extended form of hDPP III have shown that in solution protein is abundant in different conformations whose ,,semiclosed" form (Figure 1 middle) is the most presented (Tomi� c et al., 2012).
Quinazolinone nucleus is an important heterocyclic scaffold found in various biologically active natural products, synthetic compounds, pharmaceutical drugs and agrochemicals (Radwan & Alanazi, 2020).Schiff bases are a versatile family of compounds that carry imine or azomethine functional groups and exhibit a wide range of biological activities (Kajal et al., 2013).Thus, Schiff bases containing a quinazolinone core possess various biological activities such as antibacterial, antidiabetic, anti-inflammatory, and cytotoxic activity (Hricov� ıniov� a et al., 2018;Rakesh et al., 2015).Furthermore, both quinazolinone and Schiff bases derivatives have exhibited promising activity against some matrix metalloproteinases (Harney et al., 2012;Li et al., 2008).Previously, we have shown that heterocyclic compounds such as polyphenols (Agi� c et al., 2017;Blagojevi� c et al., 2021) and coumarins (Agi� c et al., 2021) have a potential to inhibit hDPP III by binding to a hydrophobic pocket and forming hydrogen bonds through OH groups with enzyme active site residues.The activities mentioned above encouraged us to experimentaly evaluate the effect of an unexplored class of heterocyclic compounds against hDPP III.
In this study we report the discovery of quinazolinone Schiff bases, synthesized and characterized earlier by Komar et al., 2021, that inhibit hDPP III at micromolar concentration.Computational techniques were used to describe the molecular basis of the observed activity of the most potent inhibitors.

Heterologous expression and purification of recombinant human DPP III
C-terminally truncated human DPP III was expressed and purified as described previously (Kumar et al., 2016).In brief, C-terminal HIS-tagged recombinant human DPP III was expressed in E. coli BL21-CodonPlus (DE3)-RIL cells (Agilent).
Overexpression was induced by 0.25 mM isopropyl-b-D-1-thiogalactopyranoside (IPTG, Fermentas) and the cells were grown 20 hours at 18 � C and 130 rpm.Obtained bacterial cells were lysed by a combination of lysozyme lysis and sonication, then treated with DNase I and centrifuged at 4000 g for 15 minutes to precipitate the cell debris.The lysate was purified by affinity chromatography on Ni-NTA column (5 mL prepacked His-trap FF, GE Healthcare).All fractions with hDPP III were pooled and incubated with TEV protease to remove the HIS-tag.hDPP III was recovered using flowthrough affinity chromatography (TEV protease is HIStagged), and additionally purified on a 16/60 Superdex-200 gel-filtration column (GE Healthcare).Fractions with hDPP III were pooled and desalted on PD-10 columns (GE Healthcare).Purified protein in 20 mM Tris HCl buffer pH ¼ 7.4 were stored at À 80 � C until use.

Human DPP III activity assay and IC 50 determination
Recombinant hDPP III (1.1 nM) was preincubated with quinazolinone-Schiff's bases (100 mM) first for 5 min at 27 � C and then for 5 min at 37 � C in 50 mM Tris-HCl buffer, pH 7.4.The enzymatic reaction was started with ARG 2 -2NA (40 lM) as a substrate, and after the 15 min incubation at 37 � C in a water bath, the reaction was stopped and the absorbance was measured using the spectrophotometric method described before ( � Spoljari� c et al., 2009).The stock solutions (8 mM) of quinazolinone-Schiff's bases were prepared daily in dimethyl sulfoxide and diluted with 50 mM Tris-HCl buffer, pH 7.4 before assay of enzymatic activity.Percentage enzyme inhibition (% inh.) was calculated by comparing the enzymatic activity without (control activity), and with inhibitor (inhibited activity) using the following formula: Compounds that showed strongest inhibition (% inh.� 90) were subjected to further IC 50 determination.The percentage of enzyme inhibition under at least 9 different concentration gradients (increasing compound concentrations from 0.05 to 100 lM) were imported into GraphPad Prism 5 software (GraphPad Prism 5.0.0.) to calculate IC 50 values of the compounds using non-linear regression analysis from the mean inhibitory value of three replicates.
Descriptor calculation for the optimized compounds was performed with Parameter Client (Virtual Computational Chemistry Laboratory, an electronic remote version of the Dragon program) (Tetko et al., 2005).Experimentally obtained values of hDPP III inhibition percentages were expressed as logarithmic values in order to normalize their distribution, and used as response values.The generation and validation of QSAR models was performed in QSARINS 2.2.4 (Gramatica et al., 2013).Constant and semi-constant descriptors, i.e., descriptors with a constant value for more than 85% of compounds and descriptors that were too intercorrelated (> 95%) were rejected by QSARINS.The final number of descriptors after the reduction was 503.A genetic algorithm (GA) was used to generate the best model, and the number of descriptors in the multiple linear regression (MLR) equation was limited to three.Concerning the limited number of compounds (28), the set was not split into the training and test, therefore the generated models are not predictive.
All generated models were validated by the internal cross-validation performed using the "leave-one-out" (LOO) and Y-scrambling method (Gramatica et al., 2013).Evaluation criteria taken into consideration were: coefficient of determination (R 2 ), adjusted coefficient of determination (R 2 adj ), cross-validated correlation coefficient (Q 2 LOO ), inter-correlation among descriptors (K xx ), the difference of the correlation among the descriptors and the descriptors plus the responses (DK), the standard deviation of regression (s), Fisher ratio (F), root-mean-square error (RMSE tr ); LOO crossvalidated root-mean-square error (RMSE cv ), concordance correlation coefficient (CCC tr ), LOO cross-validation concordance correlation coefficient (CCC cv ), mean absolute error of the training set (MAE tr ), mean absolute error of the internal validation set (MAE cv ).Y-randomization test was performed to check the robustness of models, providing R 2 Yscr , Q 2 Yscr and RMSE AV Yscr values (Gramatica, 2007).Williams plot was used to identify the possible outliers and compounds out of the warning leverage (h � ) in a model.The warning leverage is defined as 3p'/n (n is the number of training compounds, and p' the number of model adjustable parameters) (Eriksson et al., 2003).Compounds having values of standardized residuals higher than two standard deviation units were deemed as outliers.

Model setup
The initial model for docking and MD simulations was built using the semi-closed conformation of hDPP III obtained earlier (Tomi� c et al., 2012) by MD simulations of the structure available in the Protein Data Bank (PDB code: 3FVY, resolution 1.9 Å), since it has been proved that this is the most preferable enzyme form in water solution (Tomi� c et al., 2015).All ARG and LYS residues are positively charged (þ1) while all ASP and GLU residues are negatively charged (-1), as expected under the physiological conditions.To search for the best pose of the five most active compounds to the enzyme active site AutoDock Vina 1.1.2was used (Trott & Olson, 2010).The docking site for these five compounds on hDPP III was defined by cubical grid box with dimensions 80 � 80 � 80 Å 3 and the center placed on the Zn 2þ .Docking was done with the standard grid box spacing of 0.375 Å and 20 conformations were generated.The complex with the best AutoDock Vina docking score was chosen for the productive MD simulations.Parameterization of the complex structure was performed by the AMBERTools16 modules antechamber and tleap using General Amber Force Field (GAFF) (Wang et al., 2004) and ff14SB (Maier et al., 2015) force fields to parameterize the ligand and the protein, respectively.For the zinc cation, Zn 2þ , new hybrid bondednonbonded parameters were derived from our previous work (Tomi� c et al., 2019).The complex was dipped into the truncated octahedral box filled with TIP3P water molecules with a margin distance of 11 Å.To neutralize the system, added 24 Na þ ions were placed in the vicinity of charged amino acids at the protein surface.

Dynamics simulations
The resulting system, consisting of �73 000 atoms (�20 500 molecules of water), were simulated using periodic boundary conditions wherein the electrostatic interactions were calculated using the particle-mesh Ewald method.Prior to the MD simulations, the solvated complex was energy-minimized in three cycles with different restraints.In the first cycle 1 500 steps of minimization were performed, where the first 470 steps were of the steepest descent method, and the rest was the conjugate gradient.Both the protein atoms and the zinc cation were constrained using a harmonic potential of force constant 32 kcal/(mol Å 2 ), to equilibrate water molecules.In the second cycle 2 500 steps were performed and only the first 470 steps of steepest descent were used.The zinc cation and protein backbone were constrained with 32 kcal/ (mol Å 2 ).Finally, in the third cycle, the same number of minimization steps was as in the first cycle, and both protein backbone and zinc cation were constrained with 10 and 32 kcal/(mol Å 2 ), respectively.Afterward, the minimized system was heated from 0 to 300 K during 30 ps using a canonical ensemble (NVT), and then equilibrated 80 ps during which the initial constraints on the protein and the metal ion were used.This was followed by another equilibration stage of 100 ps, during which the initial constraints on the protein and the metal ion were removed and the water density was adjusted.The time step during the periods of heating and the water density adjustment was 1 fs.The equilibrated system was then subjected to 210 ns of the productive MD simulations at constant temperature and pressure (300 K and 1 atm) using the NPT ensemble, without any constraints.The temperature was held constant using Langevin dynamics with a collision frequency of 1 ps À 1 , and the Berendsen barostat respectively.Bonds involving hydrogen atoms were constrained using the SHAKE algorithm (Ryckaert et al., 1977).Simulations of the complex were performed within the AMBER16 software package (Case et al., 2016).The time step used for the productive MD simulations was set to 2 fs and the trajectory files were collected every 10 ps for the subsequent analysis.Trajectory analysis was performed by the CPPTRAJ module from the AmberTools16 program package and examined visually using Discovery Studio Visualizer, version 20.1.0.19295 (Dassault Syst� emes BIOVIA, 2019) and VMD 1.9.3 (Humphrey et al., 1996) software.
Comparing the compounds that exhibited the strongest inhibition against the enzyme (Table 1), it is clear that they all have the OH group at position 2 of the phenyl ring.The importance of the OH group at position 2 of the phenyl ring for the inhibitory effect of the tested compounds can be concluded from the comparison of 3 and 23 having OH at position 2 (% inh.90.3 and 97.4,respectively) with 21 without substituent on phenyl group (% inh.39.0), and with 16 (% inh.70.4) and 17 (% inh.66.1) having the OH group at positions 3 and 4. In addition to OH at position 2, the presence of the nitro group at position 5 of the phenyl ring proved to be essential for 23 which exhibited the strongest inhibitory potential (% inh.97.4) of all tested compounds.Compared to 3, the additional OH group at positions 3, 4, or 5 of the phenyl ring in compounds 8, 1, and 25 had no effect on increasing the inhibitory potential toward the enzyme.However, compounds 14 (% inh.81.7) and 19 (% inh.73.9) with OH groups at positions 3 and 4 and positions 3 and 5, respectively, also show that the presence of OH groups is important for the inhibitory effect.It is evident that besides the OH group, the presence of a substituent such as chlorine or methoxy group at position 2 of the phenyl ring also contributes to the inhibitory potential, as shown by compounds 4, 20, and 6 (Table 1).
The addition of bromine at positions 3 and 5 in compound 24 and bromine or chlorine at position 5 in compounds 5 and 22, respectively, causes a decrease in the inhibitory potential relative to compound 3. Also, the addition of the ethoxy group at position 3 and the benzyloxy group at position 4 in compounds 2 and 15, respectively, resulted in a reduction in their inhibitory effect on enzyme activity.The presence of bromine as the sole substituent on the phenyl ring showed a moderate inhibitory potential at position 4 in compound 26 (% inh.53.6) and a mild inhibitory potential at position 3 in compound 11 (% inh.29.9).Substitution of bromine with chlorine at position 3 did not affect inhibitory effect of compound 9. Compound 10 having a methoxy group at position 3 showed slightly better inhibition (% inh.38.1) compared to compound 18 (% inh.23.4) with an additional OH group at position 2, however, switching positions of these two substituents in compound 27 caused increased enzyme inhibition (% inh.41.6).Mild enzyme inhibition was shown in 13, which had a dimethylamino group at position 4 of the phenyl ring, and 12 with a phenyl-allylidene group.Compound 28 used as a precursor for synthesis of the screened quinazolinone Schiff bases showed a weak inhibitory potential (% inh.5.8) toward hDPP III.

QSAR analysis
The best multiple linear regression QSAR model (1) obtained for hDPP III inhibition is: where IC1 is information content index (neighborhood symmetry of 1 st order), GATS6m is a Geary autocorrelation of lag 6 weighted by mass, and EEig15d is the 15 th eigenvalue from edge adjacency matrix weighted by dipole moment.The variables in equation are listed in order of their relative importance by the standardized regression coefficient.Model satisfied the threshold values for fitting and internal validation criteria (Gramatica, 2007).The statistical parameters for obtained model are listed in To exclude the possibility that the obtained model is overfitted, the collinearity of the descriptors has been evaluated with a correlation matrix (Table 3), which showed that all correlation coefficients R were lower than 0.7.
Low collinearity between descriptors was also verified by the low value of K xx and DK being higher than 0.05 (Todeschini et al., 1999).The model satisfied fitting criteria: both coefficient and adjusted coefficient of determination (R 2 and R 2 adj ) are higher than 0.60; the concordance correlation coefficient of the training set (CCC tr ) is higher than 0.85; the root mean square error of the training set (RMSE tr ) and mean absolute error (MAE tr ) are close to zero (Table 2).Model also satisfied all of the following internal validation criteria (Gramatica, 2007): the cross-validated squared correlation coefficient Q 2 LOO needs to be higher than 0.6; the root mean square error determined through the cross-validated method (RMSE cv ) should be higher than the RMSE tr ; mean absolute error of the internal validation set (MAE cv ) close to zero.To check the robustness of the model and to eliminate chance correlation, Yrandomization test (Y-scrambling) was performed.Values of Y-scramble correlation coefficient (R 2 Yscr ) and Y-scramble cross-validation coefficients (Q 2 Yscr ) are both lower than 0.2 and R 2 Yscr is higher than Q 2 Yscr .Moreover, the root-mean-square error of Y-randomization (RMSE AV Yscr ) is higher than RMSE cv , so the chance correlation of model was excluded (Kiralj & Ferreira, 2009).Williams plot revealed that there are no compounds outside of warning leverage (h � ¼ 0.429) but there is only one outlier, compound 7 (Figure 3).
The first descriptor in model equation, IC1, is an information content index, which describes neighborhood symmetry of 1 st order.The information content descriptors evaluate and describe the effect of size, shape, and branching of molecules (Todeschini & Consonni, 2009).Therefore, a positive coefficient of IC1 means that compounds with relatively higher values of this descriptor tended to exhibit higher inhibition of hDPP III.Since higher values of this descriptors mean more complex 1 st order neighborhood, compounds with more diverse atom pairs may exhibit enhanced inhibition.The second descriptor, GATS6m, represents the Geary autocorrelation of lag 6 weighted by mass, and belongs to 2D autocorrelations.These descriptors reflect the contribution of pairs of atoms at the defined topological distance or spatial lag.The atoms represent the set of discrete points in space and the atomic property is the function evaluated at those points (Consonni & Todeschini, 2011).
GATS6m describes dependence of one atom on value of mass of other atoms at the topological distance 6. Positive regression coefficient of this descriptor in model suggests a favorable effect of heavier atoms at the 6 Å distance (Figure 4).The third descriptor in the equation is EEig15d, 15 th eigenvalue from edge adjacency matrix, weighted by dipole moment.It belongs to the edge adjacency indices, derived from the edge adjacency matrix from a hydrogen-depleted molecular graph of molecules.Descriptor EEig15d is related to the molecular polarity, describing the distribution of electronic charge in molecule (Estrada & Ramirez, 1996;Todeschini & Consonni, 2009).Since the active compounds have negative values of this descriptor, and this descriptor has a negative regression coefficient in the model equation, it means that the presence of atoms higher polarity, such as oxygen atoms from hydroxyl groups, is related to the enhanced inhibition.The importance of hydroxyl groups has already been established in earlier studies (Agi� c et al., 2017;2021).Although the obtained model is not predictive, it can provide some guidelines for the future compounds.According to the QSAR model, following structural features could improve the inhibition effect of quinazolinone Schiff bases against hDPP III: a. substitution of hydrogen atoms from the quinazolinone core with various atomic groups that contribute to the molecular polarity (e.g., hydroxyl groups) b. position these substituents at the 6 Å distance from the nitrogen and oxygen atoms of quinazolinone core.
The inspection of Williams plot (Figure 3) detected compound 7 as an outlier (compound with standardized residuals greater than 2 standard deviation units).The three methoxy groups of the phenyl ring decrease the overall neighborhood diversity and polarity of this compound.Also, these groups are sterically-crowded because of their vicinity, which is probably unfavorable for the hDPP III inhibition.

Molecular docking
The study was continued using molecular docking to obtain information on the possible mechanism of the five most active compounds (1, 3, 8, 23 and 25) in terms of binding affinity and necessary intermolecular bonds formation responsible for observed enzyme inhibition.AutoDock Vina was used to find for the best pose of these compounds in the semi-closed form of hDPP III.The best pose predicted for compounds 23, 8, 25, 1, and 3 were located deep in the enzyme binding pocket near the lower b-sheet (residues 389-393) (Figure 5A) with docking scores of À 8.9, À 8.1, À 8.0, À 7.9, and À 7.2 kcal mol À 1 , respectively.The minimum distance between the catalytic Zn 2þ ion and the screened compounds (C6 atom of the quinazolinone core) is � 4 Å and is less than the distance between the zinc cation and the nearest ligand atom (� 7 Å) of docked inhibitor molecules 3-benzoyl-7-hydroxy-2H-chromen-2-one (Agi� c et al., 2021) and luteolin (Agi� c et al., 2017), suggesting that these compounds are positioned closer to the enzyme active site relative to these previously investigated inhibitors.Also, the position of the quinazolinone core of the docked compounds is close to the substrate position in the hDPP III active site (Tomi� c & Tomi� c, 2014;Zhang et al., 2021) which is accommodated similarly to the tynorphin in the peptide binding cleft (Bezerra et al., 2012).
According to the results given in Table 4 most compounds interact with a larger number of amino acids residues that are constituents of hDPP III S1 and S2 substrate binding subsites compared to the number of residues that are part of the S1', S2', and S3' subsites (Bezerra et al., 2012).
Detailed analysis of the intermolecular bond formation for tested complexes obtained by Discovery studio software showed that the orientation of all docked compounds is conserved through interactions with the eleven identical hDPP III residues (Figure 5A).Most of these are van der Waals interactions between phenyl ring of the tested compounds and residues SER108, GLU316, SER317, ARG319, GLU327, and GLU329, as well as interactions between quinazolinone core and TYR318, PRO387, GLY389, ILE390, and ASN391 (Figure 5B-F).The exception is in the complex with compound 1, where ARG319 formed hydrogen bond (H), GLU329 formed p-cation, and GLU387 and ILE390 formed p-alkyl interactions (Figure 5B right); compound 3 where GLU316, ARG319 and PRO387 formed H bond, pi-donor H bond, and p-alkyl interaction, respectively (Figure 5C right); compound 8 where GLU329 formed H bond while PRO387 and ILE390 formed p-alkyl interactions (Figure 5D right); compound 23 where SER317, TYR318 and GLU329 formed H bonds, and PRO387 and ILE390 with p-alkyl interactions (Figure 5E right); compound 25 where GLU327 formed H bond (Figure 5F right).In addition to the latter eleven residues, conformational stability of 1, 8, and 23 were associated with six binding residues, while in 3 and 25 with three and four residues, respectively (Figure 5B-F left).Among them, notable are van der Waals interactions of important catalytic residues: HIS568 with compound 25 and GLU451 with 1, 8, and 23.The results of the above analysis of intermolecular bond are in accordance with the part of discussion in Chapter 3.1.,that the presence of the OH group at position 2 of the phenyl ring of the tested compounds has the effect of increasing their inhibitory potential.Indeed, this OH group in all five most active compounds participates in hydrogen bonding with hDPP III residues as illustrated in Figure 5B-F right.

MD Simulations
It is well known that molecular dynamics (MD) simulation provide powerful tools to study the dynamic molecular mechanisms of ligand (inhibitor) -receptor (enzyme) complexes (Chen et al., 2020;Shi et al., 2020).Therefore, MD simulations were carried out on the best docked pose of the most potent inhibitor (compound 23) in complex with semiclosed form of hDPP III to further explore the ligand-receptor interactions.Productive 210 ns-long MD simulation of the complex was performed using the programme package AMBER16, and intermolecular interactions, including native contacts and hydrogen bonding were analysed throughout.Root mean square deviation (RMSD) was used to analyse dynamic behavior and protein stability during simulation.The optimized structure of the complex was used to describe types of interactions in more detail.

RMSD of the protein backbone atoms
Figure 6 presents an analysis for the RMSD as a function of simulation time of the protein backbone atoms for free hDPP III and in complex with 23.According to the RMSD profiles, the complexed protein fluctuates less in the first 100 ns of MD simulations compared to the free protein.The average RMSD value calculated during the whole simulation time for complex and free enzyme is 1.43 ± 0.17 Å and 1.61 ± 0.22 Å, respectively.Stabilization of complexed protein is slightly better from the 60th to the 100th ns (average RMSD ± SD 1.39 ± 0.08 Å) and in the last 50 ns of MD simulations (average RMSD ± SD 1.57 ± 0.09 Å).

Native contacts
Furthermore, we studied the interactions between receptor (protein) and ligand (compound 23) by calculating the relative occupancy of native contacts during the simulation time.
The relative occupancy of native contacts is the sum of fractions of native contacts during the entire trajectory for each residue pair relative to the total number of native contacts involved with that pair.Native contacts were defined as a distance between the atoms of receptor residues and ligand atoms within a distance cutoff of 5 Å.The nativecontacts command in the CPPTRAJ module was used to calculate fractions of native contacts during simulations.that form native contacts in the complex with compound 23 (Figure 7) with the enzyme residues that form intermolecular bonds in the best docked pose of 23 (Figure 5E), it can be seen that the largest number of identical residues (GLU316, SER317, TYR318, ARG319, ASP320, GLU327, GLU329, and SER384) converged with the phenyl ring and only one residue (PRO387) with the quinazolinone core of compound 23.The reason for the weak convergence of protein residues and the quinazolinone core obtained by molecular docking and MD simulations is the movement of the quinazolinone core during system equilibration.Namely, in the first step of equilibration, the quinazolinone core of compound 23 moves from its starting position, and in the second step, it assumes an orientation that remains preserved until the system stabilization (Figure S2) and throughout the remaining time of the MD simulation (Figure 8).Of the 14 residues mentioned above that form native contacts with compound 23, seven of them (PRO387, GLU316, SER317, TYR318, ARG319, GLU327,    and GLU329) form intermolecular bonds with all other tested compounds (Figure 5A-F).A possible explanation for the differences in IC 50 observed between compounds with lower IC 50 and compounds 1, 3 and 8 is an additional intermolecular bond of compounds 23 and 25 with PHE109, which was subsequently shown to be relevant protein residue involved in native contact formation with 23.

Types of intermolecular interactions
To perform a detailed analysis of different types of intermolecular interactions, the extracted structure of hDPP III in complex with 23 obtained from the trajectory after 210 ns of MD simulations was optimized and used as a representative.
As can be seen from Figure 9 and Figure S2, more interactions are formed between the phenyl ring of compound 23 and hDPP III residues relative to the quinazolinone core.The phenyl ring of 23 formed van der Waals interactions with GLU316, SER317, GLU327, and SER384 whereas the quinazolinone core formed this type of interaction only with MET569.Furthermore, the methyl group of the quinazolinone core forms an alkyl bond with PRO387, while the core itself forms three p-type interactions (Figure S3).The results of trajectory H-bonds analysis calculated by hbond command in CPPTRAJ module showed there were four H-bonds formed during the MD process with the occupation time >10% (Table S2 and Figure 10).After system stabilization (last 110 ns of the simulation time) the first H-bond is formed by the OE2 atom of GLU329 and hydroxyl group on the phenyl ring of compound 23 (Figure 10 A) with an occupation time of 83%.The second and third H-bond was formed between the H-N of ARG319 and the O2 and O3 atoms of the nitro group of 23 with the occupation time of 11% and 12%, respectively (Figure 10 C and D).These results are in accordance with the experimental analyses which have shown that the presence of the OH group at position 2 and the nitro group at position 5 of the phenyl ring is crucial to achieving a substantial inhibitory potential of the test compounds.Only one H bond with the occupation time of 68% was formed during the simulations by the carbonyl oxygen of the quinazolinone core and H-N of ILE386 (Figure 10 B).The involvement of GLU329 has been reported in the formation of a H-bond with the OH group at C7 of coumarin derivative (Agi� c et al., 2021) suggesting that this residue as well as the carbonyl oxygen and OH group plays an important role in stabilizing quinazolinone core and the phenyl ring of 23 in their binding to hDPP III.

Conclusions
In this paper, we report for the first time the inhibition profiles of quinazolinone-Schiff's bases against the hDPP III combining an in vitro experiment with an in silico molecular modeling study.Most of the 28 tested compounds exhibited inhibitory activity at the 100 mM concentration.QSAR revealed that the higher diversity of neighboring atom pairs, distribution of heavier atoms at the 6 Å distance from each other, and the presence of higher polarity atoms contribute to the enhanced enzyme inhibition.The best docking pose predicts that the five most active compounds bind to the interdomain cleft, near the enzyme active site, while molecular dynamics simulations of the best inhibitor 23 (IC 50 ¼0.96mM) indicate that the hydroxyl group on phenyl ring and carbonyl group on quinazolinone core of 23, as well as interactions with GLU329 and ILE386, PHE109, ARG319, and GLU327 are important for the inhibitor binding and stabilization to the hDPP III active site.

Figure 1 .
Figure 1.Different forms of human DPP III structures: (left) ,,open" form (PDB code: 3FVY), (middle) ,,semi-closed" form, and (right) ,,closed" form in complex with pentapeptide tynorphin (PDB code: 3T6B).Two domains of the DPP III protein are shown in cyan and green, the zinc cation is represented as a purple sphere, while the red stick model shows tynorphin.Color figure can be viewed in the online issue.

Figure 2 .
Figure 2. Graphs used for IC 50 value determinations of compounds 1, 3, 8, 23 and 25 against hDPP III.Dose-response data points represent the mean value ± SD of three replicates.Color figure can be viewed in the online issue.

Figure 3 .
Figure 3. Williams plot (plot of standardized residuals vs. leverages (h) for each compound) of applicability domain of the Model 1.The warning leverage (h � ) is defined as 3p'/n (n is the number of training compounds, and p' the number of model adjustable parameters).Color figure can be viewed in the online issue.

Figure 4 .
Figure 4. Interatomic distances in Å (colored black) of atom pairs shown on the optimized structure of the experimentaly most active compound (23).Color figure can be viewed in the online issue.

Figure 7
Figure7and FigureS1illustrates the protein residues that were involved in forming of native contacts with relative occupancy above 50% during last 110 ns of simulation time, i.e. after system stabilization.Most native contacts (11) are formed with relative occupancy from 50% to 80% while the remaining three are involved in forming native contacts with

Figure 5 .
Figure 5. Docked pose of the five most active compounds at the hDPP III binding site.A Surface representation of the protein showing the binding site of the same enzyme residues involved in the interactions with all docked compounds.The lower b sheet is colored in red, and zinc cation is represented as a purple sphere.B-F (Left) 3D representations of the rest of the interacting enzyme residues, not listed on A (for clarity), and 2D diagram (right) of all enzyme interactions with compounds 1 (B), 3 (C), 8 (D), 23 (E) and 25 (F), colored brown, cyan, black, orange and pink, respectively.Color figure can be viewed in the online issue.

Figure 6 .
Figure 6.RMSD of the hDPP III backbone atoms obtained during 210 ns of MD simulation.Color figure can be viewed in the online issue.

Figure 7 .
Figure 7. Hydrophobic surfaces of protein residues involved in native contacts formation with 23 (orange) during last 110 ns of MD simulations.The amino acid interacting residues are shown in stick representation.Color figure can be viewed in the online issue.

Figure 8 .
Figure 8. Overlay of the hDPP III-compound 23 complexes obtained after 110 ns (cyan) and 210 ns (orange).Color figure can be viewed in the online issue.
The interplanar distances and angles between quinazolinone core of 23 and the phenyl ring of PHE109, imidazole ring of HIS568 as well as distances and angle between quinazolinone core of 23 and the amide group of ILE386 (measured between intermolecular atoms as shown in Figures S4, S5 and S6, respectively) indicate p-p stacking, p-p T-shape and amide-p stacking interactions in the most time of MD simulations.According to the 2D depiction in FigureS2the OH group at position 2 of the phenyl ring of compound 23 forms an H-bond with GLU329, the oxygen atom of the quinazolinone core forms an H-bond with ILE386, and oxygen from the nitro group forms H bonds with TYR318 and ARG319.Additionally, another oxygen atom of the nitro group forms a carbon-H bond with SER108.

Figure 9 .
Figure 9. 3D interactions of compound 23 with amino acid residues of hDPP III obtained after 210 ns of MD simulations as presented in the 2D scheme (Figure S2).Color figure can be viewed in the online issue.

Figure 10 .
Figure 10.Formation of the hydrogen bonds between compound 23 (LIG) and GLU329, ILE386 and ARG319 during 210 ns of MD simulations.

Table 2 .
The values of the descriptors included in the model are given in the Supplementary Materials Table S1.The logarithmic values of hDPP III inhibition percentage, experimentally obtained and calculated by the model Eq.(2), are given in Table 1.Logarithmic values calculated by model equation that are slightly higher than 2 indicate that these compounds could be potent inhibitors at lower concentrations.DK (global correlation among descriptors); RMSE tr (root mean square error of the training set); MAE tr (mean absolute error of the training set); CCC tr (concordance correlation coefficient of the training set); Q 2 LOO (cross validated explained variance); RMSE cv (root mean square error of the training set determined through the cross validated method); MAE cv (mean absolute error of the internal validation set); CCC cv (concordance correlation coefficient test set using cross validation); R 2

Table 1 .
Structures of analysed compounds, values of experimentally determined inhibition of hDPP III (at 100 mM concentration of compounds) and calculated logarithmic values of the % inhibition of hDPP III.

Table 2 .
The statistical parameters for the best QSAR model.(standard deviation of regression); F (Fisher ratio); K xx (multivariate correlation index); DK (global correlation among descriptors); RMSE tr (root mean square error of the training set); MAE tr (mean absolute error of the training set); CCC tr (concordance correlation coefficient of the training set); Q 2 LOO (cross validated explained variance); RMSE cv (root mean square error of the training set determined through the cross validated method); MAE cv (mean absolute error of the internal validation set); CCC cv (concordance correlation coefficient test set using cross validation); R 2 Yscr (Yscramble correlation coefficient); Q 2 Yscr (Y-scramble cross-validation coefficient); RMSE AV Yscr (root mean square error of Y randomization).

Table 3 .
Correlation matrix for the descriptors included in model 1.

Table 4 .
Substrate binding subsites of hDPP III in complex with docked compounds.