Computational Simulations to Guide Enzyme-Mediated Prodrug Activation

Prodrugs are designed to improve pharmaceutical/biopharmaceutical characteristics, pharmacokinetic/pharmacodynamic properties, site-specificity, and more. A crucial step in successful prodrug is its activation, which releases the active parent drug, exerting a therapeutic effect. Prodrug activation can be based on oxidation/reduction processes, or through enzyme-mediated hydrolysis, from oxidoreductases (i.e., Cytochrome P450) to hydrolytic enzymes (i.e., carboxylesterase). This study provides an overview of the novel in silico methods for the optimization of enzyme-mediated prodrug activation. Computational methods simulating enzyme-substrate binding can be simpler like molecular docking, or more complex, such as quantum mechanics (QM), molecular mechanics (MM), and free energy perturbation (FEP) methods such as molecular dynamics (MD). Examples for MD simulations used for elucidating the mechanism of prodrug (losartan, paclitaxel derivatives) metabolism via CYP450 enzyme are presented, as well as an MD simulation for optimizing linker length in phospholipid-based prodrugs. Molecular docking investigating quinazolinone prodrugs as substrates for alkaline phosphatase is also presented, as well as QM and MD simulations used for optimal fit of different prodrugs within the human carboxylesterase 1 catalytical site. Overall, high quality computational simulations may show good agreement with experimental results, and should be used early in the prodrug development process.


Introduction
Prodrug is a drug derivative intended to undergo chemical/enzymatic activation and thereby release an active parent drug, which is then free to achieve its pharmacological effect in the body [1]. It is designed to improve stability, biopharmaceutical, pharmacokinetic and pharmacodynamic features of the parent drug [2]. The prodrug approach can enhance drug formulation strategy and accomplish improved administration; it can also be used in achieving site specificity, improving the therapeutic drug effect and drug safety [3,4]. Prodrugs constitute a particularly useful approach for altering physicochemical drug properties (e.g., solubility) and for optimizing drug-like features of active compounds (i.e., absorption, distribution, metabolism and excretion (ADME)) [5,6]. The prodrug approach has been proven efficacious for several therapeutic areas, including anticancer, anti-inflammatory, antiviral agents and angiotensin-converting enzyme inhibitors (ACEI) [7].
The crucial step in the successful prodrug design is the inclusion of the activation mechanism that releases an active parent drug from a prodrug molecule in an efficient/controlled way, in order to achieve a successful therapeutic effect. Prodrug activation can be achieved through enzyme-mediated hydrolysis or oxidation/reduction processes or through chemical degradation within the body, triggered via a particular stimulus. The main principle of the prodrug approach is illustrated in Figure 1. This work focuses on the enzyme-mediated prodrug activation and optimization of this process through in silico methods. Numerous enzymes can be involved in the activation of prodrugs; some examples include oxidoreductases like CYP450, and hydrolytic enzymes such as carboxylesterase, butyrylcholinesterase, acetylcholinesterase, paraoxonase, β-glucuronidase, matrix metalloproteinase, alkaline phosphatase (ALP), phospholipase A 2 (PLA 2 ), human valacyclovirase and others [8,9]. Prodrugs often contain ester/amide functional groups derived from the hydroxyl, carboxyl or amine group of the parent drug; ester/amide is then activated through hydrolysis or oxidation and the active drug moiety is liberated. Such prodrugs were often designed to enhance oral drug absorption [10]. The main obstacle with these molecules is the lack of ability to foresee the bioconversion of the prodrug to the active parent drug, and thus its therapeutic effect. This can be overcome by using novel computational modeling techniques, which can improve the prodrug design of drug molecules with hydroxyl/phenol/amine functional groups [11]. In addition, several in silico techniques are employed in the optimization of the prodrug structure or prodrug-enzyme complex, in order to facilitate enzymatic cleavage. They can be based on simple empirical methods (i.e., molecular docking), or rather complex methods (i.e., quantum/molecular mechanics or free energy perturbation). Computational simulations provide an important insight into the prodrug design and overall fitting of the prodrug to the enzyme catalytical site, therefore they should be conducted prior to experimental studies, in order to reduce cost and improve the efficiency of the prodrug development process. The following section provides a brief overview of the calculations and methods used for the optimization of the prodrug activation process.

In Silico Methods for Predicting Prodrug Activation
In silico methods that can predict binding affinity between an enzyme and a substrate can be simple, empirical methods, such as molecular docking, or more complex methods, based on the laws of physics, such as quantum mechanics/molecular mechanics (QM/MM) or free energy perturbation (FEP) [12,13]. Docking studies are a rapid, routine way for molecule selection and design, however their predicting precision is rather low; whereas QM/MM, as well as FEP methods are time/effort/money consuming, but are more reliable [14,15]. Docking studies have a lower level of accuracy due to limited molecule databases, poor docking position choice, unfitting target binding site etc. However, it is the most commonly used method in fast preliminary analysis in prodrug design. QM calculations, including ab initio, semi-empirical and density functional theory (DFT), as well as MM, are being gradually employed as tools to offer structure-energy calculations for the prediction of potential drugs/prodrugs [16]. Ab initio molecular orbital method is a QM method founded on Schrödinger equation with certain estimates; it is used in limited systems (not more than 30 atoms) due to costly computing time; it is a method to calculate the electronic distribution and other features of molecules [17,18]. Semi-empirical calculations, such as AM1, PM3, MINDO, MNDO, MINDO/3 and SAM1, are grounded on the Schrödinger equation, and they include certain features that fit the experimental results; these methods can be used for studying compounds with more than 50 atoms and thus are extensively utilized in QM. DFT is also a semi-empirical QM method that can provide structural/energy calculations for intermediate-size biological systems; in the last 30 years, it has been widely used in predicting and explaining the features of various materials [19]. MM is used for the computation of various physical properties (structure, energy, dipole moment) and is utilized in evaluating large biological systems such as proteins, but is restricted by numerous torsion angle computations in structurally different molecules. Free energy perturbation methods, such as molecular dynamics (MD), are a form of in silico simulation in which atoms/molecules are permitted to interact for a certain time period through estimates of known physics, giving an outlook of the motion of the atoms; nowadays, it is an extensively used method for drug/prodrug design [20]. MD is widely used method for gathering structural/functional information about various macromolecules. It can model proteins and biological macromolecules motion, usually using the laws of classical dynamics. MD simulations can be performed together with MM or with QM/MM techniques [21]. Some potential uses include the prediction of the Michaelis-Menten complex formation between the enzyme and the prodrug, entry of water molecule into the active site of a particular enzyme, and evaluation of the steric hindrance between the prodrug and the enzyme [20]. All these computational approaches may be used for an intelligent design of novel prodrugs.

Computational Optimization of Enzyme-Mediated Prodrug Activation
Prodrug activation can be mediated through different enzymes, from oxidoreductases like Cytochrome P450 to hydrolytic enzymes such as carboxylesterase, phospholipase A 2 , β-glucuronidase, alkaline phosphatase and others. Numerous enzymes exhibit significant conformational modifications during the reaction cycle, and the function/relationship of such modifications to the chemical reaction steps should be investigated. In addition, enzymes have complex dynamics, with a range of internal motions, often including ones crucial to their function. This section provides an overview of particular enzyme-mediated prodrug activation and various in silico modeling techniques used for the optimization of the activation process.

Cytochrome P450 (CYP450)
There is substantial proof that the genetic polymorphisms of Cytochrome P450 (CYP450) enzymes contribute to the inconsistency in the activation of various prodrugs and hence to the effectiveness/safety of the drugs undergoing this bioactivation pathway [22]. The CYP450 enzymes mainly involved in the drug metabolism are CYP1A2, 2B6, 2C8, 2C9, 2C19, 2D6, 2E1, 3A4 and 3A5, which makes them of key interest for the activation of the prodrugs as well [23]. CYP3A4 is present in the hepatic tissues in the highest concentration, and is accountable for the oxidation of nearly two-thirds of all known drugs; nevertheless, other P450 are significant for drug metabolism and prodrug activation as well [24], and their expression is different in various tissues; therefore, they can be used to selectively target prodrugs to particular tissues.
CYP2C9 is responsible for the metabolism of approximately 15% of drugs on the market. Some variants of CYP2C9 show decreased enzymatic activity due to polymorphism, which can affect the clinical response of certain drugs. The uncommon allelic variant, CYP2C9.30 (which demonstrates < 50% activity when compared to wild type enzymes) is found in Japanese people and is linked to a decreased therapeutic response of the prodrug losartan [25]. MD simulations alongside energetic analyses were shown to be successful in elucidating the mechanism of amino acid (AA) replacement that influences drug metabolism/regioselectivity of CYP450 [26]. The CYP2C9 and wildtype/A477T models were solvated with TIP3P water molecules, and systems with~60,000 atoms were used. The simulations included CHARMM36 force field, and for losartan parameters, the DFT QM calculations were used. Structures of the substrate binding pocket, substrate docking, and channel dynamics were also determined using suitable software. The mutant enzyme showed higher rigidity of enzyme fragments involved in the recognition of the substrate, which could in turn affect the entrance and fitting of the substrate. MD simulations demonstrated that the A477T mutation present in CYP2C9.30 variant caused bigger rigidity of the key substrate recognition sites (SRS1 and SRS5) and the shifting of the β turn 4 of SRS6 toward helix F, thus decreasing the substrate access to some protein channels. The behavior of substrates, like losartan, was observed, and it was seen that the A477T mutant can modify the metabolism of some substrates, when compared to wildtype. The dynamics of substrate binding showed changes in substrate channel access and fitting to the active site, which could be responsible for altered catalytic enzymatic activity. Identified conformations of enzyme and substrates could be used for predicting suitable substrates of drug-drug interactions [27].
Paclitaxel has limited solubility, which restricts its widespread use. Many paclitaxel prodrugs were designed in order to improve its aqueous solubility. A recent study presents a quantitative structure property relationship (QSPR) model, which could serve as a basis for novel paclitaxel prodrug design with improved solubility. Molecular docking and molecular dynamics simulation were then used for the metabolic study with CYP1A2 enzyme, which revealed that the substituent groups in the prodrug could indeed be metabolized by the CYP1A2 enzyme [29].

Carboxylesterase
Prodrugs with ester bonds are frequently hydrolyzed via human carboxylesterase (hCE). Two major hCE enzyme classes are present in humans, hCE1 and hCE2. hCE1 has smaller catalytic site, and thereby cleaves smaller substrates, whereas hCE2 is reported to have a bigger catalytic domain [30]. However, both forms of the enzyme contain a catalytical site that is characterized by serine, histidine and glutamine; hydroxyl group in the serine can attack the carbonyl group of the ester prodrug, and the remaining 2 AA stabilize the complex. These molecular targets were employed in studying novel prodrugs of epalrestat and natural antioxidant products for diabetes complications [31]. MD and QM techniques were used in the process of prodrug design leading to the optimal accommodation within the catalytical site of hCE1. Most of the given prodrugs are monocyclic phenols, thus, due to the size of its catalytical site, the hCE1 enzyme was selected for analysis.
The main process in the hCE1 catalytical activation is the nucleophilic attack of serine hydroxyl group on the carbonyl ester group of the prodrug (based on frontier molecular orbital approach). Calculations including geometrical parameters (Burgi-Dunitz angle governing the nucleophilic attack and distance) and QM factors, such as highest occupied molecular orbital-lowest unoccupied molecular orbital (HOMO-LUMO) energy gap, were performed to determine the prospect of hydrolytic cleavage of the ester prodrugs via esterase ( Figure 3A). HOMO-LUMO difference can provide an insight into the mechanism of transfer of electrons from protein to the prodrug within the catalytic domain. The HOMO and LUMO calculations were performed using Schrödinger software, and DFT analyses (with local density approximation functional) were engaged to calculate the HOMO-LUMO energy for all structures. A gap in the energy difference within the HOMO and LUMO can be a certain stability index; high gap means high stability ( Figure 3B), and low chemical reactivity, and vice versa [32]. The MD simulation of designed prodrugs showed that the monocyclic prodrugs are more suitable substrate for hCE than bicyclic antioxidants, which was confirmed via QM as well. Governed by in silico predictions, prodrugs were synthesized and evaluated in the in vivo system; bearing in mind the in silico predictions, the most promising prodrug was evaluated for its in vivo antioxidant activity [31].

Phospholipase A 2 (PLA 2 )
Phospholipase A 2 (PLA 2 ) is overexpressed in many inflammatory conditions and cancer [33][34][35]. PLA 2 is an enzyme responsible for hydrolyzing the sn-2 bond of the phospholipid (PL) [36]. Our group studies PL-prodrugs consisting of the drug covalently bound to the sn-2 position of the PL. We aim to exploit the PLA 2 enzyme in order to hydrolyze the sn-2 bond of the prodrug, thereby releasing the drug specifically at the site of action, where the enzyme PLA 2 is overexpressed ( Figure 4A) [37]. Direct conjugation between the PL and the drug demonstrated lack of activation by PLA 2 [38], therefore, our PL-prodrugs contain different linker lengths between the PL backbone and the drug in the sn-2 position of the PL; different linker length resulted in different rate and extent of activation [8]. We developed an MD simulation in order to determine energetic changes of the PL-drug conjugates in the active site of PLA 2 [39][40][41][42]. An example of the PL-diclofenac prodrug structure is presented in Figure 4B. PLA 2 -mediated activation will occur only if the prodrug takes on a well-defined transition state geometry of the enzyme active site, characterized by interactions of the sn-2 carbonyl oxygen with the calcium atom from the PLA 2 active site, together with the protein His residue that activates a water molecule, used for the nucleophilic attack on the acyl group of the prodrug. The PL-drug concentration in the PLA 2 transition state geometry and kinetics of the cleavage are evaluated by the binding free energy of the prodrug within the enzyme active site. Lower binding free energy correlates with a higher degree of binding between the PL-drug conjugate and the enzyme. The differences in the transition state binding energy in the PLA 2 was calculated through decreasing/increasing the linker length. We used the thermodynamic integration (TI) method to compute the free energy of eliminating the -CH 2 units from the linker ( Figure 4C). TI calculates the free energy difference among the two states by ensemble-averaging the enthalpy variation along the path that connects the two states. The free energy of attaching the drug to the shorter linker was calculated using the umbrella sampling (US) and weighted histogram analysis methods (WHAM). Prodrugs containing diclofenac [43], indomethacin [44][45][46], and methotrexate [42] were evaluated through this simulation [39,40,42]. Optimal linker length was evaluated, and confirmed through in vitro [43,44] and in vivo [44] studies, with good correlation to in silico studies ( Figure 4D). Furthermore, MD and DFT simulations were performed on thio-ester pro anticancer ether lipid (S-ProAEL) and the mechanisms behind the PLA 2 -mediated activation were explained. It was shown that the presence of the water molecule has a significant part in the activation of different derivatives of PL-prodrugs (as well as for thio-esters). Further studies showed the distinction between the hydrolysis rate of the S-ester vs. the natural O-ester; significantly different ester conformation and a longer space between the carbonyl and the water molecule in the S-ester when compared to O-ester [47].

Cholinesterase
In a recent study, novel prodrugs for Alzheimer's disease targeting 2 specific targets, dual-specificity tyrosine phosphorylation-regulated kinase 1A (DYRK1A) and cholinesterase (ChE), were designed and synthesized [48]. The ChE and DYRK1A inhibitors were connected through a carbonate link. It was shown that a prodrug is inactive against both DYRK1A and ChEs and exhibits potent inhibition of ChEs only, through the oxidized form. The evaluation of promising binding interaction that contributed to the inhibition of ChE was performed via molecular docking (using Autodock Vina software). Supporting the likely interaction of prodrug with the peripheral anionic site of the enzyme, molecular docking suggests that the benzothiazole can be π-stacked with the tryptophan 286 residue of the peripheral anionic site, with a possible π-sulfur interaction among the thiazole moiety of the ligand and the same residue. The inhibition of ChE by compound of the prodrug could be clarified via a pseudo-irreversible mechanism. The carbonate link installed between the two inhibitors permitted transient masking of the inhibitory activity of the drug inhibitor against DYRK1A [48].

Alkaline Phosphatase
In an attempt to concentrate radioactive drugs inside the solid tumors, a novel technology is utilized: enzyme-mediated cancer imaging and therapy.
A radioactive quinazolinone prodrug, ammonium 2-(2 -phosphoryloxyphenyl)-6-[ 125 I]iodo-4-(3H)-quinazolinone ( 125 IQ 2-P ) is water-soluble, and can be hydrolyzed by alkaline phosphatase to a water-insoluble, 2-(2 -hydroxyphenyl)-6-[ 125 I]iodo-4-(3H)-quinazolinone ( 125 IQ 2-OH ). Upon enzymatic activation, the water-soluble molecule liberates its prosthetic group and the consequential molecule is water-insoluble and precipitates; this allows it to be permanently trapped within the extracellular spaces of targeted solid tumors. Two isoforms of the prodrug are present (IQ 2-P(I) and IQ 2-P ), and their structural alterations were studied via molecular modeling and docking methods, in order to elucidate the interaction among the prodrugs and the enzyme human placental alkaline phosphatase (PLAP), a membrane-bound hydrolytic enzyme of the tumor cell. The PLAP active site was well-defined using the AutoGrid module, and metal ions were modeled by Amber force field potentials (AutoDock). IQ 2-P(I) and IQ 2-P were docked into PLAP, and both binding free energy (∆G binding ) and inhibition constant (K i ) were calculated. Molecular docking results, following energy minimization and geometric optimization (until there was no resistance among ligand, metal ion, water molecules, and PLAP), revealed that only IQ 2-P accommodated the PLAP active binding site and interacts with the catalytic AA Ser92 that plays an important role in the hydrolytic process. The fact that IQ 2-P is a more suitable substrate was shown by the binding free energy (∆G binding ) of the isoforms to PLAP; the binding affinity of IQ 2-P was~170 fold stronger than that of IQ 2-P(I) . This was also confirmed through in vitro/in vivo studies [49].

Human Valacyclovirase
Valacyclovir prodrug was developed to overcome the low permeability of acyclovir. It was discovered that the conversion of valacyclovir to acyclovir occurs through enzyme valacyclovirase (VACVase), and that this enzyme presents an important therapeutic target for prodrug design [50]. For instance, VACVase is a suitable target for AA ester prodrugs; it is also known to hydrolytically activate antiviral nucleoside prodrugs, valganciclovir, producing AA, L-valine and corresponding parent drugs ganciclovir [51].
A homology model of VACVase was developed (via molecular operating environment software) to investigate the substrate specificity of different nucleoside analogues which offer a deeper understanding of catalytic and active site residues [52]. The docking of substrates and VACVase was performed with the MMFF94 force field with the solvation term by MOE-DOCK. It was shown that the model contains residues S122, H255 and D227 close to the catalytic triad, and charge-charge interaction sites, M52 or D123 for the R-amino group. Structural affinity of VACVase for hydrophobic amino acyl groups and restricted preference for the secondary alcohol substrates was shown. The reason for this might be attributable to the hydrophobic acyl-binding site consisting of I158, G161, I162, and L229 residues and the spatial limitation surrounding the active site by a loop from one side, serine and histidine on the other side, and L53 and L179 on the top. An additional study showed that VACVase has a specific binding mode and affinity for amino acid esters [53]; it contains a distinctive active site characterized by hydrophobic acyl pocket, a confined negative electrostatic potential, large open leaving group accommodating groove, and AA residue, Asp-123 ( Figure 5) [53] This is the first instance that a residue following nucleophile has a side chain heading for the substrate binding pocket and has an important part in substrate discrimination in serine hydrolases [53]. This study concluded that the enzyme functions as a specific α-amino acid ester hydrolase. The structure activity profile of the VACVase model is a useful tool for the future design of nucleoside prodrugs [54,55]. Another study showed that the rate of acyclovir prodrug hydrolysis could be evaluated based on the structural properties of the linker used [56]. It was shown that the rate-limiting step may be the creation of tetrahedral intermediate depending on the substitution choice. This insight could allow development of less hydrophilic prodrug system and controlled parent drug release. MM methods such as ab initio and DFT could provide better understanding of intramolecular processes, with the objective to elucidate the underlying parameters that influence the rate-determining step in prodrug activation. These methods have a crucial role in leading the intelligent design of novel prodrugs.

Discussion
Advancement in computing power and algorithms allows for the rapid development and innovation of the in silico tools that may advance our understanding of the biological systems. Molecular docking is a valuable tool for predicting binding modes for ligands, in particular protein targets, and can be used for target-based design of the prodrug (e.g., Sections 3.4 and 3.6); it is simple, relatively fast and has relatively low reliability [57]. Nevertheless, it cannot recognize relevant binding modes inside large and flexible binding sites, in cases when the enzyme dynamics should not be disregarded. This is where classical MD simulations, as well as accelerated MD and metadynamics, can help, accounting for the atomic descriptions of enzyme and prodrug dynamics, as well as the free energy of the entire conformational ensembles. These methods offer both structural and dynamic features of the enzyme (in a ≤ milisecond), which could predict likely binding sites, and propose the optimal enzyme-prodrug complex. The MD simulations are widely used and have short simulation times, however the absence of optimized tools and representation standards remains an obstacle. MM computation is focused on structures, energy, dipole moment and other physical properties, and does not take into account chemical reaction. QM methods calculate the distribution of electrons in molecules using the laws of quantum mechanics; they are extremely important in order to understand what happens in the chemical reaction (the change in the electron distribution in a reaction when bonds are broken and/or formed). For instance, nuclear tunneling is very significant in enzyme reactions that involve hydrogen transfer [58]. Nonetheless, it includes a relatively small number of atoms, which limits the understanding of the whole biological systems [59]. A potent way of overcoming these limitations is the use of QM/MM calculations, in which the system is divided into QM and MM regions, where QM regions deal with the electrons in a small area (i.e., the active site of an enzyme), and a larger part of the protein, the solvent and any other groups connected to the protein are modeled by MM. In order to include the environmental effects on the enzyme reaction, the QM and MM regions ought to interact, meaning that the QM region should be under the influence of the MM atoms. For instance, this method was shown to be successful in elucidating the catalytic function of active-site residues such as conserved proline in flavin-dependent monooxygenases [60]. It is sometimes said the QM/MM are hybrid simulations, in the sense that they are a combination of methods [61]. MM methods are also used in MD simulations and docking. Altogether, through the combination of several modeling methods, better prediction of enzyme-prodrug interaction may be obtained than the use of one technique alone. An overview and limitations of the key computational approaches used in the enzyme-mediated prodrug activation are presented in Table 1. Table 1. Comparison of computational methods used in enzyme/prodrug simulations.

Mechanics Molecular Mechanics Free Energy Perturbation
Purpose Predicting ligand-target interactions at a molecular level [62] Study reactions in molecular systems [63] Protein structure and dynamics (not applicable to chemical reactions) [59] Computing free energy differences between an enzyme and a substrate [64] Calculation

Limitations
Low reliability, cannot be used when protein dynamics should not be ignored [62] Difficulty in understanding large biological system (used for up to a few hundred atoms); high computational time and resource [63] Large number of unique torsion angles present in different molecules [65] High computational costs; absence of optimized tools and representation standards [20] Numerous prodrugs have been designed and developed aiming to improve the effectiveness of certain drugs and to overcome pharmaceutical/pharmacokinetic obstructions that occur upon drug administration (e.g., low stability, poor absorption, low solubility, reduced patient compliance) [5,66]. Many endogenous enzymes can be used for prodrug activation. In addition, recent advancement in the delivering of prodrugs to the specific tissue/organ employed non-endogenous enzymes for the site-specific activation of prodrugs: directed enzyme prodrug therapy (DEPT) [67]. However, prodrug development process in both cases is costly in terms of time, money and effort. The above mentioned in silico methods may provide trustworthy predictions, reducing the time and experimentation attempts, and making the overall process more affordable.
One of the main uses of improving the prodrug activation process is the optimization of the prodrug structure. Alongside PL-diclofenac, indomethacin and methotrexate, the use of different linkers/spacers/side chains was studied through different in silico methods; some examples include acyclovir prodrugs [56], PL-retinoid prodrug [68], anticancer ether lipids [47,69] and atovaquone [70]. The use of such groups allows the liberation of the parent drug at different, adjustable rates, depending upon the structural properties of the group itself. This approach, alongside predicting the molecular target (enzyme) for the prodrug, highlights the importance of the modern computational approached in the prodrug design. These computational methods reduce the need for numerous chemical synthesis and a long prodrug development process. In many cases, there was high correlation with experimental results (i.e., Figure 4D) [40,47,71]. A mechanism of prodrug intramolecular processes could also be used for the design of novel prodrug linkers. MM, DFT and ab initio methods were used, aiming to allocate the parameters that are influencing the rate-determining step and leading the process rate [16,70].
The molecular revolution has changed the process of prodrug development significantly. While the traditional prodrug approach focused on physicochemical drug features (e.g., charge, lipophilicity), the modern prodrug approach takes into account the molecular and cellular parameters, such as protein (e.g., enzymes, transporters) expression and distribution. Reducing the empirical features by using the targeted prodrug approach, while optimizing the activation process through in silico studies, endorses an effective overall process, as the outcomes may be considerably more predictable. Additionally, future progress should be made in the field of organic chemistry reactions that may be used in the development of certain prodrugs. Excellent understanding of the intramolecular reactions of prodrugs is essential for the advancement of this field as well. Since molecular docking, computational simulations and physical methods have their own limitations, a combination of physical methods and computational simulations yields more precise, biologically relevant results [72]. For instance, molecular docking and MD simulations come across some difficulties when it comes to predicting transient states in the enzyme reactions; the combination with QM methods (such as DFT) can help to overcome these limitations. In the coming years, an even wider use of computational approach is expected; different molecular orbital methods such as DFT and ab initio, as well as MM and MD simulations, provide a steppingstone for modern prodrug design.

Conclusions
The recent development of computational power allows powerful in silico methods to be used in the evaluation of the enzyme activity, catalysis and substrate specificity. Excellent knowledge of enzyme structure can help in evaluating enzyme-mediated prodrug activation, by offering new information regarding the structural requirements of prodrugs in order to accomplish optimal enzyme activation. QM, MM and MD simulations may provide an exceptionally thorough analysis of the main enzyme barriers, and understanding of the parameters needed for successful prodrug discovery/development. Prospective modeling studies should be performed prior to the prodrug synthesis and in vitro evaluation, in order to decrease the number of experimentations and save resources. Combined approaches, including novel in silico simulations of prodrug/enzyme complex, with rational prodrug design, present an important direction for the future development of the prodrug approach.

Conflicts of Interest:
The authors declare no conflict of interest.