Lipophilicity, Pharmacokinetic Properties, and Molecular Docking Study on SARS-CoV-2 Target for Betulin Triazole Derivatives with Attached 1,4-Quinone

A key parameter in the design of new active compounds is lipophilicity, which influences the solubility and permeability through membranes. Lipophilicity affects the pharmacodynamic and toxicological profiles of compounds. These parameters can be determined experimentally or by using different calculation methods. The aim of the research was to determine the lipophilicity of betulin triazole derivatives with attached 1,4-quinone using thin layer chromatography in a reverse phase system and a computer program to calculate its theoretical model. The physiochemical and pharmacokinetic properties were also determined by computer programs. For all obtained parameters, the similarity analysis and multilinear regression were determined. The analyses showed that there is a relationship between structure and properties under study. The molecular docking study showed that betulin triazole derivatives with attached 1,4-quinone could inhibit selected SARS-CoV-2 proteins. The MLR regression showed that there is a correlation between affinity scoring values (ΔG) and the physicochemical properties of the tested compounds.


Introduction
Nowadays, the search for new drugs requires the application of computational chemistry, including experimental and in silico analysis of physicochemical properties, pharmacokinetic features, ADMET analysis, and quantitative structure-activity relationship (QSAR) research [1,2]. The therapeutic potential of a drug depends on its distribution in the body. Often, substances with high biological activity show low bioavailability and high toxicity to healthy tissues. Nanocarriers, including liposomes, micelles, and polymer nanoparticles, enable the improvement of the biodistribution of substances. Incorporation of the hydrophobic core into hydrophilic nanotransporters can effectively influence the solubility of the drug, which contributes to the improvement of permeability through cell membranes [3,4].
A key parameter in the design of new active compounds is lipophilicity, which influences the solubility and permeability through membranes. Lipophilicity affects the pharmacodynamic and toxicological profiles of compounds. These parameters can be determined experimentally or by using different calculation methods. Experimental lipophilicity is usually determined by chromatographic methods, such as reversed phase thin layer chromatography (RP-TLC), normal phase thin layer chromatography (NP-TLC), or reversed phase high performance liquid chromatography (RP-HPLC). Theoretical methods are the second way to determine lipophilicity. However, the calculated value of lipophilicity Studies on the lipophilicity of 5,8-quinolinedione derivatives were carried out applying the RP-TLC and RP-HPLC methods. It was found that there is a relationship between lipophilicity and in silico pharmacokinetic properties of synthetic 5,8-quinolinedione derivatives [20][21][22][23]. As a continuation of our previous research, we determined the lipophilicity experimentally as well as by using computational methods for the betulin triazole derivatives with attached 1,4-quinone. The physicochemical and pharmacokinetic properties influence the bioavailability and biological activity of compounds when they are determined by computer programs. The correlation between lipophilicity and in silico determined structural parameters have been analyzed in this study.
The tested compounds contain two active moieties, i.e., 1,4-quinone and betulin. The main mechanism of activity for 1,4-quinone derivatives is the interaction with NQO1 protein, for which the overexpression is observed in many types of cancer cell lines [8,10,24]. The betulin derivatives show a wide spectrum of activities including anticancer, antiviral, antibacterial, and anti-inflammatory effects. The use of betulin and its derivatives in treatment is limited due to low bioavailability and low water solubility. For these reasons, the use of new drug nanoencapsulation procedures and/or nanoparticle-based delivery methods are a key issue. These are issues of nanomedicine interest [25,26].
The tested compounds contain two active moieties, i.e., 1,4-quinone and betulin. The main mechanism of activity for 1,4-quinone derivatives is the interaction with NQO1 protein, for which the overexpression is observed in many types of cancer cell lines [8,10,24]. The betulin derivatives show a wide spectrum of activities including anticancer, antiviral, antibacterial, and anti-inflammatory effects. The use of betulin and its derivatives in treatment is limited due to low bioavailability and low water solubility. For these reasons, the use of new drug nanoencapsulation procedures and/or nanoparticle-based delivery methods are a key issue. These are issues of nanomedicine interest [25,26].

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (trishydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 °C.
The obtained values of retardation factor (Rf) were converted to RM parameters according to Equation (1): The RM parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (RM0) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (φ0) was determined according to Equation (3):

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (trishydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 °C.
The obtained values of retardation factor (Rf) were converted to RM parameters according to Equation (1): The RM parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (RM0) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (φ0) was determined according to Equation (3):

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (trishydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 °C.
The obtained values of retardation factor (Rf) were converted to RM parameters according to Equation (1): The RM parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (RM0) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (φ0) was determined according to Equation (3):

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (trishydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 °C.
The obtained values of retardation factor (Rf) were converted to RM parameters according to Equation (1): The RM parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (RM0) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (φ0) was determined according to Equation (3):

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (trishydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 °C.
The obtained values of retardation factor (Rf) were converted to RM parameters according to Equation (1): The RM parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (RM0) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (φ0) was determined according to Equation (3):

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (trishydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 °C.
The obtained values of retardation factor (Rf) were converted to RM parameters according to Equation (1): The RM parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (RM0) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (φ0) was determined according to Equation (3):

Experimental Lipophilicity
The experimental lipophilicity was determined using the RP-TLC method according to the literature [20,21,37]. We used the modified silica gel as a stationary phase and (tris-hydroxymethyl)aminomethane (0.2 M, pH = 7.4) with acetone as a mobile phase. The percent of organic solvent volume was varied within the range of 60-90% in 5% increments.
The compounds 1-20 were dissolved in chloroform (1.0 mg/mL), then 5 µL of sample solution was applied to the chromatographic plates with a micropipette. Spots were visualized by spraying with 10% ethanol solution of sulfuric acid (VI) and then heated up to 110 • C.
The obtained values of retardation factor (R f ) were converted to R M parameters according to Equation (1): The R M parameter was calculated for every concentration of acetone and extrapolated to zero concentration of organic component in the mobile phase. The chromatographic parameter of lipophilicity (R M0 ) was calculated using Equation (2): where C is the concentration of acetone in the mobile phase, while b is the slope of the regression plot. The hydrophobic index (ϕ 0 ) was determined according to Equation (3):

Structure Optimization
The optimized chemical structures of compounds 1-16 were calculated using the DFT (B3LYP/6-311G+(d.p)) method implemented in the Gaussian 16 program package [44]. All local minima of energy were confirmed by the absence of an imaginary mode in the vibrational calculations. In our calculations, we applied the basis set of the diffuse function to heavy atoms (+) to obtain a better description of lone pair electrons with orbitals occupying a larger region of space. The obtained optimized structure is presented in Figure  S1. The geometries of compounds 1-16 were used to determine the HOMO-LUMO energy, quantum chemical descriptor, the molecular electrostatic potential, and the molecular docking study [45]. All obtained results were visualized using the GaussView, Version 5 software package (Gaussian Inc., Wallingford, CT, USA) [46].

Molecular Docking Study
The three-dimensional (3D) structures (in mol2 format) of the studied compounds were generated using the ChemOffice package (version 19.1, PerkinElmer, Waltham, MA, USA) [47]. Their low-energy conformations were calculated using Gaussian 16 (revision A.03) computer code [44] at the density functional theory (DFT, B3LYP) and 6-311+G(d,p) basis sets. Calculations were performed using the X-ray coordinates of chloroquine as the input structure obtained from the Cambridge Crystallographic Data Centre (CCDC ID: CDMQUI).
Target macromolecule for molecular docking study was obtained from the Protein Data Bank (https://www.rcsb.org/, accessed on 1 April 2020). We used the 3D crystal structures of COVID-19 main protein and papain-like protease of SARS CoV-2 (PDB ID: 5R7Z and 6W9C, respectively).
Ligands and proteins used in the calculations were prepared for docking using the AutoDockTools program (Molecular Graphics Laboratory, La Jolla, CA, USA) [48]. In this study, AutoDock Vina [49] tool compiled in PyRx [50] was employed to perform molecular docking. Implementing in AutoDock Vina scoring function was mostly inspired by X-score and combined an empirical free-energy force field with a Lamarckian Genetic Algorithm. The volume was set as 25

Correlation and Cluster Analysis
Based on the experimental and theoretical values of lipophilicity and molecular descriptor values, the correlation and cluster analysis were performed. All data used for the cluster analysis were standardized and the cluster analysis was based on the Euclidean distance. The analysis was carried out using the Statistica 13.1 software (TIBCO Software Inc., Palo Alto, CA, USA).

Experimental and Theoretical Lipophilicity
The hybrids 1-16 and substrates 17-20 were examined under the RP-TLC method. For each compound, the parameter R M was calculated from the retardation factor (R f ) according to Equation (1). Equation (2) was used to determine the R M0 and b values and the results are presented in Table 2. The high correlation coefficients (r = 0.976-0.997) for all compounds show good correlation between acetone concentration and R M value. The calibration curve is required for the conversion of R M0 to logP TLC . The following compounds were used as standard substances: acetanilide, prednisone, 4-bromoacetophenone, benzophenone, anthracene, dibenzyl, 9-phenylanthracene, and dichlorodiphenyltrichloroethane (DDT), for which the literature values of logP lit are in the range 1.21-6.38 [52,53]. The R M0 for the standard substance was determined in the same conditions as for compounds 1-20 (Table 3). b is the slope, r is the correlation coefficient and SD is the standard deviation for the linear relationship R M = R M0 + bC.
The calibration curve Equation (4) For the standard substances the logP TLC was calculated according to the calibration curve (Equation (4)). Figure S2 shows a good agreement between experimental and literature values of lipophilicity (r = 0.996).
Equations (3) and (4) were used to determine the logP TLC and the hydrophobic index (ϕ 0 ), respectively. The results are presented in Figure 2 and Table S1. Equations (3) and (4) were used to determine the logPTLC and the hydrophobic index (φ0), respectively. The results are presented in Figure 2 and Table S1. Comparing the logPTLC for hybrids 1- 16 and derivatives 17-20 shows that introduction of the 1,4-quinone fragment to the betulin moiety reduces the lipophilicity. The lipophilicity depends on the type of 1,4-quinone moiety and the order is as follows: (13)(14)(15)(16).
The hydrophobic indexes for hybrids 1-16 and compounds 17-20 are in the ranges of 79.81-88.21 and 82.62-87.13, respectively. A higher value of φ0 means that compounds 1-20 are less soluble in water. Comparing the hydrophobicity of derivatives 1-16 and 17-20 shows that the introduction of 1,4-quinone moiety slightly affects the solubility in water. In the series of 1-16, no relationship between the type of 1,4-qinone moiety and hydrophobicity index has been observed.
The hydrophobic indexes for hybrids 1-16 and compounds 17-20 are in the ranges of 79.81-88.21 and 82.62-87.13, respectively. A higher value of ϕ 0 means that compounds 1-20 are less soluble in water. Comparing the hydrophobicity of derivatives 1-16 and 17-20 shows that the introduction of 1,4-quinone moiety slightly affects the solubility in water. In the series of 1-16, no relationship between the type of 1,4-qinone moiety and hydrophobicity index has been observed.
Theoretical values of lipophilicity were determined with the available programs [38][39][40][41]. The calculated values of logP cover a wide range from 4.11 to 12.19 depending on the mathematical module used by programs. The theoretical results are presented in Figure 3 and Table S2.
For 1-20, all programs show that lipophilicity depends on the type of substituent at the C-3 position of the betulin moiety, and this correlation is consistent with the experimental results. As shown in Figure 3, compounds 17-20 have lower lipophilicity than hybrids 1-16, which is not in agreement with the experimental results. The exception is the MLOGP program, for which the theoretical lipophilicity for derivatives 17-20 are similar to the experimental values (Table S2). Comparing the logP values for 1-16 shows that the theoretical lipophilicity slightly depends on the type of 1,4-quinone moiety. Moreover, hybrids with 5,8-quinolinedione (1-4) and 5,8-isoquinolinedione (5-8) moiety have the same molecular formula. In this case, atomistic (WLOGP) and topological (MOLGP and SILCOS-IT) methods exhibit the same value of lipophilicity, but the LogP TLC is different. For this reason, to better predict the lipophilicity, it is important to find a correlation between its experimental and theoretical values. Due to the structural differences of hybrids 1-16 and betulin derivatives 17-20, separate equations were developed for these two groups of compounds (Table 4). Pharmaceutics 2021, 13, x 8 of 22 For 1-20, all programs show that lipophilicity depends on the type of substituent at the C-3 position of the betulin moiety, and this correlation is consistent with the experimental results. As shown in Figure 3, compounds 17-20 have lower lipophilicity than hybrids 1-16, which is not in agreement with the experimental results. The exception is the MLOGP program, for which the theoretical lipophilicity for derivatives 17-20 are similar to the experimental values (Table S2). Comparing the logP values for 1-16 shows that the theoretical lipophilicity slightly depends on the type of 1,4-quinone moiety. Moreover, hybrids with 5,8-quinolinedione (1-4) and 5,8-isoquinolinedione (5-8) moiety have the same molecular formula. In this case, atomistic (WLOGP) and topological (MOLGP and SILCOS-IT) methods exhibit the same value of lipophilicity, but the LogPTLC is different. For this reason, to better predict the lipophilicity, it is important to find a correlation between its experimental and theoretical values. Due to the structural differences of hybrids 1-16 and betulin derivatives 17-20, separate equations were developed for these two groups of compounds (Table 4).   The relationship between the lipophilicity and structure of compounds 1-20 was analyzed by cluster analysis (Figure 4).
As shown in Figure 4, compounds 1-20 are arranged in two main clusters. The first cluster consists of betulin derivatives 17-20 and the second one of hybrids 1-16. The high value of Euclidean distance between two clusters suggests that there is a low correlation between structures of these two groups of compounds.
The second cluster can be divided into four subclusters. In the first subcluster, the hybrids are arranged according to the type of betulin moiety, which means that they have oxo group at the C-3 position of the betulin moiety. The second and third subclusters consist of 1,4-naphthoquinone (13, 15, and 16)   The relationship between the lipophilicity and structure of compounds 1-20 was analyzed by cluster analysis (Figure 4). As shown in Figure 4, compounds 1-20 are arranged in two main clusters. The first cluster consists of betulin derivatives 17-20 and the second one of hybrids 1-16. The high value of Euclidean distance between two clusters suggests that there is a low correlation between structures of these two groups of compounds.

Physicochemical and Pharmacokinetic Properties
In the early drug discovery process, one of the important stages is to optimize the physicochemical properties of a potential drug. Based on the observation of the physicochemical properties of oral drugs, Lipinski formulated the Rule of Five. According to this

Physicochemical and Pharmacokinetic Properties
In the early drug discovery process, one of the important stages is to optimize the physicochemical properties of a potential drug. Based on the observation of the physicochemical properties of oral drugs, Lipinski formulated the Rule of Five. According to this rule, the lipophilicity (logP), molecular weight (MW), number of acceptors (HA), and donors (HD) of hydrogen bond should be a multiple of five [54]. In an attempt to improve the prediction of bioavailability, the rules were expanded by Veber to the number of rotatable bonds (RB) and topological polar surface area (TPSA) are also determined [55,56]. The molecular descriptors for tested hybrids 1-16 are presented in Table 5.
The tested hybrids 1-16 have a molecular mass above 500 g/mol, which means they do not meet the mass criterion. All compounds have less than five hydrogen bond donors (HD = 0-1). The 5,8-quinolinedione and 5,8-isoquinolinedione derivatives with hydroxyl and oxo groups at C-3 position of betulin (1-2, 5-6, and 9-10) have less than 10 hydrogen bond acceptors (HA = 9-10) and the logP TLC less than 5, which means that these compounds meet three of four Lipinski rules. Only compounds with propanoyl substituents at the C-3 position of the betulin moiety (4, 8, 12, and 16) do not meet the Veber's rules concerning the number of rotatable bonds (RB < 10). The TPSA for 1-16 is less than 140 Å, which determines the oral bioavailability.
For compounds 1-16, the relationship between the type of 1,4-quinone moiety, the physicochemical parameters and experimental lipophilicity was analyzed by means of a dendrogram ( Figure 5).  The similarity analysis shows two main clusters. The first cluster consists of 1,4-naphthoquinone derivatives (13 and 14) and the second one is divided into four subclusters. Analysis of subclusters show that compounds with the 5,8-quinolinedione (1-4) and 2methyl-5,8-quinolinedione (5-8) exhibit similar properties. Comparing the properties of 5,8-quinolinedione (1-4) and 5,8-isoquinolinedione (9)(10)(11) shows that the position of nitrogen atom in heterocyclic ring influences their physicochemical parameters. The similarity analysis shows that nitrogen atom significantly affects properties of the tested derivatives.
The multilinear regression (MLR) Equation (5) The correlation coefficient shows good agreement between physicochemical parameters and experimental lipophilicity. A comparison of the experimental and calculated parameters for the investigated compounds are shown in Table S3. The result shows that lipophilicity could be determined by in silico parameters. The similarity analysis shows two main clusters. The first cluster consists of 1,4naphthoquinone derivatives (13 and 14) and the second one is divided into four subclusters. Analysis of subclusters show that compounds with the 5,8-quinolinedione (1-4) and 2methyl-5,8-quinolinedione (5-8) exhibit similar properties. Comparing the properties of 5,8-quinolinedione (1-4) and 5,8-isoquinolinedione (9)(10)(11) shows that the position of nitrogen atom in heterocyclic ring influences their physicochemical parameters. The similarity analysis shows that nitrogen atom significantly affects properties of the tested derivatives.
The multilinear regression (MLR) Equation (5) The correlation coefficient shows good agreement between physicochemical parameters and experimental lipophilicity. A comparison of the experimental and calculated parameters for the investigated compounds are shown in Table S3. The result shows that lipophilicity could be determined by in silico parameters.
For 1-16, the logBB are lower than −1, which means these hybrids poorly cross the brain-blood barrier. Compounds do not pass into the central nervous system when the logPS is less than −2. This result shows that the tested derivatives 1-16 are not neurotoxic.
The similarity analysis showed no correlation between pharmacokinetic parameters and the structure of compounds 1-16 ( Figure S3). However, the cluster analysis showed that only betulin moiety influences the pharmacokinetic parameters of hybrids.

Molecular Properties
The energy of HOMO and LUMO orbitals determines the ability of the molecule to donated or receive an electron. The energy of orbitals can be used to calculate the global reactivity descriptors, such as: ionization potential (I), electron affinity (A), hardness (η), chemical potential (µ), electronegativity (χ), and electrophilicity index (ω) [61].
The HOMO and LUMO orbitals of 1-16 are presented in Figure S4. Table 7 shows the HOMO and LUMO energy, the energy gap (∆E = E HOMO − E LUMO ), the global reactivity descriptors, and dipole moment (DM). For all compounds 1-16, the LUMO orbitals are mainly delocalized at the 1,4-quinone moiety. Localization of HOMO orbitals depends on the type of substituent at the C-3 position of the betulin moiety. The HOMO orbitals of hybrids with hydroxy and acyl groups (1, 3-5, 7-9, 11-13, and 15-16) at the C-3 position are mainly delocalized at the isopropenyl group and A ring of betulin. For hybrids with the oxo group (2, 6, 10, and 14), these orbitals are delocalized at the D ring and oxo group ( Figure S4). The molecular properties influence the reactivity and stability of molecules. The high value of E HOMO and low value of E LUMO show that the compounds are highly reactive with nucleophilic molecules. The energy gap ∆E depends on the van der Waals interaction between molecules and its value could be correlated with the biological activity of the compounds [62,63]. The hybrids 1, 5, 9, and 13 are characterized by the lowest chemical hardness and the highest negative value of chemical potential, which means they are comparatively soft with high polarizability compared with other compounds in this series.
The arrangement of nucleophilic and electrophilic regions of a molecule influences its interaction with a biological target through hydrogen bonds and hydrophobic interactions. The distribution of the positive and negative charges of molecule is described by the molecular electrostatic potential map (MEP) [64,65]. The different charges are represented by different colors, which means that the red, blue, and green areas are negative, positive, or neutral, respectively [66]. The MEP's are sketched for an order of (−47,935) to 47,935 kcal/mol. The maps for hybrids 1, 5, 9, and 13 present in Figure 6, while for the rest of the compounds, they are shown in Figure S5.
For all hybrids 1-16, the nucleophilic regions (red color) are localized in four main area. The first and second areas are localized on 1,4-quinone moiety. The first area contains the carbonyl atom at C-5 and the second area the carbonyl group at C-8, and the nitrogen atom. The third area includes the triazole linker, and the fourth, the substituent at the C-3 position of the betulin moiety. The electrophilic regions (blue color) are localized near the methine group at 1,4-quinone moiety and the methylene group at the triazole linker, while the betulin regions are neutral ( Figure 6 and Figure S5). its interaction with a biological target through hydrogen bonds and hydrophobic interactions. The distribution of the positive and negative charges of molecule is described by the molecular electrostatic potential map (MEP) [64,65]. The different charges are represented by different colors, which means that the red, blue, and green areas are negative, positive, or neutral, respectively [66]. The MEP's are sketched for an order of (−47,935) to 47,935 kcal/mol. The maps for hybrids 1, 5, 9, and 13 present in Figure 6, while for the rest of the compounds, they are shown in Figure S5. For all hybrids 1-16, the nucleophilic regions (red color) are localized in four main area. The first and second areas are localized on 1,4-quinone moiety. The first area contains the carbonyl atom at C-5 and the second area the carbonyl group at C-8, and the nitrogen atom. The third area includes the triazole linker, and the fourth, the substituent at the C-3 position of the betulin moiety. The electrophilic regions (blue color) are localized near the methine group at 1,4-quinone moiety and the methylene group at the triazole linker, while the betulin regions are neutral (Figures 6 and S5).
In each area, the local minima have been determined for 1-16, and they are collected in Table S4. The analysis of the local minima in the first, third, and fourth areas show that the arrangement of charges does not depend on the type of 1,4-quinone moiety. However, in the second area, the amount of potential minima depends on the type of 1,4-quinone moiety, i.e., the hybrids with 5,8-quinolinedione or 5,8-isoquinolinedione moiety have two potential minima, but compounds with 1,4-naphthoquinone moiety have only one potential minima in this area (Table S4).

Molecular Docking Study
Since the end of 2019, coronavirus disease 2019 (COVID-19), caused by SARS-CoV-2, has infected many people around the word. The World Health Organization (WHO) estimated that up to April 2021 more than 142 million people tested positive for COVID-19 and more than 3 million people died from this virus [67]. At the beginning of the pandemic, attempts were made to treat with chloroquine, but this treatment was halted and now chloroquine is not recommended (Figure 6) [68]. Contemporary reports indicate the effective use of amantadine and remdesivir in the treatment of patients with SARS-CoV-2 ( Figure 7) [69,70].
Continuing our research on molecular docking to SARS-CoV-2 targets, we examined the interaction between hybrids 1-16 and Mpro and PLpro proteins using the AutoDock Vina program (referred to as Vina). As reference ligands, chloroquine, amantadine, and remdesivir were used (Figure 7). In each area, the local minima have been determined for 1-16, and they are collected in Table S4. The analysis of the local minima in the first, third, and fourth areas show that the arrangement of charges does not depend on the type of 1,4-quinone moiety. However, in the second area, the amount of potential minima depends on the type of 1,4-quinone moiety, i.e., the hybrids with 5,8-quinolinedione or 5,8-isoquinolinedione moiety have two potential minima, but compounds with 1,4-naphthoquinone moiety have only one potential minima in this area (Table S4).

Molecular Docking Study
Since the end of 2019, coronavirus disease 2019 (COVID-19), caused by SARS-CoV-2, has infected many people around the word. The World Health Organization (WHO) estimated that up to April 2021 more than 142 million people tested positive for COVID-19 and more than 3 million people died from this virus [67]. At the beginning of the pandemic, attempts were made to treat with chloroquine, but this treatment was halted and now chloroquine is not recommended (Figure 6) [68]. Contemporary reports indicate the effective use of amantadine and remdesivir in the treatment of patients with SARS-CoV-2 ( Figure 7) [69,70]. For all hybrids 1-16, the nucleophilic regions (red color) are localized in four main area. The first and second areas are localized on 1,4-quinone moiety. The first area contains the carbonyl atom at C-5 and the second area the carbonyl group at C-8, and the nitrogen atom. The third area includes the triazole linker, and the fourth, the substituent at the C-3 position of the betulin moiety. The electrophilic regions (blue color) are localized near the methine group at 1,4-quinone moiety and the methylene group at the triazole linker, while the betulin regions are neutral (Figures 6 and S5).
In each area, the local minima have been determined for 1-16, and they are collected in Table S4. The analysis of the local minima in the first, third, and fourth areas show that the arrangement of charges does not depend on the type of 1,4-quinone moiety. However, in the second area, the amount of potential minima depends on the type of 1,4-quinone moiety, i.e., the hybrids with 5,8-quinolinedione or 5,8-isoquinolinedione moiety have two potential minima, but compounds with 1,4-naphthoquinone moiety have only one potential minima in this area (Table S4).

Molecular Docking Study
Since the end of 2019, coronavirus disease 2019 (COVID-19), caused by SARS-CoV-2, has infected many people around the word. The World Health Organization (WHO) estimated that up to April 2021 more than 142 million people tested positive for COVID-19 and more than 3 million people died from this virus [67]. At the beginning of the pandemic, attempts were made to treat with chloroquine, but this treatment was halted and now chloroquine is not recommended ( Figure 6) [68]. Contemporary reports indicate the effective use of amantadine and remdesivir in the treatment of patients with SARS-CoV-2 ( Figure 7) [69,70].
Continuing our research on molecular docking to SARS-CoV-2 targets, we examined the interaction between hybrids 1-16 and Mpro and PLpro proteins using the AutoDock Vina program (referred to as Vina). As reference ligands, chloroquine, amantadine, and remdesivir were used (Figure 7). Continuing our research on molecular docking to SARS-CoV-2 targets, we examined the interaction between hybrids 1-16 and Mpro and PLpro proteins using the AutoDock Vina program (referred to as Vina). As reference ligands, chloroquine, amantadine, and remdesivir were used (Figure 7). The results obtained with the use of the Vina program are presented in Table 8. The lower the ∆G energy, the better the affinity of the tested ligands for the receptor. Calculations of the K i from the binding energy of the pose generated by Vina were performed using Equations (6) and (7). Based on the obtained results, it can be concluded that the tested derivatives show docking scores in the range of −9.3 to −7.8 and −8.3 to −6.4 for Mpro and PLpro, respectively. The obtained ∆G are lower compared to the reference compounds (−7.4 to −4.5 and −5.7 to −4.1 for Mpro and PLpro, respectively). This indicates that in preliminary in silico studies, hybrids 1-16 show a higher affinity for the proteins used than the reference drugs. Comparing the score values of hybrids 1-16 with betulin-5,8-quinolinedione derivatives [18] showed that the introduction of triazole linker between betulin and 1,4-quinone increased the ∆G value for Mpro and PL pro proteins.
The main protease Mpro, also known as 3CLpro, is one of the coronavirus nonstructural proteins (Nsp5). Inhibition of Mpro would prevent the virus from replication and as a consequence, Mpro is one of the coronavirus proteins designated as potential targets for drug development. According to the crystallographic data, amino acids His41, Met49, and residues 164-168 of the long strand play an important role in stabilizing the ligand-Mpro complexes [71].
As seen in Figure 8a, the hybrids are localized in the hydrophobic matrix of the protein, but the arrangement depends on the type of the 1,4-quinone moiety. Comparing the score values for Mpro protein shows that in the series of 5,8-quinolinedione (1-4) and 5,8-isoquinolinedione (5)(6)(7)(8), the highest ∆G values exhibit derivatives with oxo group at the C-3 position of betulin moiety. In the group of 2-methyl-5,8-quinolinedione, better docking scores are obtained for hybrids 10 and 11. However, in the series of 1,4-naphthoquinone derivatives (13)(14)(15)(16), the group at the C-3 position of the betulin moiety does not influence the ∆G.
Ligands 2, 6, and 14, which have the best score values, contain the same betulin moiety but different 1,4-quinone fragments. In the group of 2-methyl-5,8-quinolinedione hybrids, the lowest ∆G has hybrid 11. For detailed analysis of molecular docking, complexes of the Mpro with 2, 6, 10, 11, and 14 ligands have been chosen. Complete models of the possible interaction in 3D and 2D views are presented in Figure 9a,e and Figure S6A-E. The detailed data about the type and length of the binding interactions between these ligands and the protein residues are summarized in Table S5.  Table S5.  The detailed data about the type and length of the binding interactions between these ligands and the protein residues are summarized in Table S5. Comparing the arrangement of hybrids 2, 6, 10 in the active site of the protein s that the position of nitrogen atoms in the 5,8-quinolinedione moiety influences the action between the ligand and the hydrophobic matrix of protein. As seen in Figur and Table S5, ligands 2 and 10 create the hydrogen bond and hydrophobic intera between 1,4-quinone moiety and His41, Asn142, Gly143, and Cys145. In these comp the betulin moiety is stabilized by hydrophobic interaction with Pro168. Changing th quinolinedione ring to 5,8-isoquinolinedione moiety causes different arrangement o ligand in the hydrophobic matrix. In complex Mpro-hybrid 6, the nitrogen atom and zole ring create the hydrogen bond and hydrophobic interaction with Pro168 and As Glu166, Leu141, respectively (Figures 9b and S76B). In this case, betulin is bound to M Met165, and His41 by hydrophobic interaction. The 1,4-naphthoquione moiety (hybr interacts with His41, Cys145, Met165, and Gln189 by hydrophobic interaction. More the triazole unit and betulin moiety are stabilized by hydrophobic interaction with C and Pro168, respectively (Figures 9e and S6E).
Figures 9d and S6D present the possible interaction of the best fitted compou inside the binding pocket of Mpro. Corresponding amino acids that are significant volved in the hydrophobic interactions between ligand and the betulin unit are as fol His41, Met49, and Met165. Moreover, the whole complex is additionally stabilized b hydrophobic interactions of the newly introduced substituents: the triazole unit Glu166 and the 1,4-quinone moiety with Leu167 and Pro168.
It should be emphasized that the 1,4-quinone moiety plays an important role i stabilization of Mpro hybrid complexes.
As seen in Figure 8b, the hybrids are localized in the hydrophobic matrix of the tein, but the type of 1,4-quinone slightly influences the arrangement of the hybrid i pocket site of proteins. Comparing the score values for PLpro protein shows that i group of 5,8-quinolinedione hybrids (1-4 and 9-12) the highest values have deriva with hydroxy group at the C-3 position of betulin moiety, while in the series of 5-13-16 the best results are obtained for hybrids with the oxo group at this position. Comparing the arrangement of hybrids 2, 6, 10 in the active site of the protein shows that the position of nitrogen atoms in the 5,8-quinolinedione moiety influences the interaction between the ligand and the hydrophobic matrix of protein. As seen in Figure 9a,c and Table S5, ligands 2 and 10 create the hydrogen bond and hydrophobic interaction between 1,4-quinone moiety and His41, Asn142, Gly143, and Cys145. In these complexes, the betulin moiety is stabilized by hydrophobic interaction with Pro168. Changing the 5,8-quinolinedione ring to 5,8-isoquinolinedione moiety causes different arrangement of the ligand in the hydrophobic matrix. In complex Mpro-hybrid 6, the nitrogen atom and triazole ring create the hydrogen bond and hydrophobic interaction with Pro168 and Asn142, Glu166, Leu141, respectively (Figure 9b and Figure S7B). In this case, betulin is bound to Met49, Met165, and His41 by hydrophobic interaction. The 1,4-naphthoquione moiety (hybrid 14) interacts with His41, Cys145, Met165, and Gln189 by hydrophobic interaction. Moreover, the triazole unit and betulin moiety are stabilized by hydrophobic interaction with Cys145 and Pro168, respectively (Figure 9e and Figure S6E). Figure 9d and Figure S6D present the possible interaction of the best fitted compound 11 inside the binding pocket of Mpro. Corresponding amino acids that are significantly involved in the hydrophobic interactions between ligand and the betulin unit are as follows: His41, Met49, and Met165. Moreover, the whole complex is additionally stabilized by the hydrophobic interactions of the newly introduced substituents: the triazole unit with Glu166 and the 1,4-quinone moiety with Leu167 and Pro168.
It should be emphasized that the 1,4-quinone moiety plays an important role in the stabilization of Mpro hybrid complexes.
As seen in Figure 8b, the hybrids are localized in the hydrophobic matrix of the protein, but the type of 1,4-quinone slightly influences the arrangement of the hybrid in the pocket site of proteins. Comparing the score values for PLpro protein shows that in the group of 5,8-quinolinedione hybrids (1-4 and 9-12) the highest values have derivatives with hydroxy group at the C-3 position of betulin moiety, while in the series of 5-8 and 13-16 the best results are obtained for hybrids with the oxo group at this position.
In each group of 1,4-quinone, we chose one of the best ∆G, and performed detailed analysis for it. Complete models of the possible interaction in 3D and 2D views are presented in Figure 10a-d and Figure S7A-D. The detailed data about the type and length of the binding interactions between these ligands and the protein residues are summarized in Table S5. In each group of 1,4-quinone, we chose one of the best ΔG, and performed detailed analysis for it. Complete models of the possible interaction in 3D and 2D views are presented in Figures 10a-d and S7A-D. The detailed data about the type and length of the binding interactions between these ligands and the protein residues are summarized in Table S5. Comparing the optimal docking post of compounds 1, 6, and 9 shows a general trend that the betulin moiety is bound to Pro248 and Lue162 by hydrophobic interaction (Figures 10a-c and S7A-C). The hydroxyl or carbonyl group at the C-3 position creates a hydrogen bond with Thr301 (1) or Arg166 (6 and 9). The 1,4-quinone moiety and triazole ring are involved in hydrophobic interaction with Pro248 and Tyr268, respectively.
Hybrid with 1,4-naphtoquinone 14 has a different arrangement in the active site of PLpro than these with 5,8-quinolinedione or 5,8-isoquinolinedione ligands (Figures 10d  and S7D). The 1,4-quinone moiety creates two hydrogen bonds with Lys157 and one with Gly163. The complex is stabilized by hydrophobic interaction with Leu162, Pro248, Met208, and Tyr264. In this case, the triazole ring does not interact with the hydrophobic matrix of enzyme.
As indicated by literature date, the molecular, physicochemical, and pharmacokinetic properties can be used to determine the biological activity of the compounds [73,74]. The MLR analysis allows to find correlation between score values (∆G), experimental lipophilicity and in silico determined parameters. The MLR equation (8) shows the correlation between the score value for Mpro, experimental lipophilicity (logPTLC), number of Comparing the optimal docking post of compounds 1, 6, and 9 shows a general trend that the betulin moiety is bound to Pro248 and Lue162 by hydrophobic interaction (Figure 10a-c and Figure S7A-C). The hydroxyl or carbonyl group at the C-3 position creates a hydrogen bond with Thr301 (1) or Arg166 (6 and 9). The 1,4-quinone moiety and triazole ring are involved in hydrophobic interaction with Pro248 and Tyr268, respectively.
Hybrid with 1,4-naphtoquinone 14 has a different arrangement in the active site of PLpro than these with 5,8-quinolinedione or 5,8-isoquinolinedione ligands (Figure 10d and Figure S7D). The 1,4-quinone moiety creates two hydrogen bonds with Lys157 and one with Gly163. The complex is stabilized by hydrophobic interaction with Leu162, Pro248, Met208, and Tyr264. In this case, the triazole ring does not interact with the hydrophobic matrix of enzyme.
As indicated by literature date, the molecular, physicochemical, and pharmacokinetic properties can be used to determine the biological activity of the compounds [73,74]. The MLR analysis allows to find correlation between score values (∆G), experimental lipophilicity and in silico determined parameters. The MLR equation (8) shows the correlation between the score value for Mpro, experimental lipophilicity (logP TLC ), number of rotat-able bonds (RB), and the energy of LUMO (E LUMO ) orbitals. The ∆G for PLpro correlates with molecular mass (M), number of acceptors of hydrogen bond (HA), and energy of HOMO (E HOMO ) orbitals (Equation (9)).

∆G
The obtained results show that the physicochemical and molecular properties of tested hybrids influence their interaction with the SARS-CoV-2 proteins.

Conclusions
In the presented study, the lipophilicity of betulin triazole derivatives with attached 1,4quinone moiety was determined and analyzed in terms of structure, physicochemical, and pharmacokinetic parameters. Comparing the lipophilicity of betulin triazole derivatives and hybrids with attached 1,4-quinone, the introduction of 1,4-quinone moiety reduces the lipophilicity. The cluster analysis showed a correlation between the molecular structure of hybrids and their experimental and calculated lipophilicity that can be used to predict these values when designing new compounds. The bioavailability of the tested compound was described by pharmacokinetic and physicochemical properties using the Lipinski and Veber rules. The obtained in silico parameters showed that most of the hybrids could be applied orally and they did not exhibit the neurotoxic activity.
For compounds with attached 1,4-quinone, the LUMO orbitals were mainly delocalized at the 1,4-quinone moiety. Localization of HOMO orbitals depended on the type of substituent at the C-3 position of the betulin moiety. The molecular properties depended on the type of 1,4-quinone moiety. The molecular electrostatic potential showed that the negative potential sites are present in the nucleophilic atoms, i.e., nitrogen and oxygen atoms. The arrangement of the charge depends on the type of 1,4-quinone moiety.
The molecular docking study showed that betulin triazole derivatives with attached 1,4-quinone can inhibit selected SARS-CoV-2 proteins. The high affinity score values of Mpro and PLpro proteins resulted from the introduction of the triazole ring as a linker connecting the betulin moiety with the 1,4-quinone fragment. The multilinear regression equation showed the correlation between the score values for Mpro and PLpro proteins, and experimental lipophilicity (logP TLC ), molecular descriptors, and global properties. The study showed that the determination of these parameters allows for the prediction of the interactions between the ligand and SARS-CoV-2 proteins.