Pharmacophore generation and atom based 3D-QSAR of quinoline derivatives as selective phosphodiesterase 4B inhibitors

Abstract

Phosphodiesterase 4B (PDE4B) hydrolyses cyclic adenosine monophosphate (cAMP) and thus regulates its intracellular levels. The enzyme has been proposed as a potential drug target against diseases like inflammation and chronic obstructive pulmonary disease. But use of current PDE4B inhibitors is limited due to dose-dependent nausea and vomiting. Adverse effects associated with current PDE4B inhibitors are possibly results of PDE4D inhibition, a highly similar homolog of PDE4B. Here we considered quinoline analogs and applied ligand-based pharmacophore and atom based 3D-QSAR modeling with structure-based docking and ADME approach. A 5-point pharmacophore model was developed and used to derive a predictive 3D-QSAR model for the studied dataset. The obtained r² and q² values were 0.96 and 0.91, respectively. The result suggested that the generated 3D-QSAR model is reliable and can be considered for PDE4B activity prediction. Further, a pharmacophore model was employed for virtual screening to identify potent PDE4B inhibitors. The selective ligands for PDE4B were identified through docking and prime binding energy analysis of ligands in both PDE4B and PDE4D. ADME analysis was performed to confirm the drug ability of the selective ligand. To validate docking results, a molecular dynamics simulation was performed for PDE4B complexed with a top scoring ligand, AQ-390/42425549. AQ-390/42425549-PDE4B interactions reported in MD analysis were consistent with the docking results. All the hit molecules were procured and biologically evaluated for percentage inhibition of PDE4B and PDE4D in in vitro enzymatic assays. Among the total of thirteen molecules that were active against PDE4B, ten were selective with little PDE4D inhibition.

---

1. Introduction

Phosphodiesterases are a superfamily of enzymes catalysing the hydrolysis of 3′,5′-cyclic adenosine monophosphate (cAMP) and 3′,5′-cyclic guanosine monophosphate (cGMP) to the inactive 5′-AMP and 5′-GMP, respectively. cAMP serves as a secondary messenger inside the cell which transmits signals from the cell surface to a target.¹ Also, cAMP regulates immune responses by suppressing release of pro-inflammatory mediators e.g. tumor necrosis factor-α (TNF-α), interleukin-7 (IL-7), interferon-γ (IFN-γ)² and promotes the release of anti-inflammatory mediators e.g. IL-10.³ As central regulators of intracellular concentrations of cyclic nucleotides, PDEs are considered pharmacological targets for various disease therapies, such as for congestive heart failure,⁴ erectile dysfunction,⁵ inflammatory bowel diseases, and rheumatoid arthritis, chronic obstructive pulmonary disease (COPD) and asthma.^6,7 PDE4 is predominantly expressing enzyme in immune, inflammatory cells, neutrophils, T-cells, macrophages and studied in diverse spectrum of inflammatory responses both in vitro and in vivo.^8–11

The PDE4 enzyme family consists of 4 members (PDE4A, PDE4B, PDE4C and PDE4D) that are majorly expressed in neutrophils, monocytes, central nervous system (CNS) and smooth muscles of the lung.^12–15 Rolipram, non-selective PDE4 inhibitor is reported with severe side effects like nausea and emesis. Second generation PDE4 inhibitors like roflumilast came out with lesser side effects but could not overcome narrow therapeutic window.^16,17 Manning et al. proposed that the side effects of PDE4 inhibitors is due to the non-selective inhibition of PDE4 isoforms¹⁸ and hence selective inhibition of PDE4 subtypes would be beneficial to improve therapeutic effect to side effect ratio. Previously, several research groups have studied that deletion of PDE4B gene resulted in significant downfall of lipopolysaccharide (LPS)-induced TNF-α production in circulating monocytes and peritoneal macrophages.^19,20 Further, the behavioral correlation of emesis in mice by deleting PDE4D encoding gene concluded that inhibition of PDE4D and not PDE4B is responsible for emetic effect of non-selective PDE inhibitors.²¹ Therefore, selective PDE4B inhibitors can provide anti-inflammatory efficacy without side-effect. However, conserved active site residues of PDE4B and PDE4D make it difficult to design selective PDE4B inhibitors.

In the line of identifying selective PDE4B inhibitors, using the computational techniques, Guariento et al. attempted a ligand-based comparative molecular fields analysis (CoMFA) study for selective PDE4B inhibitors.²² Tripuraneni et al. worked on pharmacophore and atom based 3D-QSAR studies of pyrozolo[1,5-a]pyridine/4,4-dimethylpyrazolone analogues followed by molecular docking and molecular dynamics to design novel molecules with better PDE4 inhibitory activity.²³ Dong et al. applied Hopfinger's receptor-dependent QSAR (RD-QSAR) formalism and the Lukacova–Balaz scheme and initiated studies on structure-based multimode QSAR (SBMM QSAR) method, where the structure of the target protein is used in characterizing the ligand, and the multimode issue of ligand binding is systematically treated with a modified Lukacova–Balaz scheme.²⁴ Yang et al. conducted CoMFA studies on 315,6-dihydro-(9H)-pyrazolo-[4,3-c]-1,2,4-triazolo-[4,3-α]-pyridine analogs with variable inhibition of PDE4 to develop models for establishing three-dimensional quantitative structure–activity relationships (3D-QSAR). Comparative molecular field analysis (CoMFA) was conducted on the group of 5,6-dihydro-(9H)-pyrazolo-[4,3-c]-1,2,4-triazolo-[4,3-α]-pyridine analogs to determine the structural requirements for potency in inhibiting PDE4.²⁵ Srivani et al. developed 3D-QSAR with 44 derivatives of triarylethanes to drive structural requirement for PDE4 enzyme inhibition.²⁶ Chakraborti et al. analysed a series of 29 thieno[3,2-d]pyrimidines with PDE4 activity and subjected 3D-QSAR studies using CoMFA and CoMSIA.²⁷ Paola Fossa et al. produced rational pharmacophoric model of the PDE4 enzyme active site using DISCO approach with structurally diverse compounds.²⁸

Here, pharmacophore and 3D-QSAR modeling techniques were used to identify potent and selective PDE4B inhibitors.^29,30

The basis of such techniques is that the compounds interacting with the same target could share similar structural or physicochemical properties. A series of quinolines reported as selective PDE4B inhibitors, were used in the present study for developing pharmacophore and 3D-QSAR model.^31,32 The quinoline framework was considered as it has emerged as template for the design and identification of novel anti-inflammatory agents.³³ The developed models gave useful information of lead optimization for future rational design of PDE4B inhibitors and will be helpful in development of selective inhibitors. The best quantitative model was used as a 3D search query for screening the Specs database to identify new inhibitors of PDE4B. Once identified, the candidate compounds were subsequently subjected to filtration by molecular docking. The selectivity of filtered ligands was optimised by docking in both PDE4B and PDE4D (Fig. 1). ADME analyses were performed to confirm the drug ability of selective ligand.

Fig. 1 The overall work flow progressed during identification of selective PDE4B inhibitors.

2. Computational methods

2.1 Software and workstation

All computations were performed on Fujitsu linux workstation (xeon quad-core E3-1220 processor). Phase 3.9, Glide 6.3, LigPrep 3.0, Prime 3.6, Impact 6.3, and QikProp 4.0 of Maestro 9.8 (Schrödinger, LLC, 2014) were used to generate pharmacophore hypothesis, docking studies and ADME studies. Molecular dynamics (MD) study was done with AMBER 14 on 320 processor SUN Microsystems cluster. The analyses of MD simulations were carried out by CPPTRAJ module of AMBER 14 and Visual Molecular Dynamics 1.9.1 tool (VMD). The crystal structures of human PDE4B and PDE4D were retrieved from RCSB Protein Data Bank (PDB) and were utilized for molecular modeling studies.

2.2 Dataset

The pharmacophore and 3D-QSAR model were generated based on 49 published quinoline analogs that show wide range of biological activity and were tested using similar biological method (Table 1). The molecular structures were sketched and built with Maestro 9.8 (Schrödinger, LLC, New York software package). L2 and W1 represent identical molecule, which was tested in two different studies and was annotated with different names ((L2-molecule-4)³² (W1-molecule-3)³¹). In QSAR studies, appropriate conformation (lowest energy) of the compound is required for accurate calculation of 3D descriptors. The geometry optimization was carried out using the semi-empirical OPLS_2005 force field. All the molecules were aligned to the bioactive conformation of one of the ligands submitted as co-crystal ligand (PDB ID 3GWT, co-crystallized ligand-066, resolution-1.75 Å). The pharmacophore models were generated with the “Develop Pharmacophore Model” module of Phase 3.9.³⁴ PHASE incorporates structure cleaning step (LigPrep3.0), which attaches hydrogens, converts 2D structures to 3D, generates stereoisomers and neutralizes charged structures or determines the most probable ionization state at a user-defined pH. The prepared ligands were subjected to conformational analysis using Monte-Carlo Multiple Minimum method implemented in the Schrodinger software. The ligands were assigned as active and inactive by giving an appropriate activity threshold values to develop statistically robust models. The training set (also the test set) were divided in three categories: poor, modest and highly selective PDE4B inhibitors. We categorised all available PDE4B inhibitors into three major categories: active (>5.5), inactive (<4.7) and modest activity (4.7–5.5). All the molecules were divided into training and test set maintaining the similar biological activity diversity in both the sets for development and validation of QSAR model.

Table 1 Structures of all quinoline derivatives with their reported and observed biological activity (pIC50). The L-series belongs to Lunniss et al. and the W-series belongs to Woodrow et al

Ligand name	Structure	Activity pIC50	Predicted activity	Ligand name	Structure	Activity pIC50	Predicted activity
L_1		7.0	7.3	W_1		8.4	8.3
L_2		8.4	8.3	W_2		11.1	11.5
L_3		8.8	8.9	W_3		8.8	8.9
L_4		5.5	5.7	W_4		8.5	8.7
L_5		4.7		W_5		6.7	7.1
L_6		6.3	6.9	W_6		6.8	7.3
L-7		6.6	6.8	W_7		7.5	7.7
L-8		4.9	5.5	W_8		7.6	7.9
L-9		5.5	5.6	W_9		9.1	9.0
L-10		4.7	4.8	W_10		8.7	8.9
L-11		4.8	4.9	W_11		8.0	8.3
L-12		5.1	4.8	W_12		8.5	8.4
L-13		4.5	4.5	W_13		8.2	8.2
L-14		4.5	4.3	W_14		7.4	7.6
L-15		7.8	7.9	W_15		7.4	8.0
L-16		7.7	7.9	W_16		8.6	8.9
L-17		8.6	7.9	W_17		7.9	7.9
L-18		7.0	6.9	W_18		9.9	9.7
L-19		7.1	7.2	W_19		9.2	9.3
L-20		5.3	5.9	W_20		11.0	11.1
L-21		7.5	7.8	W_21		8.8	9.2
L-22		8.2	7.9	W_22		10.0	9.7
L-23		6.0	5.9	W_23		10.7	10.5
L-24		9.5	8.8
L-25		8.3	8.3
L-26		9.4	9.2

---

2.3 Pharmacophore modeling

The chemical features of all ligands were defined by four pharmacophoric features: four H-bond acceptor (A), two H bond donor (D), one hydrophobic group (H), and three aromatic rings (R). An active analog approach was used to identify common pharmacophore hypotheses (CPHs), in which common pharmacophores were culled from the conformations of the set of active ligands using a tree-based partitioning technique that groups together similar pharmacophores according to their inter site distances.

The resulting pharmacophores were scored and ranked. The scoring was done to identify the best hypothesis, which provided an overall ranking of all the hypotheses. The scoring algorithm included the contributions from the alignment of site points and vectors, volume overlap, selectivity, number of ligands matched, relative conformational energy, and activity. After the careful analyses of the scores and alignment of the active ligands to the generated hypothesis, a best pharmacophore hypothesis was selected.

For validation of the developed hypothesis, a common dataset of 27 active^35–38 and 260 inactive molecules was used. Inactive molecules were chosen such that they consists of similar physiochemical descriptors (molecular weight, number of rotational bonds, hydrogen bond donor count, hydrogen bond acceptor count and octanol–water partition coefficient) to active molecules but deprived of any of the chemical descriptors of the active ligands.

Various statistical parameters such as accuracy, precision, sensitivity, specificity and enrichment factor (E value) were calculated for each hypothesis.

Furthermore, for the analysis of results, E value score was calculated using the following formula

where D, A, Ht and TP represent the total number of compounds of the database, total number of actives, total number of hits and total number of true actives in the hits, respectively.

2.4 3D-QSAR modeling

QSAR modeling was carried out using Phase module from Schrodinger. An atom based alignment was used to develop QSAR model, which is more useful in explaining the structure–activity relationship in comparison to the pharmacophore based alignment. In atom-based QSAR, a molecule is treated as a set of overlapping van der Waals' spheres. Each atom (and hence each sphere) is placed into one of six categories according to a simple set of rules: hydrogens attached to polar atoms are classified as hydrogen bond donors (D); carbons, halogens, and C–H hydrogens are classified as hydrophobic/non-polar (H); atoms with an explicit negative ionic charge are classified as negative ionic (N); atoms with an explicit positive ionic charge are classified as positive ionic (P); non-ionic nitrogen and oxygen are classified as electron-withdrawing (W); and all other types of atoms are classified as miscellaneous (X). The data set of 49 molecules was divided into training set and test set in a random manner.

For the purpose of QSAR development, van der Waals' models of the aligned training set molecules were placed in a regular grid of cubes, with each cube allotted zero or more ‘bits’ to account for different types of atoms in the training set that occupy the cube. This representation gives rise to binary-valued occupation patterns that can be used as independent variables to create partial least-squares (PLS) QSAR models. QSAR models were generated for the selected hypothesis using the 40-member training set using a grid spacing of 1.0 Å. The best QSAR model was validated by predicting activities of the 8 test set compounds and 1 outlier. A four component (PLS factor) model with good statistics was obtained for the dataset whereas the maximum number of PLS factors in each model can be 1/5 of the total number of training set molecules. Further increase in the number of PLS factors did not improve the model statistics or predictive ability.

For further validation of QSAR model, external r² prediction was carried out with the data set of thirty two reported PDE4B inhibitors of diverse scaffolds^35–39 having activity ranging from 5.0 to 10 pIC₅₀. The reported IC₅₀ values were converted to pIC₅₀ values. The dataset was prepared in Ligprep and then imported to search for matches panel of Phase keeping atom based QSAR model as query. The obtained hits were analysed for their predicted pIC₅₀ values.

2.5 Virtual screening

The novel and potential leads can be predicted from large chemical databases using virtual screening (VS) techniques.⁴⁰ The 3D pharmacophore model represents the chemical functionalities responsible for the bioactivities, and hence the validated 3D queries. Pharmacophore hypothesis developed in previous step was used for VS of SPECS database containing 203152 molecules. The hypothesis was used to search a 3D specs database for structures that match the pharmacophoric features of the model. Virtual screening was carried out using search for matches panel of Phase that uses the pharmacophore to efficiently search the database of fixed conformers for pharmacophore matches. The molecules were initially submitted to Phase 3.9 to pre-generate conformers to fasten screening process. These pre-selected, prepared molecules were screened through Lipinski filter. The filtered molecules were next preceded for the virtual screening considering previously generated pharmacohore hypothesis as query. The number of conformers per rotatable bond was kept as 10 and 100 set as the maximum number. Relative energy window was kept by default at 10.0 kcal mol⁻¹. The cutoff for rotatable bond was 15. Minimum matching site points were set to 2 out of 5 and inter-site distance matching tolerance was set to 2 Å. Molecules which fit well with the pharmacophoric features of the selected hypothesis were retrieved as a hit and were docked into the PDE4B (PDB ID 3G45) and PDE4D (PDB ID 3G4G) structures (next section).

2.6 Molecular docking

Preparation of receptors and small molecules. Crystal structures of PDE4B (PDB ID 3G45, resolution-2.63 Å) and PDE4D (PDB ID 3G4G, resolution-2.3 Å) were prepared using Protein Preparation Wizard (Impact 6.3, Schrodinger). Although the selected structures were not the highest resolution structures yet they were considered due to additional N-terminal residues. PDE4B and PDE4D vary in their N-terminal region (Tyr274-4B::Phe196-4D) which can be targeted to prove the selectivity of the ligands for PDE4B.⁴¹ The proteins were pre-processed with bond order assignment, hydrogen addition, formal charges and metal treatment, and deletion of all water molecules. Hydrogen bonding network and orientation of Asn, Gln, and His residues were optimized based on hydrogen bond assignment using the exhaustive sampling option. The states of histidine (HIS, HIE, or HIP) were assigned after optimization. The proteins were minimized to an RMSD limit from starting structure of 0.3 Å using Impref module of Impact with OPLS_2005 force field. Missing side chains were structured using Prime 3.6.

Low-energy conformations of ligands that were used for docking program Glide were generated via Ligprep 3.0 with Epik 2.8 of Schrodinger. The optimised structures were produced based on force field OPLS_2005, with protonation states generated at target pH 7.0 ± 2.0. Thirty-two stereoisomers computed by retaining specified chirality and keeping low energy ring conformation 1 per ligand.

Receptor grid generation and docking. Docking is based on a grid represented by physical properties in the receptor volume that is searched for ligand–receptor interaction during docking process. Grid files were prepared with the “Receptor Grid Generation” panel of Glide. Docking grids were generated with the default settings in Glide using the co-crystal ligand to define the centre of the grid box (20 Å × 20 Å × 20 Å). The default parameters for van der Waals' radius scaling factor 1.0, partial charge cutoff 0.25 and charge scale factor 1.0 was used and no constraints were included during grid generation.

Docking was carried out to increase the reliability of the pharmacophore-based screening and to discriminate between the active–inactive and selective ligands for PDE4B. Prepared small molecules were docked into the protein structure using Glide 6.3 XP.⁴² The score function of Glide or Glide score is a modified and expanded version of ChemScore, was used for binding affinity prediction and ligand ranking. The docking can be on the level of either standard (SP) or extra precision (XP). The addition of large desolvation penalties to both ligand and protein, assignment of specific structural motifs that contribute significantly to binding affinity, and expanded sampling algorithm makes XP superior over SP. In this study, extra precision docking was applied, and the rest of the parameters like ‘dock flexibly’, ‘add epik state penalties to the docking score’ were kept default. The scaling factor was 0.8, and the partial charge cutoff was 0.15. The 3D complex structures of all hits docked in PDE4B and PDE4D were analysed for Glide score and H-bonding interactions. The molecules with considerable difference in docking score and H-bonding interaction in PDE4B and PDE4D were selected and evaluated for binding free energy calculation of docked complex using Prime-MM-GBSA-3.6, employing VSGB continuum dielectric model as solvent model.

2.7 ADME

The analysis of ADME (Absorption, Distribution, Metabolism and Excretion) properties are crucial determinants for the successful development of new drugs, and are imperative for rational drug design. The failure of which might lead to rejection of a molecule in the later stages of drug development process. All the shortlisted hits were further processed for ADME properties analysis using QikProp 4.0 tool of Schrodinger. QikProp is built using experimental details of 710 compounds including 500 drugs and heterocyclic compounds. It calculates physico-chemical properties like molecular weight, molecular volume, no. of H-bond donors, no. of H-bond acceptors, polar surface area, predicted octanol/water partition coefficient (QP log P_o/w) and violations related to Lipinski's ‘Rule of 5’ and Jorgensen's ‘Rule of 3’ to filter out compounds with clear cut undesirable properties.

2.8 Molecular dynamics simulation

Among newly identified hit molecules, the stability of protein–ligand interactions of docked AQ-390/42425549 (AQ) in PDE4B (PDB ID 3G45) was evaluated by molecular dynamics simulations. “ff14SB” force field in AMBER 14 (ref. 43) was used to generate topology and parameter files for the protein residues. The amber force field (GAFF) and the Parmchk module of Antechamber were used for generating missing parameters for the small molecule, Zn and Mg. All systems were simulated in a standard protocol of energy minimization which is 2500 steps of steepest descent followed by 1000 conjugate gradient steps. Thus minimized system was gradually heated from 0 to 298 K in 200 ps. Following equilibration, MD simulation was carried out for 10 ns with periodic boundary conditions in an isothermal–isobaric (NPT) ensemble at a temperature of 298 K with Berendsen temperature coupling and a constant pressure of 1 atm with isotropic molecule-based scaling. SHAKE algorithm was applied to fix all covalent bonds containing hydrogens and particle-mesh-Ewald (PME) method was used for long-range electrostatic interactions. The MD simulation analyses were carried out by CPPTRAJ module of AMBER 14 and VMD visualizer.⁴⁴ All MD simulations were performed on 320 processor SUN Microsystems cluster.

2.9 Biological evaluation

PDE4 activity was determined for subtypes B2 and D2 using respective assay kits from BPS Biosciences, San Diego, CA (catalogue no. 60343 and 60345, respectively). We screened twenty one compounds (procured from Specs) for their inhibitory activities. The assay kits work on the principle of fluorescence polarization. In brief, fluorescently active nucleotide monophosphoate which is generated by PDE4B2 (or PDE4D2) subtype binds to a binding agent followed by detection of this adduct using fluorescence polarimetry. All testing were performed as per the manufacturer's instructions. In brief, 25 μL of FAMCyclic-3′,5′-AMP (200 nM) was incubated with 20 μL of PDE4B2 (stock: 7.5 pg μL⁻¹) at 25 °C for 1 hour containing 30 μM of test compounds (stock solution dissolved in DMSO). After 1 hour incubation, 100 μL of binding agent was added to each well and the mixture is further incubated for 1 hour with gentle shaking.

The fluorescent polarization (FP) was measured in the Microtiter plate (Black, low binding NUNC Microtiter plate) using a fluorescence reader (TECAN Infinite F200 Pro Microplate Reader) set for excitation at wavelengths ranging from 475–495 nm and detection of emitted light ranging from 518–538 nm. A blank value was subtracted from all other values. The FP value for substrate control wells set as 0% activity while positive control wells set as 100% activity. The % activity for test compound wells were calculated with respect to the positive and substrate control wells. % inhibition was determined by subtracting individual % activity values from 100.

Likewise, PDE4D inhibition was separately calculated in similar experimental setup using 20 μL of PDE4D2 (2.5 pg μL⁻¹) subtype enzyme instead.

3. Results and discussions

3.1 Pharmacophore model

For the generation of pharmacophore model, we have considered 41 compounds having activity >5.5 against PDE4B as active as they contain important structural features crucial for binding to the receptors binding site. We used five as minimum sites and maximum sites to have optimum combination of sites or features common to the most active compounds. Twenty two common pharmacophore models were generated with different combination of variants; the results are illustrated in (ESI Table 1†). Among these pharmacophores, the models which showed the superior alignment with active compounds and good survival score in comparison to inactive molecules, was identified to generate QSAR model. The survival scoring function identifies the best candidate hypothesis from the generated models and provides an overall ranking of all the hypotheses. The scoring algorithm includes contributions from the alignment of site points and vectors, volume overlap, selectivity, number of ligands matched, relative conformational energy, and activity. However, these pharmacophore models should also discriminate between the active (most active) and inactive (less active) molecules (ESI Table 2†). The hypothesis that lacks critical site explaining ligand binding is not supposed to be reasonable for discriminating between actives and inactives. The models with maximum adjusted survival score selected for generating atom-based alignment of PDE4B inhibitors. Model AADRR.1108 has been selected because it produced good alignment of actives in comparison to inactives. The special arrangement of features along with their distance present in five-featured pharmacophore, AADRR, is shown in Fig. 2A and B. As depicted in the Fig. 2A, among the two ring aromatic features, features are mapped to the both aromatic rings of all 41 active inhibitors. The both hydrogen bond acceptor and negative ionic features are mapped to the sulphone groups and carboxyl group of carboxamide substituent. For generating an atom-based 3D QSAR hypothesis, we have used a dataset of 40 (training set) compounds having inhibitory activity against PDE4B. The model was validated using 8 (test set) compounds, which cover the same range of PDE4B inhibitory activity. But 1 compound (L-5) was found as outlier. The pre aligned molecules (Fig. 2C) were used for developing generated pharmacophore model AADRR. The alignment generated by the best pharmacophore model AADRR was used for QSAR model generation. Fig. 2D presents good alignment of the active ligands and scattered alignment of inactive ligands Fig. 2E to the developed pharmacophore model.

Fig. 2 (A) Pharmacophore model generated for PDE4B. Brown rings (R9 and R10): aromatic rings. Blue sphere (D7): hydrogen-bond donor and pink sphere (A3 and A5): hydrogen-bond acceptor. (B) Common pharmacophoric sites of active ligand with distance. All distances are in Å unit. (C) Alignment of all ligands (quinoline analogues) to the bioactive conformation of crystal ligand 066 of PDB ID 3GWT. (D) Alignment of all active ligands to the pharmacophore. (E). Alignment of all inactive ligands to the pharmacophore. Atom-based PDE4B 3D-QSAR models visualized for (F) positive effect (blue cubes) of electron withdrawing group/atom position (W+) mapped on the most active ligand, (G) negative effect (red cubes) electron withdrawing group/atom mapped position (W−) on the least active ligand, (H) negative effect (red cubes) and positive effect (blue cubes) of position hydrophobic groups (H+,−), (I) positive effect (blue cubes) of position H-bond donor (D+), (J) negative effect (red cubes) of position H-bond donor (D−), (K) negative effect (red cubes) of position of positive ionic groups (P−), (L) negative effect (red cubes) of position of negative ionic groups (N−).

A four-PLS factor model with good statistics and predictive ability was generated for the dataset (ESI Table 2†). The number of PLS factor included in model development is four as incremental increase in the statistical significance and predictivity was observed for each incremental increase in the incorporated PLS factors up to four. The model expressed r² value 0.95 exhibited by quinoline based derivatives, is close to one and signifying an acceptable agreement of fitting points on the regression line for the observed and PHASE-predicted activity that is shown in Fig. 3A and B and is summarized in Table 2.

Fig. 3 (A) Observed and predicted activities of test set compounds associated with PDE4B. (B) Observed and predicted activities of training set compounds associated with PDE4B.

Table 2 PLS statistical parameters of the selected 3D-QSAR model^a

ID	# factors	SD	r^2	F	Stability	RMSE	q^2	Pearson-R	P
a SD-standard deviation of the regression; r^2-squared value of r for the regression; F-variance ratio (large values of F indicate a more statistically significant regression); P-significance level of variance ratio (smaller values indicate a greater degree of confidence); RMSE-root-mean-square error; q^2-squared value of q for the predicted activities; Pearson-R-value for the correlation between the predicted and observed activity for the test set.
AADRR.1108	1	1.29	0.47	33	0.57	0.57	0.81	0.93	1.38 × 10^−6
	2	0.86	0.77	60.21	0.67	0.46	0.88	0.95	3.29 × 10^−12
	3	0.54	0.91	120.43	0.81	0.52	0.84	0.98	1.68 × 10^−18
	4	0.39	0.95	185.31	0.93	0.39	0.91	0.99	1.42 × 10^−22

---

Generated hypotheses were assessed using dataset of active (27) and inactive (260) ligands. Dataset from literature was employed to validate the generated pharmacophore hypotheses. To each models, different statistical parameters like accuracy, precision, sensitivity, and specificity of the best pharmacophore models were calculated. Furthermore, an E value of 1 was calculated for model (Table 3). From the overall validation results, we assured that hypotheses can differentiate between the actives and inactive molecules.

Table 3 The statistical parameters obtained from pharmacophore validation test

S. No.	Parameters	Hypothesis AADRR.1108
1	Total compounds in database (D)	287
2	Total number of actives in database (A)	27
3	Total hits (Ht)	210
4	Active hits (TP)	22
5	True negative (TN)	215
6	Enrichment factor or enhancement (E)	1.1
7	False negatives (FN = A − TP)	5.0
8	False positives (FP = Ht − TP)	188
9	Accuracy = (TP + TN)/(TP + TN + FP + FN)	0.6
10	Precision = TP/(TP + FP)	0.1
11	Sensitivity = TP/(TP + FN)	0.8
12	Specificity = TN/(TN + FP)	0.5

---

3.2 Atom based 3D-QSAR model

A 3D-QSAR study was carried on the analogues of quinoline based derivatives to investigate the effect of spatial arrangement of structural features on the PDE4B inhibition. The stability of statistically significant regression model is indicated by the large value of F (185.6) and high degree of confidence supported by the small value of the variance ratio (P). Further, the data used for model generation are best for the QSAR analysis is obvious from small values of standard deviation (0.38) of the regression and RMSE (0.39). The model was validated by cross-validated correlation coefficient (q²-0.91) that was obtained by leave one out or leave N out method. The reliability and robustness of q² statistical parameter is more than r² because it is obtained by external validation method by dividing the dataset into training and test set.⁴⁵

The external r² prediction with the reported PDE4B inhibitors was carried out to validate the QSAR model. The data set of thirty two molecules was taken as an input for searching matches against the developed model in Phase panel. The output window showed twenty molecules as hits. The molecules with the reported and predicted pIC₅₀ values are given in ESI Table 3.† The predicted activity and experimental activity are in statistical limits of confidence and hence the derived model can be accepted for further studies.

The atoms visualize 3D characteristics of the ligands (atoms or pharmacophores) as that contribute positively or negatively to activity. The QSAR model displays 3D characteristics as cubes that represent the model and color according to the sign of their coefficient values, which is indicated as blue for positive coefficients and red for negative coefficients. Positive coefficients indicate an increase and negative coefficients as decrease in activity. This might provide the information that which functional groups/atoms are desirable or undesirable at certain positions in a molecule. The blue cubes in 3D plots of the 3D pharmacophore regions indicates the ligand regions in which the specific feature is important for better activity, whereas the red cubes shows that particular structural feature or functional group that is not essential for the activity or likely the reason for decreased binding potency or affinity.

Visual analysis of Fig. 2F demonstrates that the presence of the blue cubes at the A3, A5, –NH linker and –NHCO group of dimethyl carboxamide and carboxamide substituent, methoxy group and sulphone group pointing out the positive potential of electron withdrawing characteristic of the molecules and is requisite for the activity at this particular place.

It is also apparent from the literature data that the replacement of electron withdrawing groups at the A3-sulphone group (W_3, W_4, W_1) to sulphide (W_6) or amide (W_7, W_8) or removal (W_5) lead to decline in activity. The phenyl sulphones are observed to be better fit to solvent filled pocket of PDE4B active site.

If 3-methoxyphenyl (W_3) was replaced with 3-ethylphenyl (W_11) or cyclohexyl (W_14) or phenyl (W_13), the activity was decreased as methoxy phenyl group was reported to be good fit to the surface.

3-Dimethylcarboxamide group (W_20) attached to phenyl ring, if replaced with 4-methoxy group (W_18) or H (W_3) or 4-t-butyl (W_19) or 4(3-furyl, W_21), all replacements resulted in leap in potency. This was observed that the meta-CONMe₂ develops significant van der Waals' contacts with Ser454 and desolvation effect by displacing trapped water molecules.

Likewise replacement of –NH linker (L_15) with –O atom (L_4), –S atom (L_5) and –NMe (L_23) lead to decrease in potency. The role of this –NH linker is reported in developing intramolecular H-bonding with 3-carboxamide substituent at quinoline ring.

Apart from this, 3-carboxamide substituent to quinoline ring was observed to be involved in H-bonding with Asp567, Gln615 and water molecules. Hence its replacement with –CN (L_6), –CH₂NH₂ (L_7), –COOH (L_8), –CONHMe (L_9), –H (L_12) or with oxadiazoles (L_10 and L_11) produced detrimental effect on the biological activity.

The negative impact of electron withdrawing groups was mapped on inactive ligand to produce remarkable picture. The red cubes in vicinity of –N atom of quinoline ring, –NH linker and sulphone group indicates that presence of electronegative atom in close to these group may produce leap in activity Fig. 2G.

Fig. 2H demonstrates the blue cubes and hence positive impact of hydrophobic groups on PDE4B inhibition and the results are supported by the evident potencies. It can be deduced from the figure and reported as well that hydrophobicity of quinoline ring (pi–pi cloud with Phe618), methyl substituent at quinoline ring (van der Waals contact with Met603), methoxy benzene ring (interaction with protein surface) and benzene ring attached to sulphone ring (interaction with solvent filled pocket) produced increase in biological activity. The substitution of hydrophobic groups near –NH linker to methoxy benzene is unacceptable (L_23) as shown by red cubes or may hinder the binding of the molecules to the receptor active site and will result in decreased PDE4B inhibition.

Fig. 2I illustrates that H-bond donor characteristics is necessary at D7 and –NH linker as consistent with the trend of activity of L_4-12. The red cubes (Fig. 2J) in proximity to D7, –NH linker (L_23) and A5 (L_9) demonstrate negative potential of extension of H-bond on those positions.

Further, Fig. 2K showed negative effect of positive ionic group in replacement with or near to methoxy substituent (red cubes) supported by the reported potencies of W_11, W_14-17. Also Fig. 2L shows the negative impact of negative ionic group near –N atom of quinoline ring and A5 (red cubes) consistent with the decreased biological activity of L_8 and L_6.

3.3 Virtual screening

Virtual screening of the molecular libraries is considered as one of the efficient approach to drug discovery. Pharmacophore based database searching is considered as a type of ligand-based virtual screening, which can efficiently be used to find novel and potential leads for further development. Initially Lipinski filter was applied and 203 152 molecules were retrieved which were next filtered. A potent pharmacophore model possesses the chemical functionalities responsible for bioactivities of potential drugs; therefore, it can be used to perform a database search by serving as a 3D query. The pharmacophore AADRR was used as a 3D structural query for retrieving potent molecules from the specs chemical database through. To avoid possible incorrect predictions, cut off limit of 15 rotatable bonds filters were applied. Altogether, 3D structural query retrieved 500 compounds.

3.4 Molecular docking and ADME analysis

Docking study was performed using Glide program⁴⁶ on previous hit molecules to identify the key interactions with PDE4B structure. Glide computes the binding energy (BE) in terms of the Glide docking score with respect to the docked ligands. Among all hits, molecules that have Glide score of less than −8 in PDE4B were only considered for further analyses.

The active site residues are mostly conserved between PDE4B and PDE4D. Beside catalytic Gln615 (Gln535 in PDE4D), other H-bonding residues in PDE4B active site like Tyr405 (Tyr325 in PDE4D), His406 (His326 in PDE4D), His410 (His330 in PDE4D), Asp447 (Asp367 in PDE4D), Asp564 (Asp484 in PDE4D), Asn567 (Asn487 in PDE4D) are also conserved between two closely related homologs. Likewise, the hydrophobic residues are also conserved.

Besides structurally resolved catalytic domain, full PDE4 proteins (all PDE4A-4D) also consist of upstream regulatory regions (UCR1-2). Recently, Burgin et al. showed that UCR2 closing over to PDE4B active site region prevent access to substrate ATP.⁴¹ Tyr274 of PDE4B-UCR2 is occupied by Phe196 of PDE4D-UCR2 which could make H-bond with the bound inhibitor and therefore provide selective PDE4B inhibition, which can be explored to design selective PDE4B inhibitors. Therefore, in our study we used structures with these additional resolved residues (PDE4B:PDE4D::3G45:3G4G) to perform docking and molecular dynamics study.

Most of the current PDE4B inhibitors are associated with PDE4D related adverse-effects. To minimise PDE4D related toxicity, molecules obtained with glide score < −8 were next docked in PDE4D structure (with resolved N-terminus residue, PDB ID 3G4G) and molecules with similar (or high) docking score in PDE4D structures were selected for further analyses. Molecules with minimum glide score difference of ≥1.5 (PDE4D < PDE4B) were only chosen as selective PDE4B inhibitors.

The ligands exploited the catalytic cavity of PDE4B to much extent and showed H-bond interaction with key residues e.g. Gln615, Tyr274, Asp564, Asn277, His406, Met519 and Ser614. Besides the H-bond, hydrophobic interaction (pi–pi) was established with Tyr405, Tyr274, His450, His406, F586 and F618. In addition to this, metal interaction with Mg²⁺ was also observed in some of the ligands. Whereas the intrinsic ligands cAMP, and reported inhibitors such as rolipram, roflumilast and 3-isobutyl-1-methylxanthine (IBMX) showed similar glide score, H-bond and hydrophobic ambience in both PDE4B and PDE4D, the screened ligands showed better selectivity profile for PDEB in terms of glide score, H-bond and pi–pi interaction. The ligands could develop H-bond in only three ligands with Tyr325 and Asn375. Phe538 interacted with ionic bond while metal interaction was observed with Mg²⁺. Pi–pi hydrophobic interaction was associated with Phe538, Phe506, His326 and Tyr325. Above all, while most of the ligands is found to possess either H-bond or pi–pi interaction with Tyr274 of N-terminal of PDE4B, none of the ligands interacted with N-terminal of PDE4D i.e. no interaction with Phe196 of PDE4D, showing selective profile of screened ligands.

54 ligands were finalized based on considerable difference in docking score of both PDE4B and PDE4D. 27 ligands were selected for considerable difference in Prime binding energy of PDE4B and PDE4D docked complexes.

The ADME properties of all the obtained ligands were assessed using Qikprop 4.0, and these ADME properties of best hits are listed in ESI Table 4.† The 21 of the hits showed optimum ADME properties while 6 were rejected. The ADME property of these hits makes them promising candidates for as PDE4B inhibitors. The selective inhibitors, their glide scores along with the hydrogen bonding of key residues are presented in Table 4 and ESI Table 5.† The binding pose of top scoring ligand is depicted in Fig. 4A.

Table 4 Glide score (kcal mol⁻¹) of the docked molecules in both PDE4B and PDE4D, and ΔG binding (kJ mol⁻¹) of PDE4B and PDE4D docked complexes

Molecule	G score		H-Bond residues		ΔG binding of complex
	PDE4B	PDE4D	PDE4B	PDE4D	PDE4B	PDE4D
AQ-390/42425549	−11.8	−8.7	Gln615, Tyr274	No	−96.5	−77.6
AG-690/11972161	−11.7	−8.7	Gln615	No	−95.5	−82.3
AO-476/43407280	−11.3	−9.3	Gln615	No	−88.6	−74.1
AG-690/36276051	−11.3	−7.8	Gln615	No	−76.4	−65.7
AO-022/43390834	−11.3	−4.7	Tyr274, Asp564	No	−27.6	−10.4
AG-205/36564043	−11.3	−8.5	His406	No	−30.5	−4.1
AG-690/36276042	−11.1	−8.2	Gln615	No	−66.8	−46.6
AG-690/36873050	−10.9	−8.8	Tyr405, Asn455	No	−57.3	−30.9
AJ-292/42152568	−10.6	−7	Asp564	Asn375	−72.5	−63.6
AF-399/42316263	−10.5	−5.7	Gln615	No	−64.9	−51.0
AK-968/37005156	−10.4	−8.4	Gln615	No	−83.9	−74.9
AO-022/43453692	−10.3	−8.1	Tyr405, Tyr274	Tyr325	−92.2	−78.3
AF-399/41980308	−10.2	−8	Gln615	No	−55.9	−43.0
AP-964/40915318	−10.0	−8.1	Tyr274, Met519, Asn277	No	−55.9	−18.7
AG-205/33161053	−10.1	−7.6	Gln615, Tyr274, Ser614	No	−34.8	−21.7
AT-057/43469096	−10	−5.3	Tyr274	No	−76.4	−39.1
AG-690/15437723	−9.7	−5.3	Gln615, Tyr274, His450	No	−59.2	−43.9
AK-968/12386394	−9.7	−5.9	Gln615	No	−69.5	−57.8
AF-399/15335138	−9.6	−7.4	Tyr274	No	−64.3	−50.8
AG-690/10252051	−9.4	−7.7	Gln615	No	−96.6	−85.5
AP-124/43383688	−9.2	−7.5	Tyr274	No	−57.7	−41.3
Rolipram	−9.2	−9.3	Gln615	No	−81.2	−83.4
IBMX	−7.3	−7.7	Gln615	Gln535	−51.5	−46.8
Roflumilast	−9.7	−9.6	Gln615	Gln535	−71.1	−79.1
cAMP	−8.2	−8.6	Gln615	Asp438, Asn199, Asn375	−43.7	−44.4

---

Fig. 4 (A) Docked pose of AQ-390/42425549 (AQ) in PDE4B (PDB ID 3G45). The molecule showed H-bond with Tyr274 and Gln615 (blue solid line) and pi–pi interaction with Phe618. These interactions were also observed in the molecular dynamics simulation of PDE4B-AQ complex. The values represents the average and standard deviation of different distances during last 3 ns of MD simulation. (B) RMSD plot of protein and ligand for the 10 ns MD simulation run.

3.5 Molecular dynamics simulation

Protein–ligand (PDE4B: AQ) complex was simulated for 10 ns. Fig. 4B shows the temporal root mean square deviation (RMSD) plot of protein residues (backbone atoms) and ligand (excluding hydrogens) for the 10 ns molecular dynamics (MD) simulation run. The early large variation stabilizes at 3 ns and remains stable during further simulation. We performed all our analyses on last 3 ns part of the MD simulation. RMSD values of protein and ligand during last 3 ns run were 1.73 ± 0.14 Å and 1.17 ± 0.08 Å respectively. AQ-PDE4B interactions which we observed in docking analyses, remain stable during MD simulation (Fig. 4A and ESI Fig. 1 and 2†). For example, the constant H-bonds with Gln615 and Tyr274 were observed during last 3 ns simulation. Apart from H-bonds, hydrophobic interaction with Phe618 was consistent during MD simulation (Fig. 4A and ESI Fig. 3†). The stable protein–ligand interactions during MD simulation support our previous docking results and confirmed the role of new molecules as potent PDE4B-inhibitors.

3.6 Biological evaluation

The hits obtained from pharmacophore based virtual screening and molecular docking analysis were purchased from respective vendors of Specs and initially evaluated for PDE4B inhibitory activity using the PDE4B2 Assay Kit from BPS BioScience. Table 5 lists the PDE4B % inhibitory activity for twenty one hit molecules, and indicates that thirteen compounds were highly active against PDE4B (>60% inhibition at 30 μM). Further, to evaluate the selectivity of identified hits toward the PDE4B, these 13 active hits were evaluated for PDE4D inhibitory activity using the PDE4D2 Assay Kit from BPS BioScience. Out of thirteen, ten compounds were potent and selective against PDE4B. The non-selectivity of the ligands AO-022/43453692 and AJ-292/42152568 might be due to the predicted H-bond with PDE4D in molecular docking studies (Table 4). Ligand AP-124/43383688 is more selective to PDE4D that PDE4B, this behaviour justify the low difference between the glide score of ligand in PDE4B and PDE4D.

Table 5 Inhibitory activity of shortlisted molecules against PDE4B and PDE4D^a

Molecules	% inhibition of PDE4B	% inhibition of PDE4D
a n.d.: not determined.
AQ-390/42425549	73.16	30.39
AG-690/11972161	65.98	26.96
AG-690/36276051	63.93	39.95
AG-205/36564043	68.24	37.75
AG-690/36276042	70.49	36.52
AP-964/40915318	61.27	30.39
AG-205/33161053	71.72	38.73
AG-690/15437723	77.87	39.46
AK-968/12386394	65.37	39.71
AF-399/15335138	70.08	37.01
AJ-292/42152568	78.48	59.56
AO-022/43453692	60.25	56.13
AP-124/43383688	63.93	67.16
AG-690/36873050	30.94	n.d.
AF-399/42316263	29.10	n.d.
AK-968/37005156	31.15	n.d.
AF-399/41980308	40.57	n.d.
AT-057/43469096	46.72	n.d.
AG-690/10252051	39.34	n.d.
AO-476/43407280	51.43	n.d.
AO-022/43390834	44.88	n.d.

---

4. Conclusions

A ligand-based pharmacophore model was generated using quinoline derivatives with PDE4B inhibitory activity. The developed model revealed critical pharmacophoric features responsible to elicit biological activity. Further, pharmacophore-based alignment of PDE4B inhibitors was used to derive 3D-QSAR. It helped to identify spatial arrangements of various substituents which may influence PDE4B binding, and thus, inhibitory activity. Among generated models, selected model showed significant correlation coefficient of 0.95 and also validated with the external r² prediction. Further, this model also explained the extent of electron withdrawing, hydrophobic, H-donor, positive ionic and negative ionic moieties in molecular structure influencing the PDE4B inhibition which was found to be consistent with reported pattern of biological activity.

Finally, 21 potential molecules were identified with PDE4B selectivity using virtual screening, docking and ADME properties analysis. Molecular dynamics study was carried out for a top scoring ligand, AQ and interaction pattern of ligand protein complex was found to be consistent with the docking results. All the hits were purchased and biologically evaluated for percentage inhibition of PDE4B and PDE4D with in vitro enzymatic assay. Out of these, thirteen hits have shown potent and ten showed selective inhibitory activity against PDE4B. All this data provides confidence level in the present study to find its application for designing the potent and selective PDE4B inhibitors.

Acknowledgements

The authors are indebted to Department of Science and Technology (DST), India for providing financial assistance to carry out this project (Grant No. SB/FT/CS-013/2012). The authors thank Vivek Kumar, IIT Delhi, India for critical reading of this manuscript and technical assistance for molecular dynamics work.

References

1. C. Burnouf and M.-P. Pruniaux, Curr. Pharm. Des., 2002, 8, 1255–1296.
2. W. Baumer, J. Hoppmann, C. Rundfeldt and M. Kietzmann, Inflammation Allergy: Drug Targets, 2007, 6, 17–26.
3. S. Oger, C. Méhats, E. Dallot, D. Cabrol and M.-J. Leroy, J. Immunol., 2005, 174, 8082–8089.
4. A. M. Feldman, M. R. Bristow, W. W. Parmley, P. E. Carson, C. J. Pepine, E. M. Gilbert, J. E. Strobeck, G. H. Hendrix, E. R. Powers and R. P. Bain, N. Engl. J. Med., 1993, 329, 149–155.
5. M. Boolell, S. Gepi-Attee, J. Gingell and M. Allen, Br. J. Urol., 1996, 78, 257–261.
6. J. E. Souness, D. Aldous and C. Sargent, Immunopharmacology, 2000, 47, 127–162.
7. S. Jin, S.-L. Ding and S.-C. Lin, Chang Gung Medical Journal, 2012, 35, 197–210.
8. T. J. Torphy, Am. J. Respir. Crit. Care Med., 1998, 157, 351–370.
9. N. J. Press and K. H. Banner, Prog. Med. Chem., 2009, 47, 37.
10. C. Page and D. Spina, in Phosphodiesterases as Drug Targets, Springer, 2011, pp. 391–414.
11. F. E. McCann, A. C. Palfreeman, M. Andrews, D. P. Perocheau, J. J. Inglis, P. Schafer, M. Feldmann, R. O. Williams and F. M. Brennan, Arthritis Res. Ther., 2010, 12, R107.
12. M. D. Houslay, P. Schafer and K. Y. Zhang, Drug Discovery Today, 2005, 10, 1503–1519.
13. B. J. Lipworth, Lancet, 2005, 365, 167–175.
14. Z. Huang and J. A. Mancini, Curr. Med. Chem., 2006, 13, 3253–3262.
15. A. Lerner and P. M. Epstein, Biochem. J., 2006, 393, 21–41.
16. A. Hatzelmann and C. Schudt, J. Pharmacol. Exp. Ther., 2001, 297, 267–279.
17. K. F. Rabe, E. D. Bateman, D. O'Donnell, S. Witte, D. Bredenbröker and T. D. Bethke, Lancet, 2005, 366, 563–571.
18. C. D. Manning, M. Burman, S. B. Christensen, L. B. Cieslinski, D. M. Essayan, M. Grous, T. J. Torphy and M. S. Barnette, Br. J. Pharmacol., 1999, 128, 1393–1398.
19. J. Spond, R. Chapman, J. Fine, H. Jones, W. Kreutner, T. Kung and M. Minnicozzi, Pulm. Pharmacol. Ther., 2001, 14, 157–164.
20. S.-L. C. Jin and M. Conti, Proc. Natl. Acad. Sci. U. S. A., 2002, 99, 7628–7633.
21. A. Robichaud, P. B. Stamatiou, S.-L. C. Jin, N. Lachance, D. MacDonald, F. Laliberté, S. Liu, Z. Huang, M. Conti and C.-C. Chan, J. Clin. Invest., 2002, 110, 1045–1052.
22. S. Guariento, O. Bruno, P. Fossa and E. Cichero, Mol. Diversity, 2016, 20, 77–92.
23. N. S. Tripuraneni and M. A. Azam, J. Mol. Model., 2015, 21, 1–12.
24. X. Dong, J. O. Ebalunode, S. J. Cho and W. Zheng, J. Chem. Inf. Model., 2010, 50, 240–250.
25. Z. Yang and P. Sun, Int. J. Mol. Sci., 2007, 8, 714–722.
26. P. Srivani, K. Kiran and G. N. Sastry, Indian J. Chem., Sect. A: Inorg., Bio-inorg., Phys., Theor. Anal. Chem., 2006, 45, 68–76.
27. A. K. Chakraborti, B. Gopalakrishnan, M. E. Sobhia and A. Malde, Bioorg. Med. Chem. Lett., 2003, 13, 1403–1408.
28. P. Fossa, G. Menozzi and L. Mosti, Quant. Struct.-Act. Relat., 2001, 20, 17–22.
29. N. H. N. Moorthy, N. M. Cerqueira, M. J. Ramos and P. A. Fernandes, RSC Adv., 2013, 3, 23409–23422.
30. N. S. Sapre, N. Jain, S. Gupta and N. Sapre, RSC Adv., 2013, 3, 10442–10451.
31. M. D. Woodrow, S. P. Ballantine, M. D. Barker, B. J. Clarke, J. Dawson, T. W. Dean, C. J. Delves, B. Evans, S. L. Gough and S. B. Guntrip, Bioorg. Med. Chem. Lett., 2009, 19, 5261–5265.
32. C. J. Lunniss, A. W. Cooper, C. D. Eldred, M. Kranz, M. Lindvall, F. S. Lucas, M. Neu, A. G. Preston, L. E. Ranshaw and A. J. Redgrave, Bioorg. Med. Chem. Lett., 2009, 19, 1380–1385.
33. S. Mukherjee and M. Pal, Drug Discovery Today, 2013, 18, 389–398.
34. S. L. Dixon, A. M. Smondyrev, E. H. Knoll, S. N. Rao, D. E. Shaw and R. A. Friesner, J. Comput.-Aided Mol. Des., 2006, 20, 647–671.
35. R. P. Gangwal, M. V. Damre, N. R. Das, G. V. Dhoke, A. Bhadauriya, R. A. Varikoti, S. S. Sharma and A. T. Sangamwar, J. Mol. Graphics Modell., 2015, 57, 89–98.
36. K. Naganuma, A. Omura, N. Maekawara, M. Saitoh, N. Ohkawa, T. Kubota, H. Nagumo, T. Kodama, M. Takemura and Y. Ohtsuka, Bioorg. Med. Chem. Lett., 2009, 19, 3174–3176.
37. R. Adepu, D. Rambabu, B. Prasad, C. L. T. Meda, A. Kandale, G. R. Krishna, C. M. Reddy, L. N. Chennuru, K. V. Parsa and M. Pal, Org. Biomol. Chem., 2012, 10, 5554–5569.
38. C. J. Mitchell, S. P. Ballantine, D. M. Coe, C. M. Cook, C. J. Delves, M. D. Dowle, C. D. Edlin, J. N. Hamblin, S. Holman and M. R. Johnson, Bioorg. Med. Chem. Lett., 2010, 20, 5803–5806.
39. T. Goto, A. Shiina, T. Yoshino, K. Mizukami, K. Hirahara, O. Suzuki, Y. Sogawa, T. Takahashi, T. Mikkaichi and N. Nakao, Bioorg. Med. Chem. Lett., 2013, 23, 3325–3328.
40. J. O. Ebalunode, Z. Ouyang, J. Liang and W. Zheng, J. Chem. Inf. Model., 2008, 48, 889–901.
41. A. B. Burgin, O. T. Magnusson, J. Singh, P. Witte, B. L. Staker, J. M. Bjornsson, M. Thorsteinsdottir, S. Hrafnsdottir, T. Hagen and A. S. Kiselyov, Nat. Biotechnol., 2010, 28, 63–70.
42. M. P. Jacobson, D. L. Pincus, C. S. Rapp, T. J. Day, B. Honig, D. E. Shaw and R. A. Friesner, Proteins: Struct., Funct., Bioinf., 2004, 55, 351–367.
43. J. A. Maier, C. Martinez, K. Kasavajhala, L. Wickstrom, K. E. Hauser and C. Simmerling, J. Chem. Theory Comput., 2015, 11, 3696–3713.
44. W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33–38.
45. S. Pirhadi, F. Shiri and J. B. Ghasemi, RSC Adv., 2015, 5, 104635–104665.
46. T. A. Halgren, R. B. Murphy, R. A. Friesner, H. S. Beard, L. L. Frye, W. T. Pollard and J. L. Banks, J. Med. Chem., 2004, 47, 1750–1759.

---

Footnote

† Electronic supplementary information (ESI) available: List of all hypotheses produced. List of all quinoline analogues with their QSAR-model set and pharma set classification. Predicted absorption, distribution, metabolism and excretion (ADME) properties of all PDE4B selective molecules. Chemical structure of all PDE4B selective molecules and their Specs IDs. Data for molecular dynamics simulation and table for external r² prediction of QSAR model. See DOI: 10.1039/c6ra11210b

---