Ligand-based pharmacophore modeling, virtual screening, and 2D quantitative structure-activity relationship performance on anti-Hepatitis B virus flavonols

1. Introduction

Hepatitis B virus (HBV) is a significant worldwide health problem because it can lead to acute and chronic liver infections [1,2]. HBV is a DNA virus with partial double-stranded characteristics, and a member of the Hepadnaviridae family. The global health impact is significant, affecting more than 300 million individuals with chronic HBV (CHB) [3] infections and causing approximately 820,000 deaths each year from complications such as cirrhosis and hepatocellular carcinoma [4]. One major obstacle presented by HBV is its ability to persist in the hepatocytes of the host, which can result in the development of chronic infection from an initial acute infection [5,6,7,8,9]. It is also associated with miscarriage and premature birth in pregnant women who are infected [10,11]. In long-term, CHB condition raises the chance of serious liver diseases, such as liver cirrhosis, liver failure, and a higher probability of developing hepatocellular carcinoma (HCC). HBV mainly spreads through contact with infected blood or bodily fluids, posing a risk for healthcare workers, drug users, and those having unprotected sex [12]. Moreover, HBV can be passed from mother to child during childbirth [8,9]. HBV deceptively can stay symptomless for long periods, causing gradual immune-mediated damage to the liver, and its detection is limited to the sensitivity of the currently available techniques. This feature, combined with the potential for transmission by asymptomatic carriers, highlights the crucial need for effective prevention, early detection, and control of HBV infections worldwide.

Currently available HBV treatment strategies focus on inhibiting viral replication, decreasing liver inflammation, and avoiding future complications [3,13,14]. Nucleoside analogs (such as lamivudine, tenofovir) and nucleotide analogs (such as entecavir) are essential components of HBV treatment. These antiviral drugs function by blocking the reverse transcriptase enzyme, which is necessary for HBV DNA replication. Pegylated interferon-alpha, along with other interferon-based treatments, provides both antiviral and immunomodulatory benefits. Recently, there have been promising developments with direct-acting antivirals (DAAs) like Bulevirtide (Myrcludex-B) [15,16,17,18] and capsid assembly modulators (CAMs) [19,20,21,22] that are progressing in different clinical trial phases, targeting specific stages in the HBV life cycle. The development of innovative therapeutic advancements in treating CHB infection is primarily driven by the need to overcome several challenges that hinder the achievement of successful and persistent virological responses with current FDA-approved medications. Accordingly, the emergence of drug resistance is a major expected concern in viral infection, especially when nucleos(t)ide analogs are used for an extended period of time [23]. Patients may need continuous treatment throughout their lives to ensure that the virus remains suppressed, which may result in challenges with adherence and possible side effects [24]. Interferon-based treatments frequently lead to flu-like symptoms, tiredness, and depression [25], which reduce their tolerability. Therefore, there is still a crucial demand for effective treatment approaches that are strong, safe, and able to achieve long-lasting virological response without the emergence of resistance. In this regard, nature always offers potential sources of yet remained to be discovered antiviral chemicals.

Flavonoids are a varied category of natural substances that are plentifully present in fruits, vegetables, and medicinal plants [26]. These compounds are being noticed for their possible antiviral effects, displaying potential in blocking different phases of HBV life cycle as we have reviewed before [27]. Research has shown that flavonoids can disrupt viral entry, replication, and assembly, making them promising options for antiviral treatment. A variety of different flavonoid subclasses are introduced with different mechanisms of action against HBV. Accordingly, Quercetin is distinguished among flavonoids for its capacity to block the production of HBsAg and HBeAg, important indicators of HBV replication [28,29,30]. Moreover, Baccharis species like B. spicata contain 5-caffeoylquinic acid (5-CQA) and 3,5-dicaffeoylquinic acid (3,5-DCQA) which have shown antiviral effects against HBV [31]. These compounds exhibit potential in decreasing viral DNA synthesis and viral reproduction, making them possible contenders for additional investigation in creating new antiviral treatments for HBV. Ongoing research in this field reveals that the wide variety of flavonoids discovered in the natural world could offer a valuable resource of possible antiviral treatments for fighting against HBV and other viral infections.

Due to diversification of flavonoids from different herbal species they probably possess different pharmacological activities that might be preserved among flavonoids from same subclass of different plant species. In this regard, a pharmacophore-based model was generated with 100% specificity to discovery novel flavonoids with potent anti-HBV activities. Furthermore, the flavonoids with antiviral activities against HBV used for establishing a robust and reliable QSAR model for predicting the probable activities of novel compounds. The findings of the present study also can be deployed for creating the stand-alone applications, and to highlight plants with potent flavonoids from diverse herbariums.

2. Materials and Methods

Retrieving Flavonoids’ Structures

The flavonoid compounds experimentally approved for their anti HBV activities were reviewed before [27], and the 3D structures of them were obtained from chemical databases, PubChem (https://pubchem.ncbi.nlm.nih.gov/) [32] and ChEMBL (https://www.ebi.ac.uk/chembl) [33].

Establishing Pharmacophore Model

The flavonol-based pharmacophore model was created with LigandScout v4.4 [34] as described before [35]. The compounds were grouped into the flavonoid subclasses as reviewed before [27]. Accordingly, a set of nine flavonols with anti-HBV activities, including Kaempferol [29,36], Isorhamnetin [36], Icaritin [37], Hexamethoxyflavone (Hex) [38], Chrysoeriol-6-C-b-D-boivinopyranosyl-41-O-b-D-glucopyranoside [39], Hyperoside (quercetin-3-O-galactoside) [40,41], Quercetin-3-O-glucuronide (Q3G) [29], [28,36], and Myricetin 3-rhamnoside [30] was used to train the model.

Furthermore, eight flavones (Luteolin 7-O-glucuronide [42], Isovitexin [43], Isoorientin [44], [45,46], [47], Swertisin [48], 4K [49], and Robustaflavone [50]), three flavanones ((-)-Epigallocatechin-3-gallate [51,52,53,54,55,56], oolonghomobisflavan C (OHBF-C) [57], and sikokianin A [58]), one anthocyanin (proanthocyanidin [57]), one chalcones (rosmarinic acid [59]), one biflavonoid ([60]), and one isoflavone (compound 8f [61]) were used to test and validate the model. Further polyphenols and triterpenes with anti-HBV activities, including PHAP [62], Nobiletin, Linalool [63], Betulinic acid [64], Ursolic acid [65], Taraxasterol [66], and Solamargine [67] were also included as decoys.

The compounds were clustered according to the pharmacophore RDF-code similarity measure with maximum cluster distance calculation methods. To generate the flavonol-based pharmacophore model, the conformers were generated according to the software iCon best setting. Briefly, maximum number of conformers value was set to 200 with energy window of 20.0 and max pool size of 4000. The model was created based on pharmacophore fit and atom overlap scoring function. Also, pharmacophore type was based on Merged Feature Pharmacophore to assure each feature is scored and those that do not match all input molecules are removed. Further settings were kept as its default. The best developed model based on the higher score was adopted for virtual screening.

Screening Large Chemical Database

According to our previous study on screening chemicals in a large database of natural and herbal products for the discovery of novel antiviral compounds [68], the PharmIt server (https://pharmit.csb.pitt.edu/search.html) was used for high-throughput screening (HTS) of the flavonol-derived pharmacophore model. In this regard, eleven built-in libraries comprised of 1,652,702,330 conformations of 347,839,756 compounds, including CHEMBL32, ChemDiv, ChemSpace, Enamin, MCULE, MCULE_ULTIMATE, MolPort, NCI_Open, PubChem, WuXi LabNetwork, and ZINC were screened. No filter was applied on the screening process. Moreover, the hit that matched all the pharmacophore features was preserved to provide a library for assessing the precision of the model.

Creating Chemical Libraries and Assessing Model Precision

Two datasets of chemical libraries were made to evaluate the accuracy of the established flavonol-based pharmacophore model. The relatively fitted hits (with diverse RMSD values) from the PharmIt screening results was employed to establish the first library set of compounds relatively matched to the model. As reported before, the hits from the built-in datasets might composed of duplications. The duplications were removed with Open Babel command line as described before [69]. The second decoy set was derived from the previous study consisting 1.7K Lipinski's rule of five-filtered FDA approved drugs to evaluate the specificity of the model [70]. The accuracy is reported in a receiver operating characteristic (ROC) curve.

Establishing a Predictive QSAR Equation

Same combinations that used to develop the pharmacophore model, used to establish a reliable predictive 2D quantitative structure-activity relationship (QSAR). Accordingly, RDKit (http://www.rdkit.org/), Mordred [71], and Dragon [72] tools were used to calculate 202, 1348, and 2489 descriptors for each molecule, respectively, in the setting of UseGalaxy server (https://usegalaxy.eu/) [73,74]. The calculated descriptors were checked to identify and remove duplications and constant or near-constant variables. The rest of the descriptors were retained to predict the biological activity (IC₅₀ μM) of flavonol compounds and the response variable. The response biological activity variable was extracted from literature. The activity of compounds that was reported as g.L^-1 was converted to μM according to the following formula:

$$\mathrm{\mu }M=\left(\frac{\mathrm{\mu }g}{ML}/MW\right)\times 1000$$

Where µM is the concentration in micromolar, µg/mL is the concentration in micrograms per milliliter, and MW is the predicted molecular weight.

Microsoft Excel (MS-Excel) v2019 was used to perform the QSAR modeling with multiple linear regression (MLR) as we have reported before [75,76]. Accordingly, the model was built using a 95% confidence interval (CI) and a tolerance level of 0.0001. The model was validated using two separate sets of test observations: one with random samples of non-flavonol flavonoids, and another with non-flavonoid compounds. This method enabled assessing the model's replicability by conducting two separate runs of the multiple linear regression (MLR) test with the same explanatory variables.

The process of selecting the model involved finding the best model based on the adjusted R². The number of predictor variables in the models were considered within a range of 2 to 4. HC0 was utilized for applying covariance adjustments through the Newey West (adjusted) method of heteroscedasticity. Cook's distances threshold > 1 was used to identify outliers or leverage observations in the statistical model. The established model was also tested on a set of anti-HBV flavonols rather than those used for model training. The Goodness of Fit statistics, including coefficient of determination (R²), adjusted R², mean squared error (MSE), root mean squared error (RMSE), and mean absolute percentage error (MAPE). Furthermore, the consistency between the predicted and actual observations in the external flavonol set was determined by squared cross-validated coefficient (Q²).

Multicollinearity Measurement

To evaluate multicollinearity in the QSAR model, the Variance Inflation Factor (VIF) was calculated among the independent variables. The following formula was used to calculate the VIF values for each included predictor.

$${VIF}_i=\frac{1}{1-R_i^2}$$

Where ​$R_i^2$ represents the coefficient of determination achieved by performing a regression of the $i$^th independent variable on all other independent variables. If VIF was greater than 10, it was considered as a multicollinearity. This was carried out to assurance the stability and credibility of the model's predictions.

Principal Component Analysis

The training flavonols set and selected predictors was used for PCA analysis. This was carried out in order to reduce the dimensionality of data, visualize connections between variables, and identify trends in compounds through their molecular properties. Additionally, the exact molecular weight (ExactMolWt) was incorporated as a dummy variable [76]. This indicated the weight of each molecule in the data set. A PCA was carried out using Pearson correlation to investigate the relationship between variables. Additionally, the variables were not standardized before conducting PCA. This means that in the PCA analysis the original scale of the variables was maintained. The principal components (PCs) were subjected Varimax rotation with Kaiser normalization. Five PCs were identified within the dataset to offer an overview of the underlying framework. The findings consisted of a symmetrical biplot that displayed the connections between compounds and features in one graph. The coefficients of the biplot were automatically calculated for visualization. Furthermore, a Bootstrap (on 50 samples) observation chart was developed to assess the reliability of the PCA results.

Applicability Domain (AD) Assessment

The Euclidean distance ($d$) was determined for each compound in both the training and external prediction sets using MS-Excel v2019 using the formula described below. A threshold for determining the validity range was calculated by summing the mean $d$ in the training set and the standard deviation, multiplied by two. This limit represents the maximum distance at which a compound is judged to be within the model's AD. The $d$ value of each chemical in the external prediction set was also calculated. If a compound's distance surpassed the threshold, it was categorized as 'Outlier'; otherwise, it was labeled as 'Inlier'.

$$\begin{gathered} d=\sqrt{{\left({PredictedIC}_{50}-{ObservedIC}_{50}\right)}^2} \\ Threshould=Meanofd+\left(2\times Std.devofd\right) \end{gathered}$$

3. Results

Flavonol-Based Pharmacophore Models

Among the generated pharmacophore models, the best one with the score of 130.7671 was selected (Figure 1). The model was comprised of 57 features, including 1 hydrophobic center, 3 aromatic rings, 7 H bond acceptor regions (HBA), 4 H bond donors, and 42 exclusive volumes. Different range of conformer (1 to 200) was generated for each ligand to find the best matched conformation. The mean of pharmacophore score for training set was 121.01 ± 11.96 (minimum and maximum of 103.4 and 132.54, respectively). Also, the score for the validation set was 49.35 ± 31.36 (minimum of 30.23 and maximum of 106.24). The established model was used for HT screening.

Virtual Screening

The flavonol-based model was deployed for a high-throughput screening (HTS) in the PharmIt database, resulting in 1,032 hits. Figure 2 shows that more hits were found in the PubChem dataset. After removing repeated entries, 509 distinct hits were discovered and used for establishing of a library of active compounds to evaluate the model's precision. To further confirm the accuracy of the model, another set of FDA-approved chemicals was used. The validation process assumed that the library from HTS primarily had active compounds, which were then matched with the pharmacophore model using an algorithm that only chose those with a strong pharmacophore feature correspondence. On the other hand, even though the FDA database may not have enough flavonoids for specific model comparison, it was still used to test the overall performance of the model as a decoy set.

According to Figure 3, the model was assessed using ROC curve analysis, showing a sensitivity of 71% and a specificity of 100%. The discriminative ability of the model was highlighted by the Area Under the Curve (AUC) values, which were 1.00, 1.00, 1.00, and 0.85 for cutoffs of 1%, 5%, 10%, and 100%, respectively. Furthermore, the Enrichment Factor (EF) measurements highlighted the effectiveness of the model in detecting true positives, showing EF values of 4.4 at the tested thresholds. Most importantly, the model showed zero false positives (FP), validating its accuracy and reliability in identifying active compounds. Taken together, these findings show the strength and possible usefulness of the flavonol-based model in identifying compounds that match pharmacophores, providing important information for future drug discovery efforts.

The QSAR Model

The QSAR model was established to reliably predict the biological properties of the flavonols used in constructing the pharmacophore model. However, two observations were randomly separated to be used as the external validation set. Different models were created using varying numbers of random internal validation sets of the non-flavonol flavonoids or non-flavonoids. The best model with the highest adjusted R² was selected for the report. The selected model was trained using Kaempferol, Isorhamnetin, Icaritin, Hexamethoxyflavone, Chrysoeriol-6-C-b-D-boivinopyranosyl-41-O-b-D-glucopyranoside, Hyperoside, and Q3G. Also, the random validation test sets were consisted of matched numbers of flavonoids (Proanthocyanidin, Amentoflavone, Rosmarinic acid, Sakuranetin, OHBF-C, Luteolin 7-O-glucuronide and Isovitexin) or non-flavonoids (PHAP, Nobiletin, Linalool, Betulinic acid, Ursolic acid, Taraxasterol, and Solamargine). This was due to the limited number of training set, and preventing from bias in the model.

Three sets of 2D descriptors were generated using RDKit, Dragon, and Mordred. After removing duplications and constant or near-constant variables the final number of predictor variables was 759. The results of the variable selection of modeling process are summarized in Table 1. Based on the Type III sum of squares, only two variables, the average connectivity index of order 4 (X4A) and the quantitative estimate of drug-likeness (qed) had significant information to explain the variability of the biological activities. Among the explanatory variables, based on the Type III sum of squares (Table 2), X4A was more influential. No multicollinearity was observed for these predictors (VIFs < 10). It was also observed that same predictors are adequate for establishing the model either with non-flavonol flavonoids as the internal validation test set or non-flavonoids compounds. This suggests that same variables distinguish flavonols from other compounds, and determine the specificity of the model.

The fitness and performance of the selected model was also further measured in training- and test-sets by Goodness-of-fit statistics (see Table 3). Accordingly, the model achieved a R² value of 0.87 for the training set. The Adjusted R² value of 0.85 further confirmed the robustness of the model, indicating that the inclusion of descriptors appropriately balances model complexity and goodness of fit. As a findings the model’s error metrics were small, suggesting a relatively good accuracy of the trained model in predicting the observations. Also, Q² (0.9) was significantly high, indicating a good prediction power of the model on the external prediction set (see AD analysis results). Furthermore, the model error parameters in the validation sets were interestingly differ. MSE, RMSE, and MAPE were substantially lower in the model validated by the non-flavonol flavonoids than that validated by non-flavonoids. This suggests the association between all flavonoids, and yet discriminative power of the model on the flavonol compounds.

Further analysis of variance (ANOVA) results showed that the model is statistically significant (95%CI F = 36.91 df(2), p < 0.0001), indicating that the selected descriptors collectively have a significant impact on the biological activities of the compounds. The model explains a significant amount of the variability in the data, as indicated by the high sum of squares associated with the model (7589.40). The model over fitness was checked by the standardized coefficients measured for the variables in the mode (Table 4). Furthermore, the robustness of predictive strength of the model was evaluated by the analysis of residuals (Table 5 and Table 6).

The fitness plot (Figure 4a,c) also shows the best fitted linear regression of predicted observations within training set. Comparing models validated with non-flavonol flavonoids and non-flavonoids highlights unique features in the linearity and fitness of their training sets. The model tested with chemicals that are not flavonoids shows a stronger level of linearity and fitness (Figure 4b) when compared to the training set. However, the model confirmed using flavonoids also shows significant fitting ability. These results emphasize the importance of using distinct validation sets to thoroughly investigate structural variations in flavonoids. Furthermore, the Cook’s distance metrics (Figure 4b) was checked, and no outliers were observed (Cook’s distance < 1).

Applicability Domain (AD) Analysis

AD analysis was performed to determine the range of chemical space where the QSAR model can make reliable predictions, and to identify the boundaries of the training data in which predicted biological activities of the external set fall within, ensuring that the model is applicable to predict for same compounds. Unfortunately, the number of flavonols were limited. Of nine compounds, seven were randomly used to train the QSAR model, and the rest were used to perform Euclidian distance (d) calculation and determining AD of the model. Two flavonols, Quercetin and Myricetin 3-rhamnoside, used as an external validation set. The compounds' predicted biological activities were measured in terms of distance, and if it fell within the specified threshold, it was considered an inlier, indicating the model's applicability.

The mean calculated $d$ for the training set was 7.06 ± 5.56. As described above, the threshold was measured as 2x St.dev plus mean $d$. Accordingly, the threshold was set as 18.178. As showed in Figure 3, the calculated $d$ for both compounds were inlier, indicating the applicability of the model for other flavonols with anti-HBV activities.

Figure 5. AD of the predicted external validation set. The size of the circle is based on the calculated $d$

Figure 5. AD of the predicted external validation set. The size of the circle is based on the calculated $d$

PCA Measure of the Predictors

The principal component analysis (PCA) was conducted on the training set’s variables, including Biological activities, X4A, qed, and MW as a supplementary dummy variable. According to the Pearson correlation matrix, Biological activities exhibited a strong positive correlation with X4A (r = 0.70) and a moderate negative correlation with qed (r = -0.48). Bartlett's sphericity test indicated that the variables were not independent (Chi-square = 8.65, p = 0.03), supporting the significancy of the correlation between variables, and the need for PCA analysis. The eigenvalues of the components revealed that the first three components (F1 to F3) explained 59.09%, 38.89%, and 2.02% of the total variability, respectively. The cumulative percentage of variance explained by the first two components was 97.98% (Figure 6), suggesting that these components capture the majority of the dataset's variability.

It was found that the first component (F1) had a significant loading for Biological activities (0.74), X4A (0.60), and qed (-0.32), indicating their strong influence on this component. The second component (F2) exhibited significant loadings for X4A (0.55), qed (0.83), and MW (-0.80), suggesting their association with this dimension (see Figure 7a). In contrast, the third component (F3) had a minor eigenvalue of 0.06, indicating its marginal role in explaining the dataset's variance. Considering the scree plot method and the Kaiser criterion, the analysis supports retaining the first two principal components. These components, F1 and F2, cumulatively explain nearly 98% of the total variance in the dataset, making them suitable for further analysis and interpretation. Therefore, based on the proportion of variance explained and eigenvalues, these two PCs were retained to capture the essential patterns and relationships within the data.

As showed in Figure 7b, PC1 primarily measures the bioactivity of compounds. The second principal component (PC2) shows large positive associations with qed (0.90) and X4A (0.59), while negatively associated with MW (-0.80). PC2 primarily measures qed and structural complexity (X4A) of compounds.

As showed in Figure 7c, the factor scores from the PCA revealed distinct patterns among the flavonols. F1 showed a positive association with Kaempferol, Isorhamnetin, and Chrysoeriol-6-C-b-D-boivinopyranosyl-41-O-b-D-glucopyranoside, indicating these compounds have higher scores on this factor. In contrast, Factor 2 is negatively associated with Chrysoeriol-6-C-b-D-boivinopyranosyl-41-O-b-D-glucopyranoside, Hyperoside, and Q3G, suggesting these compounds have lower scores on this factor. This pattern was further confirmed with Bootstrap sampling analysis with a minor variation (Figure 7d).

4. Discussion

The development of effective therapeutic interventions against HBV infection is importance given its significant global health burden and potential for severe complications. Flavonoids, a broad range of natural substances present in many plants, are being studied for their possible antiviral effects, such as their capacity to disrupt different phases of the HBV lifecycle [27]. The variety of flavonoids offers an interesting chance to find new antiviral drugs, leading to the creation of models based on pharmacophores and QSAR to help identify and predict compounds with strong anti-HBV effects through virtual HTS. Previous studies have also highlighted the potential of utilizing computational modeling like QSAR to offer potential compounds with anti-HBV activities [75,76,77,78]. The present study not only adds to our knowledge of flavonoid pharmacology but also has implications for creating new antiviral treatments and finding plants with high levels of bioactive flavonoids.

The development and validation of pharmacophore models are essential in contemporary drug discovery endeavors [35,69,79,80,81]. In this study, a validated pharmacophore model based on flavonols was carefully developed, consisting of 57 features such as hydrophobic centers, aromatic rings, HBA and HBD regions, and exclusive volumes. The accuracy of the top model showed its potential to precisely screen compounds with pharmacophore characteristics related to flavonols. In this regard, The ROC curve analysis showed encouraging performance metrics, including a sensitivity of 71% and a specificity of 100%. The discriminative ability of the model was also confirmed by the AUC values. Accordingly, a curated collection of 359 distinct compounds was highlighted as active true positive chemicals with potential anti-HBV activities. These findings provide a great opportunity for vast number of research on these compounds in vitro. Also, since the training set of flavonols showed different anti-HBV mechanisms of action, the hits resulted of the screening potentially can display different anti-HBV activities. Also, future efforts on discovery of such compounds provide further chance of developing more sensitive model. The screening results were not focus of this study, yet the required information can be found in the supplementary material.

The QSAR model was constructed to predict the biological properties of flavonols utilized in the pharmacophore model. Noteworthy is the model's level of predictive power with no sign of collinearity or overfitting penalties, achieving an R² value of 0.87 in the training dataset. The model's explanatory capability was significantly due to the predictors X4A and qed. The model's robust performance was further confirmed by its ability to generalize well to external validation set, demonstrated by the 90% Q² value and inlier AD that was calculated for the external validation set. The developed model benefits from the biological activity of compounds in inhibiting the secretion of HBsAg. Therefore, the prediction is applicable on flavonoids with anti-HBsAg activities. This provides a further opportunity to investigate and reconstruct another model to search for variables predictive for biological activities of flavonols on other aspect of viral infection like HBeAg. In the present study, modeling different QSARs were not possible due to the limited number of flavonols that affecting other part of viral components. The presented model, was reevaluated with different number of input variables and other affirmatory statistics, and also revalidated with different number of validation sets, yet, the results were the same. In this regard, the validated models based on non-flavonol flavonoids and non-flavonoids are reported for the sake of comparison. Accordingly, one intriguing finding was the varying success of the model when tested with non-flavonol flavonoids compared to non-flavonoids. Although both models demonstrated good fit, the non-flavonoid model displayed stronger linearity and fitness with the training data. This implies that the model can accurately differentiate between flavonoids and other compounds by analyzing their structural characteristics. Nevertheless, the model that was verified with flavonoids also showed significant predictive ability, emphasizing the distinct structural differences found in flavonol substances. This might be attribute to the predictors. The average connectivity index of order 4 (X4A) is a topological index, and has been showed to have an influence in the binding affinities of Per- and polyfluoroalkyl substances (PFAS) ligand to Human serum albumin (HSA) receptor, and is significantly influential score in establishing QSAR model [77]. In a QSAR study by Hădărugă [82], X4A was identified as a significant structural parameter related to the inhibitory action of flavonoids against cytochrome P450. The study found that the biological activity increased with the decrease in the X4A parameter value. Specifically, higher X4A values were associated with lower inhibitory activity against cytochrome P450. This indicates that X4A plays a role in describing the hydrophilic-hydrophobic properties or polarizabilities of the flavonoid molecules in relation to their inhibitory activity against cytochrome P450. These suggest that X4A is a major contributor to topological structure of flavonoids that facilitate the binding of the chemical to the target. Here, the PCA results showed a significant positive correlation between biological activity of compounds with X4A, suggesting the increase in X4A will resulted in enhancing flavonols’ IC₅₀. Another effected predictor that was quantitative estimate of drug-likeness (qed). It was found that qed has orthogonal association with biological activity, and the variables are not redundant and each contributes unique information. qed is a measure used to assess the drug-likeness of chemical compounds, and combines several molecular properties associated with drug-like behavior into a single score. The mean qed calculated for flavonols was 0.41 ± 0.19, highlighting a minor trend toward druglike compounds. Accordingly, Lee et al, showed that flavonoids with qed ≥ 0.35 have a higher likelihood of exhibiting drug-like behavior, making them promising candidates for studies related to drug development or therapeutic effects [83]. Additionally, it was observed that the strategy of establishing a QSAR model with the qed score provides an efficient approach to discover flavonoid derivatives with potent anti- dipeptidyl peptidase-IV (DPP-IV) activity [84]. Taken together, features selected for developing a QSAR model in the present study is confirmed by the previous studies, suggesting both structural and drug likeness properties of flavonoids strength the predictive power of the models.

Although the pharmacophore and QSAR models have been largely successful, it is important to recognize their various limitations. The belief that the HTS library was mostly made up of active compounds may not completely represent the diversity of existing chemical libraries, particularly specialized flavonoid libraries, which could further improve the model's accuracy. Also, only flavonols were used to train both pharmacophore and QSAR models. This contributed to advancing the existing understanding of flavonols' effectiveness against HBV, however, utilizing models trained with different subclasses of flavonoids could be essential. We are currently studying models trained with flavonoids for uncovering anti-HBV drugs, particularly using different machine learning algorithms other than MLR to assess the accuracy of the models. Additionally, the results of the current research offer a way for identifying plants that contain flavonoids with anti-HBV properties by tracking the screened hits.

5. Conclusions

Druglike flavonoids with anti-HBV activities are possessed a yet to discover potentials. In the present study, an accurate flavonol-based pharmacophore model was validated for screening large chemical libraries, highlighting novel hits for further experimental studies. Also, a QSAR model was validated by two independent sets of chemicals for the sake of reproducibility of the model. In addition, the model was found to be applicable on other flavonols due to the predictive selected features, x4a and qed. Findings pave a great way of future anti-HBV drug discovery modeling, and experimental studies. Discoveries provide a promising path for future modeling and experimental research on anti-HBV drug development. This study built upon existing research on flavonoids as potential treatments for HBV, showing promising prospects for future HBV therapeutics.