**3.1. Atomistic simulations**

Within the fully atomistic treatment, only few nanocarriers could be completely simulated. This is the case, for instance, for cyclodextrins [30, 31], calixarenes [32], dendrimers of low numbers [33], and so on. Nevertheless, in most cases, in order to get information of the nanoparticles, further approaches should be done:


Is in this last approach where we will emphasis in this section due to the amount of work present in the literature [7, 38–43] and the direct connection with experimental model bilayers. Within atomistic simulations it is possible to explore, for instance, the drug partition in the different bilayer regions and go deeper into the main interactions responsible for their localization [37]. The drug localization could be guided by specific interactions, such as hydrogen bond [36], cation-π (for aromatic molecules) [44], salt bridges [45, 46], or entropic effects (hydrophobic drugs usually partitioned in the hydrophobic core [47].

Several biophysical techniques could be used to access bilayer properties. Between them, X-ray and neutron-scattering experiments allow us to obtain the electron density profiles that give us an idea on the overall organization of a lipid bilayer. Within the simulation not only it is possible to reproduce this profile but also access each system component contribution. The profiles are calculated dividing the simulation box in slices normal to the bilayer and conducting the time averaging of the net charge (equivalent to the total electron associated with each atom) per given thick slabs. Other quantities like radial distribution function could be used in order to obtain more details on specific interactions between the drug and the lipids or to extract information about water networks and compare with the structure factor of X-ray experiments [20].

In force fields, the potential energy could be divided into two parts, the intramolecular potential energy (the bonded part) and the intermolecular potential energy (the non-bonded part). Bonded atoms usually involve the stretching, the bending, and the torsion of the bonds. The non-bonded potential energy could be modeled as Lennard-Jones interactions and electrostatic interactions. These models are also applied between atoms of the same molecule, which are more than three bonds apart. Atoms that are three bonds apart (1–4 pairs) are treated dif-

The force field parameters are fitted to experimental and quantum mechanical data to match the spectroscopic, thermodynamic and crystallographic data of the molecule or, in the case of CG force field, they could also be derived from atomistic simulations. Among the classic effective potentials available are the atomistic force fields OPLS [23], CHARMM [24], AMBER [25], GROMOS [26] and CG Martini [21], SIRAH [27], and others [28]. The choice of the force field is a crucial step, and it should be carefully chosen since it strongly affects the results of an MD,

Within the fully atomistic treatment, only few nanocarriers could be completely simulated. This is the case, for instance, for cyclodextrins [30, 31], calixarenes [32], dendrimers of low numbers [33], and so on. Nevertheless, in most cases, in order to get information of the

• small models of the complete systems: small micelles, nanostructured particles [34];

Is in this last approach where we will emphasis in this section due to the amount of work present in the literature [7, 38–43] and the direct connection with experimental model bilayers. Within atomistic simulations it is possible to explore, for instance, the drug partition in the different bilayer regions and go deeper into the main interactions responsible for their localization [37]. The drug localization could be guided by specific interactions, such as hydrogen bond [36], cation-π (for aromatic molecules) [44], salt bridges [45, 46], or entropic effects

Several biophysical techniques could be used to access bilayer properties. Between them, X-ray and neutron-scattering experiments allow us to obtain the electron density profiles that give us an idea on the overall organization of a lipid bilayer. Within the simulation not only it is possible to reproduce this profile but also access each system component contribution.

• focus on single-drug interaction with specific components of the system [35];

• bilayers as a model of liposomes, polymersomes, niosomes, etc. [36, 37].

• internal structure of the nanoparticle in the presence of the drug;

(hydrophobic drugs usually partitioned in the hydrophobic core [47].

and the number of works are constantly checking the precision of them [29].

**3. What we can get from simulations?**

nanoparticles, further approaches should be done:

ferently in each force field.

78 Molecular Dynamics

**3.1. Atomistic simulations**

Local anesthetics (LA), pain-relief drugs, are among the most extensive drugs studied by simulations in their interaction with membranes. The structure and physicochemical features of each LA agent determine its potency and toxicity. Encapsulation into DDS -such as liposomes- can modulate the anesthetic effect (onset of action, duration of sensory block) and toxicity [5]. In this direction, simulations could guide the development of DDS devoted not only for short-term procedures (ambulatory and surgical interventions) but also for chronic pain treatment.

Most LAs have a pKa between 7 and 9; therefore, both neutral and protonated forms are relevant at physiological pH. Most works in literature overcome this issue studying the drug at a high and low pH, where all molecules are protonated or deprotonated, respectively. From MD simulation, most authors report neutral LA partition in the hydrophobic region of the bilayer and the protonated species in the lipid head/water interface and water phase, with differences between them (depending on hydrophilicity, steric factors, the lipid bilayer, etc.) [47–53]. For instance, the behavior described for lidocaine in 1,2-dimyristo yl-sn-glycero-3-phosphocholine (DMPC) is similar to the one reported by us, using prilocaine in 1-Palmitoyl-2-oleoylphosphatidylcholine (POPC), even using different force fields (GROMOS and CHARMM27, respectively) [48, 50]. In both cases the neutral drug presents a bimodal distribution. Also, crossing events, defined as the drug crossing from one monolayer to the other bilayer, were observed for both drugs. On the other hand, in our knowledge, no crossing events were reported for any protonated LA, even in long simulations (400ns, Slipid force field ) [53]. The differential distribution of prilocaine (PLC) depending on its ionization state is illustrated inside the bilayer as shown in **Figure 2A**. On the other hand, more hydrophobic LA at the neutral state, such as etidocaine, shows higher distribution than prilocaine in the bilayer core [50]. The localization of the etidocaine in this region promotes a higher disorder in the lipid chains. This property could be accessed through the calculation of the order parameter, defined by:

$$S\_{mol} = \frac{1}{2} \{ 3 \cos^2 \theta\_u - 1 \} \tag{3}$$

Here, *θ<sup>n</sup>* is the angle between the normal to the bilayer and the normal to the plane defined by two carbon-hydrogen bonds. The order parameter is related to the tilt angle of the chains and to the *trans-gauche* distribution of chain dihedrals. The effects on the lipid packing from guest molecules such as PLC distributed within membranes can be experimentally determined from changes in the carbon-deuterium segmental order parameter along the lipid chain. In that case, the experimental order parameter, *SCD* = -1/2*Smol*, is derived from the measured quadrupole splitting Δν = (3/4)(*e*<sup>2</sup> *qQ*/*h*)*SCD* [54, 55]. In **Figure 2B** we show the order parameter of POPC neat and containing etidocaine.

**3.2. Coarse grain**

encapsulation efficiency.

clinically effective dose.

account Henderson-Hasselbach equation:

*pH* <sup>=</sup> *pKa* <sup>+</sup> log [*neutral PLC*]

sent the molar fractions of neutral and protonated species, respectively.

Through the reduction of degrees of freedom, CG models are useful to efficiently simulate drug delivery systems, such as liposomes, polymersome, and micelles, relieving the size and time-scale limitations of atomistic simulation but losing in details. One of the crucial factors, representing the capacity of potential drug delivery systems, is the partition coefficient of a

In particular, the encapsulation of prilocaine into liposome was studied using this approach (MARTINI-based force field) [59]. Following the atomistic results, neutral PLC was fully encapsulated in the interior part of the lipid membrane where it adopts an asymmetric bimodal distribution. Our simulation results therefore suggest that although protonation leads to a structured interaction between drug and host, hydrophobicity is the major driving force of drug encapsulation. Our results also depend on the protonated PLC initial simulations conditions [59]. This observation suggested different preparation schemes of liposomedrug complexes, leading to prilocaine trapped within vesicles that could increase overall drug

On the other hand, it is demonstrated that the two LA species (neutral and protonated) are present at physiologic pH, contributing together to the anesthetic effect [50, 57] and they could be important in the development of DDS. For example, in a liposomal formulation prepared at pH 7.4, while the neutral form of the anesthetic is able to quickly and efficiently cross the bilayer, the protonated form is mainly found in the water phase. Both species depend on each other: the protonated needs the neutral one in order to efficiently cross the membrane, reaching the adjacent water compartment, while the uncharged species depends on the high solubility of the protonated LA in the water phase to reach the

The protonation/deprotonation reaction cannot be directly simulated using classical MD simulation. In order to partially overcome this, we represented the physiological pH, taking into

and the fact that PLC apparent pKa in membranes is 7.6 [60]; values around 0.4 and 0.6 repre-

In **Figure 3A**, we illustrated the results with a snapshot. In order to get an idea of how PLCs distribute inside the vesicle, we have calculated the average number of PLCs as a function of the PLC (center of mass) distance to the vesicle center. In **Figure 3**, we also show separated histograms for neutral (B) and protonated (C) PLCs. Considering that the average vesicle radius is 75Å, we can see from **Figure 3C** that all neutral PLCs are essentially found inside the vesicle and show two main peaks, in good agreement with previous results [47, 57]. On the other hand, just a small fraction (~14 molecules on an average) of the protonated

\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_\_ [*protonated PLC*] (4)

Molecular Dynamics Simulations to Study Drug Delivery Systems

http://dx.doi.org/10.5772/intechopen.75748

81

potential drug candidate between the aggregate and surrounding water.

**Figure 2.** (A) Electron density profiles (EDPs) of different bilayer components as a function of the membrane normal: POPC in black, water in blue and drug in red. Solid and dotted line are used for protonated and neutral prilocaine, respectively. (B) Calculated order parameter, −SCD, for plain bilayers (black); etidocaine at 1:6 and 1:3 drug:lipid ratio in green and violet, respectively. Circles and circles are used for saturated and unsaturated chains of POPC.

In a similar direction, simulations with the two bupivacaine (BVC) enantiomers helped to elucidate the higher cardiotoxicity observed for the R-form than the S-form [56]. The effect of R-BVCs (at 1:3 molar concentration) is to disorganize the membrane (decrease the order parameters); this effect is seen for both tails [51]. This is essentially related with the empty space in the lower part of the lipid tails due to the localization of the LA and the lateral expansion of the bilayer. On the other hand, S-BVC only promotes a soft disorganization for carbons above 10 (in both tails): in this case, the localization of the LA in the interior of the bilayer compensates the effects of the lateral expansion.

Thus, in our hands, molecular simulation of LA in model membranes was found useful to explain different aspects on the anesthesia mechanism and drug encapsulation [57]. Besides, within simulations, it was possible to explain the experimentally observed differences [49] between less (prilocaine/lidocaine) and more hydrophobic (etidocaine) isomers regarding the depth of their preferential insertion into bilayers [50], with possible implication on the increased potency and toxicity of the more hydrophobic analogs (Meyer-Overton rule) [5, 58].
