Volume 13 - Issue 4

Review Article Biomedical Science and Research Biomedical Science and Research CC by Creative Commons, CC-BY

Significance of MD Simulation in Pharmaceutical Sciences: A Review

*Corresponding author: S M Shakib, Faculty of Pharmacy, Benazir Bhutto Shaheed University Liyari, Karachi, Pakistan.

Received: June 23, 2021; Published: July 21, 2021

DOI: 10.34297/AJBSR.2021.13.001895


Experimental methods for the analysis and characterization of pharmaceutical product, drug or any chemical substance, either at research level or commercial, are always time consuming, costly and require laborious efforts. In-silico characterization by molecular dynamics simulation technique can be utilize to analyze the physicochemical characteristic of any molecular system of pharmaceutical drug rather than just molecular systems of proteins or biological molecules. MD simulation also enables us to observe and visualize the interactions within given molecular system at atomic level and sub atomic level. The interpretations of simulation run predict any desire parameter like stability, thermodynamics and other characteristic on the basis of different bonding and non-bonding energies involve within the molecular system of drug. In recent years, different in-silico methods used to analyze the antiviral resistant mutations and discovery, the molecular motions and stability profile of antiviral drugs. Furthermore these techniques could be more beneficial in the field of pharmaceutical sciences and could enhance the precision of analysis and characterization of drugs and overcome the drawbacks of time consuming experimental methods. By MD simulation, the molecular motion and interaction of drugs with its environment can be visualized as 3D models that can never be done by in-vitro characterization methods.

Keywords: MD simulation; Molecular Mechanics; Quantum Mechanics; Pharmaceutical sciences; Computer-aided drug design; Force field; In silico modeling; Molecular modeling; Physicochemical properties; Physical stability


Computational Chemistry

The two most familiar branches of science are experiment and theory. In theory, claims have been made that might be falsifiable or predictions that can either be true or false on the basis of models represented while the experiments testify these predictions by accepting or rejecting true or false theories. Now the third pillar, which is currently emerging rapidly in academic research, is simulation. As we know, theories basically are conceptual ideas or equation and usually these equations are quite complexed and sometimes they are too complicated to solve with unobvious conclusions so instead of going towards expensive time consuming experimental methods, these theories can be assist by performing simulations, so basically computational chemistry [1-3] is the use of simulation in context of chemical systems or models. By using computational methods we simulate these complexed equations of theories and also sometime the experimental results generate some complex data that can be simulated into simplified conclusions by interpreting the experimental data. Sometimes there are different parameters in these equations of theories and by using experimental methods different values set for these parameters for simulation while giving solutions for these equations [4].

Molecular Mechanics and Force Fields

Molecular mechanics can be defined as group of models which implement empirical, algebraic, atomistic energy function for any molecular system. Here the empirical and algebraic energy function represents the experimental data, formulas, and equations and the atomistic energy function describe the level of details at individual atomic level for chemical system and biological molecules [5].

Here the total energy of the system consists of two components i.e. bonded energy and non-bonded energy [6].

Bonded energy included stretching energy between atoms, angle energy and torsion.

While non-bonded energy consists of electrostatic energy and Van der Waal energy [7,8].

Force Field

The reliance of energy of a molecular system on the coordinates of its atoms or particles is mathematically expressed as force field [9]. Force fields are basically the expression of parameterized empirical energy, involved in the interactions at atomic level of chemical systems containing large number of molecules or biomolecular systems [10,11]. In these systems evaluation of energy with high level calculations are very difficult, so instead of solving Schrodinger’s equation [12] by complex quantum mechanics (QM) calculations, molecular mechanic (MM) methods using force fields are appropriate. Force field can be implemented on either large systems of macromolecules like protein receptors [13,14], lipid bilayers [15,16] and nucleic acid [17,18] or smaller molecular system having interaction between few molecules or interaction of molecules with the solvents [19,20]. Depending upon the nature of molecular system force field can be varied. For larger molecular systems we have number of force fields such as AMBER [21], CHARMM [22,23] and GROMOS [24] and for smaller molecular system we have universal force field (UFF) [25], MMFF94 or MFF94x [26]. In any force field the potential energy between the particles is described in functional form and this potential energy functional form consist of two forms bonding potential [27] and non-bonding potential [28].

MMFF94 Force Field: Merck molecular force field (MMFF) was first developed by Merck research laboratories using computational data for its parameterization. The commencing published version of MMFF is MMFF94 [26]. In computational chemistry MMFF94 force field implemented or parameterized in MD simulation studies of gas phase small organic molecules [29]. It is applied in simulating the thermodynamic properties, structure, conformation and interactions of smaller molecules with limited number of atoms [30,31]. In some studies it has been reported that the analysis of proteins in organic solvents by parameterization of MMFF94 force field showed weak non-bonded interactions as compare to other force fields like AMBER or OPLS-AA [32]. Most widely MMFF94 force field is implemented in molecular modelling, i.e the computational method of characterizing, visualization and analysis of the structure [33] of pharmaceutical drugs [34] or other pharmacologically active compounds [35] or metal complexes [36].

OPLS Force Field: Optimized potentials for liquid simulations (OPLS) first developed by Prof William L. Jorgensen [41] is similar to AMBER force field and implanted in characterizing simulation studies of biomolecules [42] and optimized for liquid-state and gas phased parameters [43]. Depending upon the nature and type of atoms or molecules there are different set of parameters like OPLSUA (united atom) [44] and OPLS-AA (all atoms) [45].

CHARMM Force Field: Chemistry at Harvard macro molecular mechanics (CHARMM) is basically a set of different force fields that are developed by group of scientists working with Martin Karplus at Harvard [22]. It is implemented for all types of atoms and parameterized for simulation of proteins [46], lipids [47], RNA, DNA [23] and other biomolecules but not favorable for smaller organic molecules. It is reported that CHARMM is compatible diversified force field having potential energy functions for bond angle, torsion, dihedrals, out of plane motions and non-bonded electrostatic and Van der Waal interactions for different molecular systems [48].

Quantum Mechanics

In this section a brief introduction of quantum mechanics is discussed including theoretical background and different approaches of quantum mechanics and how it is differentiated from the classical mechanics.

Erwin Schrodinger in 1926 [49], expressed mathematical equation known as Schrodinger’s equation of wave function ᴪ or time dependent Schrodinger’s equation (TDSE) which described any system with respect to time.

This equation gives the idea about the wave nature of electrons, and it replaces the classical mechanics Newton’s equation F= ma which applies to Newtonian systems while the Schrodinger’s equation applies to quantum systems and describes the system three-dimensional wave function [50]. Another simplified form of Schrodinger’s equation is time independent Schrodinger’s equation (TISE) or eigenvalue equation.

In this linear partial differential equation ‘E’ is a constant represents energy level of system or in algebraic term it is the eigenvalue of corresponding Hamiltonian and ‘H’ is the Hamiltonian operator, define as a set of mathematical operators explaining all the interactions dictating the state of system which can be interpreted as the total energy of particle [50]. Total energy is associated with kinetic energy T and potential energy V.

Both the kinetic energy T and potential energy V further mathematically expressed as;

Here, in the system of electrons and neutrons ‘m’ represents the mass of ith particle and ‘h’ is Planck’s constant.

Here the potential energy is the columbic interaction of particles within the single atom, without considering the external field.

Interpretations and Approximations of Schrodinger’s Equation: Schrodinger expressed this mathematical equation of wave function for hydrogen atom and this equation is applicable on single electron system like hydrogen but in terms of molecular level or different chemical systems it is too complicated to analyze the wave function by the derivatization of Schrodinger’s equation. To overcome this complication, different scientist proposed different approximations and interpretations in the Schrodinger’s equation and simplify the methods for quantum mechanics simulation of different chemical systems [51]. Some of these approximations include Born-Oppenheimer approximation [52], which neglects the electron- neutron attractive potential energy term and only consider the neutron-neutron repulsive potential energy term as a constant in a given molecular system. The Hatree-Fock method is another approach to solve the Schrodinger’s equation [53] and calculate the energy of molecular system considering spin orbital of electron as wave function. Many other methods and different approaches have been introduced by different scientists to give the solution of Schrodinger’s equation for different molecular systems at different levels. Some of these methods include multideterminant methods [54], pair methods, semi empirical methods, density function theory [55] etc.

Molecular Dynamics Simulation

Techniques of Simulation: The fundamental purpose for simulation is to analyze the molecular system with respect to the structure and interaction of concerned molecules. There is a wide range of simulation techniques but the two basic techniques include molecular dynamics (MD) and Monte Carlo (MC) simulation and all other techniques are hybrids of these two.

Molecular Dynamics: Molecular dynamics is simulation technique that interprets classical Newton’s equation of motion to analyze the trajectories, movements and interactions in a given molecular system modelled at atomic level [3,56-58]. These atoms and molecules are allowed to simulate for a fixed period of run time and their potential energies and forces within the molecules analyzed by interpretations of molecular mechanics force field. Molecular dynamics simulation gradually analyzes the solution of equation of motion and for any molecular system it can be written as;

By solving equation of motion we can calculate the forces f acting on each atom of molecular system and acceleration experienced by these atoms and these forces are derived from the potential energy U (rN) of these atoms represented as set of 3N coordinates.

MD Integration Algorithms: In any molecular system of MD simulation, the potential energy is the function of atomic position considering all three coordinates (3N). The Newton’s equation of motion cannot be solved analytically because of complexed potential energy function. In order to get the numerical solution of this equation, different algorithms have been proposed by integrating the equation of motion. Depending upon the nature of molecular system, the appropriate algorithm was chosen for MD simulation [59]. These algorithms consider appropriate for the simulation on the basis of its computational efficiency, the duration of run time it can provide for simulation and its ability to conserve energy and momentum of molecular system. Every integration algorithm suggests the position, velocities and acceleration can be derived from Taylor’s series [60] expansion:

In these mathematical expressions ‘v’ represents velocity, ‘r’ is position and ‘a’ is acceleration.

Verlet’s Algorithm: Verlet’s algorithm is derived as;

By the combination of above two equations we obtained;

In this algorithm, the position and acceleration at time ‘t’ is used to estimate the new positions at time (t+δt), without any use of explicit velocities. Verlet’s algorithm is beneficial in straight forward molecular systems with adequate need of energy. But the precision of this algorithm is not highly accurate [61,62].

The Leap-frog Algorithm:

in leap-frog algorithm both velocities and position can be estimated by leaping over each other, first the velocity ‘v’ is calculated at time (t+1/2 δt) then by using this velocity, position ‘r’ can be estimated at time (t+δt) [63]. This is the benefit of this algorithm that velocities can be estimated explicitly for any particle or atom in the MD simulation but the drawback is that the velocity and position cannot be calculated explicitly [64]. Further, the velocity at any time ‘t’ can be evaluated by this mathematical expression;

Beeman’s Algorithm: This algorithm is derived from Verlet’s algorithm;

Beeman’s algorithm evaluated the velocity of particles more precisely by providing more accurate expression and conserve the energy as well. The only drawback for this algorithm is the computational cost for complex calculations [65].

Ensembles in MD Simulation: To perform MD simulation for different molecular systems we have different ensembles that represent different microscopic states for identical macroscopic or thermodynamic systems [66]. Ensembles of different characteristics are discussed here.

Micro-canonical Ensemble (NVE): NVE represents an isolated system having specified number of molecules N, stable volume V and fixed energy E of the system. Here the trajectories represent the exchange of potential energy and kinetic energy keeping the total energy constant. As the system in NVE ensemble is isolated, the temperature of the system is subjected to vary during the simulation and the system can be referred as adiabatic system with no exchange of heat so in NVE macromolecular systems of proteins, the temperature usually raised during simulation of exothermic conformation variations or interactions within the molecular system [67,68].

Canonical Ensemble (NVT): NVT is a set of molecular system which characterized the thermodynamic properties with constant no of atoms N, fixed volume V and fixed temperature T as well. Here an external heat reservoir is used to control the temperature of molecular system. There are set of different algorithms that are used to control the temperature of NVT ensemble system [68]. One of the commonly used algorithms is Berendsen algorithm [69], in this algorithm the temperature is controlled by scaling the velocity of each particle by a factor λ at every point where the exchange of heat energy is done by heat reservoir. Mathematically the factor λ is expressed as;

in this equation τ is the relaxation time, t is integration time step and T1 and T0 represents the target temperature and self-point temperature respectively. This algorithm is highly efficient incontrolling temperature of heat reservoir. There is also an alternative algorithm that used the coefficient fraction, which influence the forces acting on all particles, to control the temperature of heat reservoir. This algorithm is known as Nose- Hoover chain algorithm [70].

Isobaric-isothermal Ensemble (NPT): NPT ensemble characterized by constant number of atoms N, constant pressure P and constant temperature T. In this ensemble to control the pressure of the molecular system, a barostat is also needed in addition to the thermostat. NPT ensemble systems are usually significant in simulation of systems containing biological membranes [68,71].

Grand canonical ensemble (μVT): μVT ensemble characterized the thermodynamics of molecular system by considering constant chemical potential μ, constant volume V and constant temperature T [68].

Literature Review of MD Simulation of Pharmaceutical Products

Characterization of Metallic Drug Nanoparticles by MD Simulation

Theoretical approach of MD simulation, in companionship with experimental methods, is significant in the characterization of metallic nanoparticles at atomic level [72]. The potential of silver nanoparticles in the field of medical sciences is known for centuries, MD simulation is one of the advanced cost-effective techniques for researchers to characterize important features of silver nanoparticles. On March, 2020, Sohraby performed MD simulation of silver nanoparticles in order to theoretically estimate the potential of rosemary active compounds in coating of nanoparticles and analyze the dynamic behavior of these compounds with silver nanoparticles [73]. In 2018 and 2019, M.M. Blazhynska and his associates study the polarization effects of different force fields in MD simulation of silver nanoparticles and determine the morphological stability of different shapes of silver nanoparticle by analyzing the effect of temperature and binding energy by MD simulation [74,75]. In another study, complete melting, sintering and surface pre-melting point of silver nanoparticles was analyzed by MD simulation [76,77]. L Chen in 2017 perform the MD simulation of silver nanoparticles to characterized the physicochemical properties, stability of lattice structure and surface energy of NP [78]. In the research work of Paul Martin, he studied the nanotoxicology by the MD simulation of silver nanoparticles and analyzed that the reactivity of nanoparticles can be influence by the structure and size variations and by controlling the size, toxicity profile of nanoparticles can be managed [79].

Characterization of Antiviral Drugs by MD Simulation

In 2013, the antiviral resistance mutation was studied in HBeAg positive patients, administering the entecavir therapy. Molecular docking and MD simulation was performed on 3D model of HBV polymerase to analyze the antiviral resistant mutations [80]. In 2015, Dante Morgnanesi presented his work to explain the importance of computational chemistry and molecular dynamics simulation in the discovery of hepatitis B antiviral drugs. This study was focused on three computational techniques including homology modeling to analyze the structure and conformation of protein HBV polymerase, molecular docking to understand the antagonizing mechanisms of hepatitis B antiviral drugs and molecular dynamics simulation to study the molecular motions and stability profile of antiviral drugs [81]. Another computational study was done for the design and pharmacological analysis of antihepatitis B drugs, by molecular docking and molecular dynamics simulation. 500 picoseconds simulation at constant pressure was run to get the relative stability of antiviral drugs [82]. In 2014, insilico computational techniques were studied for the identification of novel HBV target NTCP, sodium taurocholate co‐transporting polypeptide. Molecular modelling and screening characterization of antiviral drugs was done by molecular dynamics simulation and experimental validation and proposed a significant target for anti HBV drugs [83].


In modern research, simulation is the bridge between the theoretical and experimental science. This technique made it possible to visualize and characterize the individual matter, its internal or free energies and its motion as a function of time. Now molecular dynamics simulation techniques has become a crucial tool in analysis, modelling and characterization of basic physical, structural and functional properties of biological molecules specially proteins, amino-acids or enzymes. This study elaborates the significance of molecular dynamics simulation in analysis of antiviral drugs, metallic nanoparticles and other pharmaceutical products and it helps to understand the nature of these drugs at the atomic level and this analysis could be used to establish a better pharmacokinetic and stability profile and to get the enhanced and safer drug delivery to its target site.


Sign up for Newsletter

Sign up for our newsletter to receive the latest updates. We respect your privacy and will never share your email address with anyone else.