Permeability of TB drugs through the mycolic acid monolayer: a tale of two force fields

Abstract

Tuberculosis (TB) treatment becomes challenging due to the unique cell wall structure of Mycobacterium tuberculosis (M. tb). Among various components of the M.tb cell wall, mycolic acid (MA) is of particular interest because it is speculated to exhibit extremely low permeability for most of the drug molecules, thus helping M.tb to survive against medical treatment. However, no quantitative assessment of the thermodynamic barrier encountered by various well-known TB drugs in the mycolic acid monolayer has been performed so far using computational tools. On this premise, our present work aims to probe the permeability of some first and second line TB drugs, namely ethambutol, ethionamide, and isoniazid, through the modelled mycolic acid monolayer, using molecular dynamics (MD) simulation with two sets of force field (FF) parameters, namely GROMOS 54A7-ATB (GROMOS) and CHARMM36 (CHARMM) FFs. Our findings indicate that both FFs provide consistent results in terms of the mode of drug–monolayer interactions but significantly differ in the drug permeability through the monolayer. The mycolic acid monolayer generally exhibited a higher free energy barrier of crossing with CHARMM FF, while with GROMOS FF, better stability of drug molecules on the monolayer surface was observed, which can be attributed to the greater electrostatic potential at the monolayer–water interface, found for the later. Although both the FF parameters predicted the highest resistance against ethambutol (permeability values of 8.40 × 10⁻³⁴ cm s⁻¹ and 9.61 × 10⁻³¹ cm s⁻¹ for the CHARMM FF and the GROMOS FF, respectively), results obtained using GROMOS were found to be consistent with the water solubility of drugs, suggesting it to be a slightly better FF for modelling drug–mycolic acid interactions. Therefore, this study enhances our understanding of TB drug permeability and highlights the potential of the GROMOS FF in simulating drug–mycolic acid interactions.

---

1. Introduction

The efficacy of many drugs against infectious diseases, whose targets reside within the bacterial cells, depends on their penetrative capacity through the bacterial membrane. Many infectious pathogens exhibit significant resistance against drug penetration, leading to drug resistance in conventional therapeutic routes. One prime example of such infections/diseases is tuberculosis (TB), the world's top infectious killer. TB is responsible for ∼10.6 million infections and around 1.5 million deaths annually.¹ The majority of TB infections and related deaths occurred in poor and developing countries like India, China, Bangladesh, Nigeria, and South Africa.² Among these countries, India holds the largest share of global TB infections.¹ In 2019, more than 2.4 million TB cases were found in India.³ Although TB mainly infects the lungs, infections can be found in other parts as well. TB is contagious and spreads through air when a patient with lung TB sneezes, coughs, or spits. Co-infections like HIV, emergence of multi-drug resistance TB, lack of healthcare infrastructure, and discontinuation of the treatment midway are some of the root causes behind the drastic global impact of TB.²

Mycobacterium tuberculosis (M.tb) is the causative agent of the infectious disease tuberculosis, and it belongs to the Gram-positive bacteria family, Mycobacteriaceae.⁴ It is distinct from most other Gram-positive bacteria because of the presence of an extra thick cell wall structure located outside the peptidoglycan layer found in Gram-positive bacteria.⁵ This extraordinarily thick cell wall structure is divided into several layers. The outermost layer consists mainly of proteins and gluten, with a small amount of lipids.⁶ The next layer is named the mycomembrane, which acts as the main permeability barrier. The inner leaflet of the mycomembrane contains tightly packed, parallelly arranged long chain fatty acids, namely mycolic acids (MA), covalently bonded to a polysaccharide called arabinogalactan (AG).⁷ The mycolic acid leaflet (monolayer) bestows M.tb with unique properties that defy medical treatment by lowering the efficacy of antibiotics/biocides by acting as a hydrophobic permeability barrier mainly for hydrophilic molecules, making the organism more resistant to chemical damage and dehydration.⁸ Mycolic acids are also responsible for the bacterium growth inside macrophages, effectively hiding it from the host immune system.

Mycobacterial mycolic acids possess some distinct characteristics. They are longer and contain 70–90 total carbon atoms, among which typically 24 carbon atoms belong to a fully saturated R-branch and the rest are part of the unsaturated mero chain.⁹ Usually two positions are there in the mero chain, which may be occupied by functional groups.⁹ The proximal position (near the α-hydroxy acid) contains exclusively cis- or trans-olefin or cyclopropane.¹⁰ However, the distal position may be the same as the proximal position or contain one of a variety of oxygen moieties such as R-methyl ketone, R-methyl methyl ether, methyl-branched ester, or R-methyl epoxide.⁹ The M.tb cell wall contains three kinds of mycolates: α-mycolates, methoxymycolates, and ketomycolates, with relative abundances of 51%, 36%, and 13%, respectively.⁹ The exact packing of the mycolic acids in the bacterium cell wall remains elusive. Previously, it was reported that the mycolic acids can fold into three distinct “W”, “U” and “Z” conformers,¹¹ although additional conformations have also been listed in the literature.¹²

Because of the important role of mycolic acid in the drug resistance of M.tb, it is of utmost importance to study the interactions of existing standard TB drugs or any potential TB drug candidate with the mycolic acid monolayers and to obtain a quantitative idea about the permeability of drug molecules through them. Machine learning based approaches have been utilized to determine the permeability probability of drug-like compounds through the M.tb cell wall, but those models do not provide the thermodynamic details of drug diffusion through the cell wall.^13–15 Previously, Hong et al. explored drug–mycolic acid interactions by computational means and calculated the diffusivity of several drug molecules inside the mycolic acid monolayer; they did not consider the anomalous diffusion of the drugs within the mycolic acid clusters, which usually takes place inside a tightly packed biomembrane.¹⁰ Moreover, although the permeability of water, oxygen, small organic molecules, and potential drug candidates through lipid bilayers has been widely explored through computational means together with the development of newer methods of permeability calculation,^16–25 no such investigation has been carried out on the mycolic acid monolayer. Against this background, we tried to compute the effective resistance and permeability of the mycolic acid monolayer against some of the well-known first- and second-line hydrophilic TB drugs, namely ethambutol, ethionamide, and isoniazid. For simplicity, we have modelled only α-mycolates, the abundant type found in the M.tb bacterial cell wall. It is noteworthy that the success of such modelling crucially depends on many factors like the generation of correct molecular structures, the level of structural detailing, and the choice of atom–atom interaction parameters or force field parameters. Apart from quantitative assessment of the permeability of the mycolic acid monolayer, the present study also aims to dissect the influences of some of the governing factors stated earlier on the outcomes of our in silico modelling. To achieve this goal, we have modelled the mycolic acid chain and the TB drugs using two levels of molecular detailing. In the first approach, we have explicitly considered all the atoms and implemented the widely used CHARMM36 all-atom (CHARMM) FF to describe the atomic interactions. On the other hand, we have also designed mycolic acid and drug molecules using the united atom approach, where hydrogen atoms bonded to carbons are modelled as a single bead. The GROMOS 54A7-ATB (GROMOS) FF was employed to describe the intra- and inter-atomic/molecular interactions. We have compared the outcomes of these two approaches and commented on the advantages and disadvantages of each approach. Summing up, the present study not only explores the permeability of the mycolic acid monolayer in a quantitative manner, but also delves into the effects of the choice of force field parameters and level of detailing on such studies.

2. Simulation details

2.1 Modelling of mycolic acid–drug systems

α-Mycolic acid was modelled using both all-atom and united atom representations. Packing of 100 α-mycolic acid molecules in a rectangular single monolayer was performed using MEMGEN and PACKMOL 20.3.1.^26,27 100 water molecules per mycolic acid chain were added to the simulation box, and the water molecules were placed above and below the monolayer. For all-atom representations, the systems were modelled using CHARMM all-atom FF parameters and a CHARMM-modified TIP3P water model.^28–30 The parameters were generated using CHARMM-GUI and CGENFF web interfaces.^31–35 For the united atom representations, GROMOS FF parameters, generated using an Automated Topology Builder and the SPC/E water model, were used to compute the molecular interactions.^36–38 The initial configurations of the systems together with the structures of the drug molecules used for the study are depicted in Fig. 1.

Fig. 1 Molecular structures of drug molecules and initial configurations of the drug–monolayer systems. Structures of tuberculosis drug molecules used in this study are shown in the upper-middle panel of the figure. The structure of α-mycolic acid is depicted in the middle. Initial configurations of the ethambutol–monolayer system is presented for (a) CHARMM36 and (b) GROMOS 54A7-ATB FFs. Starting configurations of the ethionamide–monolayer system is shown in (c) and (d) for the CHARMM36 FF and the GROMOS 54A7 FF, respectively. Snapshots of isoniazid–monolayer systems at 0 ns are shown for (e) CHARMM36 and (f) GROMOS 54A7 FFs. Colour code for mycolic acid: C, pale blue; O, red; H, white. Colour code for the drug molecule is based on atomic mass. For ethambutol/isoniazid: H, red; O, blue; N, light blue; C, whitish blue. For ethionamide: S, deep blue; H, red; C, reddish white; N, white. Water is presented using an ice-blue transparent surface.

The systems were then subjected to energy minimization in 200 000 steps using the steepest descent algorithm to remove any close contacts between atoms. The minimized system was subsequently subjected to a 10 ns equilibrium MD run in the NVT ensemble at 300 K, followed by a 100 ns run in the NPT ensemble at 1 bar pressure. The drug molecule (ethambutol/ethionamide/isoniazid) was then added to the equilibrated system, 1.0 nm above the surface of the MA monolayer, replacing the necessary number of water molecules. The drug molecules were placed on the side of the monolayer containing the –COOH groups of the α-mycolic acids (head region). The drug–monolayer systems were again energy minimized using the steepest descent algorithm to remove any close contacts. The minimized systems were then again equilibrated in the NVT ensemble for 1 ns (at 300 K), followed by a second phase of equilibration in the NPT ensemble (at 1 bar pressure) for 10 ns. Position restraints were put on the drug molecules during the equilibration phases. The equilibrated systems were then subjected to a production run of 250 ns and steered molecular dynamics simulation for either 1 ns (for all-atom representation) or 750 ps (for united atom representation). A modified Berendsen thermostat with a coupling constant of 0.1 ps was used for temperature coupling, while the pressure was maintained using a Parrinello–Rahman barostat with a coupling constant of 0.5 ps.^39,40 A semi-isotropic coupling scheme (pressure can vary independently in the z-direction and in the x/y-direction) was employed to account for the anisotropic structure of the mycolic acid molecules. Throughout the simulation time, all bonds involving hydrogen atoms were restrained using the LINCS algorithm.⁴¹ A simulation box of 6.32 × 6.32 × 12.51 nm³ was used for the systems with the CHARMM FF, while the dimension of the simulation box was 6.58 × 6.51 × 11.95 nm³ for the GROMOS FF. Periodic boundary conditions were employed in all three directions with a cut-off distance of 1.2 nm to compute the short-range LJ interactions and short-range part of the Coulomb interactions. The particle mesh Ewald (PME) summation technique was used to calculate long-range coulombic interactions.⁴² The integration time step used was 2 fs and the trajectories of the systems were saved at an interval of 10 ps for subsequent analysis. All simulations were performed using the GROMACS-2022.3 package and the visualization was done using VMD 1.9.3 software.^43,44 All of the analyses were carried out using gmx modules and in-house python scripts. We have studied a total of six drug–monolayer systems.

2.2 Steered MD simulations and umbrella sampling

Steered MD simulations were performed to pull the drug molecules through the monolayer with a constant velocity of 0.01 nm ps⁻¹. During pulling, a harmonic spring with a spring constant of 1000 kJ mol⁻¹ nm⁻² was applied on the drug molecule along the pulling direction, while the monolayer was taken as the static reference. The trajectory was recorded for every 0.1 ps. In the case of isoniazid using all-atom representation (CHARMM36 force field), a spring constant of 1500 kJ mol⁻¹ nm⁻² was used to keep the pulling motion steady and to ensure overlapping among the probability distributions of various windows during umbrella sampling.

Windows/configurations for the umbrella sampling were generated from the steered MD simulations at an interval of 0.1 nm along the z-direction. A 10 ns run was performed on each window to obtain the probability distribution of the systems for every configuration. For isoniazid and CHARMM FF parameters, an uneven sampling scheme (windows with 0.05 nm spacing for up to 1.5 nm distance from the COM of the monolayer and after that windows with 0.1 nm spacing were taken) was opted to ensure the overlap between the probability distributions of two consecutive sampling windows. The PMF profile was computed using the weighted histogram analysis method (WHAM) embedded in GROMACS, from the probability distribution obtained for each window.⁴⁵

We have used the umbrella sampling technique (US) to further compute the local resistance (R(z)), effective resistance (R_eff), and effective permeability (P_eff) of the membrane in the following manner. First, the position dependent diffusion coefficient (D(z)) was calculated for every umbrella sampling window using the following equation:^46,47

whereHere, var(z) corresponds to the variance of positional z-coordinates of the COM of the drug molecule and c_zz(t) is the positional autocorrelation function. The local resistance value R(z) for each umbrella sampling window was obtained using eqn (3) given below.^46,47Here, ΔG or the free energy difference was obtained from the PMF curve at different values of z.

Effective resistance (R_eff) and effective permeability (P_eff) of the monolayer were calculated by integrating eqn (3) and are given by eqn (4).

3. Results and discussion

In this section, thermodynamic free energy barriers of crossing the monolayer will be illustrated using both the CHARMM36 all-atom FF and the GROMOS 54A7-ATB united atom FF together with a description of the drug–MA monolayer interactions during the passage of drugs through the monolayer.

3.1 Drug–monolayer interactions during steered MD simulations

Different drugs interacted differently with the MA monolayer while being pulled through the latter. The various modes of the drug membrane interactions during pulling for two different force fields are shown in Fig. 2. With the CHARMM FF, the primary mode of drug–membrane interaction is the vdW interactions (Fig. 2(a)). The strength of vdW interactions was found to be the highest for ethambutol, followed by ethionamide and isoniazid (Fig. 2(a)), respectively. For all three drug molecules, the highest vdW interaction strength was recorded between distances of ∼−1 nm and ∼2 nm from the COM of the monolayer, i.e. around the middle portions of the mycolic acid assembly (Fig. 2(a)). The magnitude of vdW interaction energy decreased toward both ends of the mycolic acid layer (Fig. 2(a)).

Fig. 2 Different modes of drug–monolayer interaction during the pulling of drug molecules through the monolayer. Spatial variation of drug–membrane vdW interaction energy during the passage of drugs through the MA monolayer for (a) CHARMM and (b) GROMOS FFs. Positional dependency of electrostatic interactions for (c) CHARMM and (d) GROMOS FFs. In the figures, vertical black dashed lines represent positions of water–monolayer interfaces.

The nature of coulombic (electrostatic) interactions were noticeably different for different drug molecules (Fig. 2(c)). Electrostatic interaction of isoniazid was the highest in the head group region on the mycolic acid and became zero afterward (Fig. 2(c)) while pulling through the MA monolayer. Ethionamide showed a similar trend, but the strength of interaction was lower (Fig. 2(c)). Ethambutol interacted with the membrane more strongly than ethionamide and unlike the other two molecules, it electrostatically interacted with the middle region of the membrane as well (Fig. 2(c)). The hydrogen bond formation also followed a similar pattern with the CHARMM FF. Isoniazid formed the highest number of H-bonds with the head region of the MA monolayer but did not form any after that (Fig. S4(a) in the ESI†). Ethambutol formed H-bonds with the head and middle portion of the monolayer, while ethionamide formed H-bonds with every region of the membrane (Fig. S4(a) in the ESI†).

A similar trend in vdW interactions was observed with GROMOS parameters with the highest strength for ethambutol, followed by ethionamide and isoniazid (Fig. 2(b)). The coulombic interactions are also similar to CHARMM force fields, but for isoniazid, non-zero coulombic interaction is observed for every region of the MA monolayer (Fig. 2(d)). Moreover, isoniazid formed hydrogen bonds with both head and tail regions of MA with the GROMOS FF (Fig. S4(b) in the ESI†).

The strength of vdW interactions for each drug molecule was found to be very close for CHARMM and GROMOS FF parameters (Fig. 2(a) and (b)). Significant differences have been observed in terms of electrostatic interactions. The CHARMM FF offered stronger coulombic interactions for ethambutol for the head region of the monolayer and slightly less electrostatic interactions for ethionamide (Fig. 2(c)). The most significant difference is noticed for isoniazid, which has been mentioned earlier. No significant change in the number of hydrogen bonds was probed for ethambutol and ethionamide with the change of the force field parameters (Fig. S4 in ESI†). Notable changes in the number of hydrogen bonds appeared for isoniazid, as mentioned earlier.

3.2 Free energy barrier and membrane permeability

We have employed a 1-dimensional umbrella sampling (1D-US) technique to calculate the potential of mean force (PMF). It should be noted here that although 1D-US is one of the most widely used methods for free-energy calculation, it is reported that for some cases, it does not generate an appropriate description of the free energy of the systems.⁴⁸ Computationally challenging procedures like 2D Hamiltonian replica exchange umbrella sampling or 2D metadynamics can determine the free energy landscape more accurately and with respect to more degrees of freedom.⁴⁸ However, 1D-US is sufficient for the permeability calculation, and thus, we have stuck to 1D-US. The PMF obtained from 1D-US for all three drugs using the CHARMM FF is plotted in Fig. 3(a). From the figure, it is evident that the free energy barrier is the highest for ethambutol (41.30 kcal mol⁻¹), followed by ethionamide (39.39 kcal mol⁻¹) and isoniazid (20.07 kcal mol⁻¹). The free energy reached its minimum on the surface of the monolayer for all of the drug molecules, which, in corroboration with the other results, indicates that the drug molecules prefer to remain on the surface of the monolayer. The maxima of the free energy barrier for all three drugs were observed to be around the distance of 2 nm (Fig. 3(a)). It is worth stating that the order of the free energy barrier also followed the order of vdW interaction strength (Fig. 2(a) and 3(a)). Two peaks appeared in the free energy profile for every drug molecule, specifically for isoniazid because of the spatial change in the coulombic interactions and H-bond formation (Fig. 2(c) and Fig. S4(a) in ESI†), which are most prominently seen for isoniazid. Local minima can be seen at an ∼3 nm distance (Fig. 3(a)), which is at the vicinity of the lower surface of the monolayer. These minima appeared because of the preference of the drug molecules to attach to the membrane surface.

Fig. 3 Mycolic acid monolayer offers high free energy barrier for drug molecules. Potential of mean force (PMF) for various pharmaceutical molecules for (a) CHARMM and (b) GROMOS FFs. Free energy barrier was found to be higher for ethambutol and ethionamide for the CHARMM FF, compared to the GROMOS FF. Snapshots of different drug–monolayer systems have been depicted at various points of the PMF curve, with the colour code for mycolic acid: C, pale blue; O, red; H, white. Colour code for the drug molecule is based on atomic mass. For ethambutol: H, red; O, blue; N, light blue; C, whitish blue. For ethionamide: S, deep blue; H, red; C, reddish white; N, white. In the figures, vertical black dashed lines represent positions of water–monolayer interfaces.

The shape of the PMF profiles of the different drug molecules with GROMOS FF parameters is plotted in Fig. 3(b). The major differences between the GROMOS PMF profiles and the CHRAMM PMF profiles are the following:

(1) With GROMOS parameters, ethambutol exhibited the highest free energy barrier (37.79 kcal mol⁻¹) like CHARMM; it is followed by isoniazid (27.99 kcal mol⁻¹) and ethionamide (21.78 kcal mol⁻¹) (Fig. 4(c)).

Fig. 4 Mycolic acid membrane is responsible for low drug penetration properties of the TB bacterium. Resistance of the mycolic acid monolayer against different drug molecules for (a) CHARMM and (b) GROMOS FFs. Higher resistance recorded for ethambutol and ethionamide for the CHARMM FF, compared to the GROMOS FF. In the figures, vertical black dashed lines represent positions of water–monolayer interfaces.

(2) Two distinct peaks that appeared in the CHARMM PMF profile is not noticed for GROMOS parameters.

The other characteristics of the GROMOS PMF profiles are similar to those of CHARMM. For instance, the affinity of drug molecules for the membrane surfaces can be clearly seen from GROMOS PMF profiles. Besides this, the height of the free energy generally followed the spatial strength of the drug–membrane vdW interactions (Fig. 2(b)). However, although the strength of isoniazid–mycolic acid vdW interactions was slightly lower than that between ethionamide and the membrane, the free energy barrier for the former is higher (Fig. 3(b)).

The disappearance of the dual peaks in the GROMACS PMF profiles can be attributed to the significantly altered spatial profile of coulombic interactions and the number of H-bonds, compared to CHARMM parameters (Fig. 2(c), (d) and Fig. S4 in the ESI†). More importantly, from Fig. 2(d), it is clear that the electrostatic interaction between isoniazid and mycolic acid is stronger than that between ethionamide and the mycolic acid monolayer in the tail region of the molecules. This significantly affected the PMF profile and lowered the free energy barrier of ethionamide compared to isoniazid. This feature is absent for CHARMM parameters, and thus, the free energy barrier of ethionamide was higher than that of isoniazid.

The position-dependent resistance of the mycolic acid membrane against all three drug molecules, which follows the same pattern as the free energy barrier, has been represented using a log scale in Fig. 4, which followed a similar pattern like the PMF profile (Fig. 3), for both the FFs. The diffusion constant is presented in Fig. 5, and it shows that the diffusivity is extremely low for the mycolic acid monolayer, in agreement with the existing literature.⁶ The effective permeability of the membrane for different drug molecules is tabulated in Table 1 for both CHARMM and GROMOS FFs. For CHARMM force field parameters, isoniazid was found to be the most effective in penetrating the mycolic acid monolayer, and ethambutol is the least effective. However, ethionamide was found to be more effective in terms of permeability under the GROMOS FF (Table 1). It should be noted that the effective permeability was found to be higher for ethambutol and ethionamide for the GROMOS FF (Table 1). Besides this, the position dependent diffusion constant was always calculated to be higher for the GROMOS FF for all of the drug molecules (Fig. S5 in the ESI†), which is also a reason for the higher permeability of ethambutol with GROMOS parameters (Table 1).

Fig. 5 Mycolic acid monolayer acts as a diffusion barrier for drug molecules. Position-dependent diffusion coefficients of various pharmaceutical molecules in the mycolic acid monolayer for (a) CHARMM and (b) GROMOS FFs. In the figures, vertical black dashed lines represent positions of water–monolayer interfaces.

Table 1 Effective permeability (cm s⁻¹) of different drug molecules using two different FFs: the CHARMM36 all atom FF (CHARMM) and the GROMOS 54A7-ATB (GROMOS) united-atom FF

Name of the drugs	GROMOS	CHARMM
Ethambutol	9.61 × 10^−31	8.40 × 10^−34
Ethionamide	3.57 × 10^−19	1.45 × 10^−31
Isoniazid	1.68 × 10^−22	2.94 × 10^−17

---

3.3 Insight into the in silico behaviours

We will start the current sub-section with a brief description of the physicochemical properties of the TB drugs used for the study, followed by a summary of the characteristics of each FF employed for modelling. We will then summarize the previously described observed behaviors of different drug molecules with the underlying probable causes. Lastly, we shall comment on the advantages and disadvantages of each FF.

Along with isoniazid, ethambutol is a first-line drug used to treat TB. Ethambutol diffuses into the M.tb cell and inhibits arabinosyltransferases, hindering the process of cell wall formation and cell division.^49–51 Isoniazid, on the other side, is a prodrug that is activated by bacterial catalase. Once activated, it inhibits mycolic acid synthesis.⁵² Ethionamide is a second-line drug for the treatment of TB and is used in combination with other drugs. Like isoniazid, ethionamide is also a prodrug that requires activation.⁵³ Once activated, it disrupts the fatty acid synthesis necessary for the cell wall.⁵³ It is noteworthy that all of the drug molecules used in the present study are hydrophilic in nature and have less affinity toward the lipid phase. Hence, it is expected that the mycolic acid monolayer will act as a permeability barrier for all of these drugs, as depicted in the above-mentioned observations. Some essential characteristics, such as solubility and lipophilicity of the drug molecules, are tabulated in Table S1 in the ESI.†

As far as the FF parameters are concerned, both CHARMM and GROMOS FFs are well-known for the simulation of the lipid and lipid-like molecules. Historically, CHARMM27 was the first all-atom force field, which was widely used for different kinds of lipids.^54,55 Further modifications of the lipid parameter were incorporated in the later version, CHARMM36, and specific focus was given to the description of the lipid head group structure.⁵⁶ The CHARMM force field also includes a general version, the CHARMM general force field (CGENFF), for small biomolecules.^34,35 The functional form includes various internal degrees of freedom like the bond, proper and improper dihedral terms, including the Urey-Bradley term for the covalent bonds, determined mainly using quantum mechanical calculations.⁵⁷ CHARMM parameters are optimized to work within a 1.0–1.2 nm cut-off distance for the LJ interactions and they use Lorentz–Berthelot combination rules for LJ parameters.⁵⁴ Moreover, CHARMM parameters should be used in combination with a modified TIP3P water model (TIP3PM).⁵⁶ It is worth mentioning that the CHARMM FF family sometimes underestimates the area per lipid for constant pressure simulations.⁵⁷

In contrast, the GROMOS FF is a united atom force field, where non-polar CH, CH₂, and CH₃ groups are represented as a single bead.⁵⁷ Bonded parameters and charge distributions in the GROMOS FF are optimized based on electronic structure computations, while the development of non-bonded parameters is more focused on the reproduction of thermodynamic properties like enthalpy and free energies of solvation of model compounds.⁵⁷ Both CHARMM and GROMOS families of FFs exclude nearest and next-nearest non-bonded (1–2 and 1–3) interactions between bonded atoms from the sum over non-bonded terms. For the GROMOS FF, 1–4 interactions are also excluded for atoms within an aromatic ring.^58–60 The GROMOS family uses a geometric combination rule for LJ interactions and can be used with SPC or SPC/E water models.⁶⁰ GROMOS parameters, in many instances, overestimated the surface area/lipid for various kinds of lipids like DOPC, DLPC, DMPC and POPC, and such overestimation led to enhanced water penetration in the membrane.⁶¹

In the current study, we observed that the drug molecules differed in terms of their affinity toward the monolayer surface as a function of FF parameters. In the case of ethionamide and isoniazid, the GROMOS FF provided better stability of the drug molecules on the monolayer surface, whereas the affinity of ethambutol toward the monolayer surface was found to be similar for both the FFs (Fig. S11, ESI†). The underlying reasons can be explained in the following manner.

From Fig. S13 (ESI†), it is clear that the vdW interaction was almost similar for both sets of FFs in the case of ethambutol, but weaker for the CHARMM FF for ethionamide for approximately the first 125 ns. For isoniazid, vdW interaction was found to be weaker for the CHARMM FF almost throughout the simulation window (Fig. S13(c), ESI†). This weaker vdW interaction may contribute to the less stable dynamics, observed for ethionamide and isoniazid under the CHARMM FF. One interesting observation is that ethionamide electrostatically interacted with the monolayer (with the CHARMM FF), while the magnitude of coulombic interaction for the other two drugs (ethambutol and isoniazid) is not that significant (Fig. 5). This strong electrostatic interaction helped ethionamide to have more stability on the monolayer surface, compared to isoniazid under the CHARMM FF, as seen from the evolution of the minimum distance (Fig. S11(b) and (c), ESI†).

While comparing the dynamics of ethionamide using two different FFs, one interesting aspect that appears is that despite having greater coulombic interaction, ethionamide was found to have less affinity toward the monolayer under the CHARMM FF. To figure out the underlying reason, we have computed the radial distribution function (g(r)) of water oxygen atoms surrounding the drug molecules to assess the hydration effect on the dynamics of the drug molecules, and the calculated g(r) values for both sets of FFs are presented in Fig. S5 of the ESI.† From Fig. S5 (ESI†), it can be clearly seen that under the CHARMM FF, ethionamide and isoniazid had been surrounded by a greater number of water molecules compared to the GROMOS FF. Higher solvation contributed to the less stable behavior of ethionamide and isoniazid on the mycolic acid monolayer surface under the CHARMM FF, despite having greater drug–monolayer electrostatic interaction (for ethionamide). The solvation behavior was quite similar for ethambutol for both sets of FFs (Fig. S5(a) in ESI†), and it exhibited a similar affinity toward the monolayer surface for CHARMM and GROMOS.

It is a well-known fact that the GROMOS FF, in many instances, overestimated the surface area/lipid for various kinds of lipids like DOPC, DLPC, DMPC, and POPC, and such overestimation led to enhanced water penetration in the membrane.⁶¹ In the present study, the area per lipid head group was found to be higher for the GROMOS FF (Fig. S6 in ESI†). This gives rise to higher diffusivity for ethambutol, when the GROMOS FF is used. For the other two molecules, several other factors, like the drug–monolayer interactions, played an important role in determining the diffusivity, and this is reflected in the observed trend in permeability, as presented in Table 1.

Apart from the above-mentioned facts, we have computed the electrostatic potential and electric field across the membrane for both the FF in the absence of any drug molecule in order to obtain a better insight into the behavior of the drug molecules inside the membrane, and they are presented in Fig. S7 of the ESI.† The electrostatic potential was found to be higher at the water–monolayer interface for GROMOS (Fig. S7(a) in the ESI†), which, we believe, helped to better stabilize the drug molecules on the surface as well (Fig. S11, ESI†). Also, for the GROMOS FF, the electrostatic potential changed more smoothly inside the membrane, compared to that for the CHARMM FF. As per a previous report, a higher electrostatic potential compared to the other FF is a characteristic of the GROMOS FF.⁶² This difference also may have played a role in the notable difference in the PMF profile of the drug molecules discussed earlier. Another important fact is the drastic difference in the effective membrane permeability of different drug molecules, which arises from the exponential dependency of effective permeability on the free energy profile (see Section 2), because of which a slight change in the PMF profile results in a drastic change in effective permeability.⁶² In the current study, the possibility of passive transport of some first and second-line tuberculosis drugs has been explored using two different well-known FFs. Now, the general question arises about the suitability of the force field between GROMOS and CHARMM for modelling mycolic acid and drug interactions. We can pen a few points to answer this query. Firstly, the experimentally observed area per molecule for mycolic acid is ≥50 Ų with a limiting molecular area of 48 Ų.^63,64 Both FF sets estimate this quantity closer to the limiting molecular area (Fig. S6 in the ESI†), but the area per lipid is slightly close to the experimental values for the GROMOS FF. Secondly, the better stability of drugs on the monolayer surface while using the GROMOS FF probably represents a more realistic scenario as it is well known that the drugs used in this study can diffuse through the M.tb cell wall, and the first step of such phenomena is a steady attachment to the surface.

Lastly and most importantly comes the reliability of the permeability values we obtained from both sets of FFs. It is a well-known fact that the transport across membranes takes place in two major ways: active (proteins are involved) and passive (diffusion of small molecules across the membrane).⁶⁵ Relatively more hydrophobic/less soluble molecules can passively transport across the membrane comparatively easily than more soluble molecules.⁶⁵ From Table S1 in the ESI,† it is clear that among the drug molecules, ethambutol has the highest water solubility, followed by isoniazid and ethionamide. Therefore, the general expectation is that ethambutol will experience the highest free energy barrier, followed by isoniazid and ethionamide. This corroborates well with the GROMOS FF results. It should be iterated at this point that the permeability of drug molecules complexly depends on many physico-chemical characteristics other than solubility like the number of H-bond acceptors/donors, lipophilicity, molecular weight, size, shape, polar surface area, pK_a, minimum inhibitory concentration (MIC) etc.¹⁵ Less permeability of ethambutol also gets reflected in its higher MIC value (Table S1 in the ESI†).

Generally, the majority of the drug compounds face an energy barrier in the range of 5–50 kcal mol⁻¹ in the lipid assembly, and our results agree with this range for both FFs.^66–68 Permeability of isoniazid has been previously measured both computationally and experimentally against the POPC bilayer and artificial membrane.^69,70 Experimental results suggest that isoniazid is a low permeable drug (permeability <10 × 10⁻⁶ cm s⁻¹). Also, the permeability of isoniazid varies significantly with the change in membrane properties (from 0.3 ± 0.1 cm s⁻¹ to 5.4 ± 0.07 × 10⁻⁶ cm s⁻¹).^69,70 In the current study, we obtained even lower effective permeability values for all the drugs against the MA monolayer for both sets of FFs (Table 1), which is expected because mycolic acid is known to act as a very efficient permeability barrier.⁶ Unfortunately, to the best of our knowledge, there are no experimental observations about permeability of drugs against the MA monolayer to directly compare our results. Hence, considering all the points stated above, we can conclude that although GROMOS appears to be a slightly better FF for modelling drug–MA interactions, it is impossible to draw a final conclusion as to which FF is better to study the drug–MA interactions at the present state because of the deficit of experimental studies.

4. Conclusions

In the current study, we have successfully employed the atomistic MD simulation method to obtain a molecular level insight into the experimentally reported resistance of the MA monolayer toward conventional drug molecules. Our study delved into the phenomenological causes behind the observed low permeability of the MA monolayer and illustrated the importance of FF parameters to model such systems. The following key conclusions can be drawn from the present work.

(a) The MA monolayer offers a high free energy barrier of penetration toward drug molecules. The MA molecules significantly interacted with the drug molecules during their passage through the monolayer, which led to such high free energy barrier.

(b) For both the FFs, the drug molecules preferred to stay on the monolayer surface and interacted with the monolayer mainly via vdW interactions. However, better stability of drug molecules on the surface of the monolayer was observed for GROMOS FF, specifically for ethionamide and isoniazid. This can be attributed to the lower solvation of drug molecules together with a higher electrostatic potential at the water–membrane interface found for GROMOS FF. On the other hand, ethionamide exhibited affinity toward the tail of the mycolic acid monolayer for adsorption, which was not observed for simulations using the GROMOS FF.

(c) The PMF profile of different drug molecules looked significantly different for GROMOS and CHARMM FFs because of the notable difference in drug–mycolic acid coulombic interactions within the monolayer (during the passage of the drugs through the monolayer). The free energy barrier height also differed significantly under the two FFs because of the same reason. However, ethambutol was found to experience the highest free energy barrier for monolayer crossing for both the FFs.

(d) Apart from the drug–monolayer interaction, the height of the free energy barrier and diffusivity of the drug molecules were also found to depend on the area per lipid. A higher diffusivity was observed for the GROMOS FF, compared to that for the CHARMM FF, because of the higher area per head group of lipid obtained for the former FF set. Our simulation results indicate that the drugs will typically take time on the order of milliseconds to cross the monolayer through passive diffusion.

(e) GROMOS appeared to be a better FF for similar studies. But, lack of experimental evidence prohibits us from drawing a deterministic conclusion.

To summarize, the current work provides an atomistic insight into the origin of the thermodynamic barrier and low permeability of the MA monolayer for some of the well-known TB drugs. We also present a detailed comparison of two well-known FFs and their influences on the free energy barrier and permeability. From the current work, it has been clear that the interactions between MA molecules and drugs play an immense role in the permeability of the former and this knowledge provides a guideline for modeling drug–membrane systems with more details to design more potent TB drug candidates, as similar studies can be performed with potent drug-like molecules to assess their interactions with the M.tb cell wall and in turn their ability to diffuse within bacterial cells.

Data availability

The data validating the findings of our simulation study are available upon request to the corresponding author.

Conflicts of interest

The authors declare no conflicts of interest.

Acknowledgements

The authors thank the Department of Biotechnology (DBT), Government of India, for funding this work (Grant number: BT/PR33123/MED/29/1497/2020). The authors also thank Prof. Narayanaswamy Jayaraman, Prof. Dipankar Chatterji, and Dr Anirban Ghosh, IISc, Bangalore, for critical reading of the manuscript and useful suggestions.

References

1. J. Chakaya, M. Khan, F. Ntoumi, E. Aklillu, R. Fatima, P. Mwaba, N. Kapata, S. Mfinanga, S. E. Hasnain and P. D. M. C. Katoto, Global Tuberculosis Report 2020–Reflections on the Global TB burden, treatment and prevention efforts, Int. J. Infect. Dis., 2021, 113, S7–S12.
2. W. H. Organization, Global tuberculosis report 2013, World health organization, 2013.
3. W. H. Organization, The WHO global task force on TB impact measurement, World Health Organization, 2019.
4. S. V. Gordon and T. Parish, Microbe Profile: Mycobacterium tuberculosis: humanity's deadly microbial foe, Microbiology, 2018, 164, 437–439.
5. H. Nikaido, Prevention of drug access to bacterial targets: permeability barriers and active efflux, Science, 1994, 264, 382–388.
6. H. Marrakchi, M.-A. Lanéelle and M. Daffe, Mycolic acids: structures, biosynthesis, and beyond, Chem. Biol., 2014, 21, 67–85.
7. D. E. Minnikin, C. Ratledge and J. Stanford, The biology of mycobacteria, Biol. Mycobact., 1982, 1, 94–184.
8. P. A. Lambert, Cellular impermeability and uptake of biocides and antibiotics in Gram-positive bacteria and mycobacteria, J. Appl. Microbiol., 2002, 92, 46S–54S.
9. X. Hong and A. J. Hopfinger, Construction, molecular modeling, and simulation of Mycobacterium tuberculosis cell walls, Biomacromolecules, 2004, 5, 1052–1065.
10. X. Hong and A. J. Hopfinger, Molecular modeling and simulation of Mycobacterium tuberculosis cell wall permeability, Biomacromolecules, 2004, 5, 1066–1077.
11. W. Groenewald, M. S. Baird, J. A. Verschoor, D. E. Minnikin and A. K. Croft, Differential spontaneous folding of mycolic acids from Mycobacterium tuberculosis, Chem. Phys. Lipids, 2014, 180, 15–22.
12. W. Groenewald, R. A. Parra-Cruz, C. M. Jäger and A. K. Croft, Revealing solvent-dependent folding behavior of mycolic acids from Mycobacterium tuberculosis by advanced simulation analysis, J. Mol. Model., 2019, 25, 1–12.
13. B. Merget, D. Zilian, T. Müller and C. A. Sotriffer, MycPermCheck: the Mycobacterium tuberculosis permeability prediction tool for small molecules, Bioinformatics, 2013, 29, 62–68.
14. S. Nagamani and G. N. Sastry, Mycobacterium tuberculosis Cell Wall Permeability Model Generation Using Chemoinformatics and Machine Learning Approaches, ACS Omega, 2021, 6, 17472–17482.
15. S. Janardhan, M. Ram Vivek and G. Narahari Sastry, Modeling the permeability of drug-like molecules through the cell wall of Mycobacterium tuberculosis: an analogue based approach, Mol. BioSyst., 2016, 12, 3377–3384.
16. M. Lundborg, C. L. Wennberg, A. Narangifard, E. Lindahl and L. Norlén, Predicting drug permeability through skin using molecular dynamics simulation, J. Controlled Release, 2018, 283, 269–279.
17. W. Shinoda, Permeability across lipid membranes, Biochim. Biophys. Acta, Biomembr., 2016, 1858, 2254–2265.
18. R. Paul, S. Bera, M. Devi and S. Paul, Inhibition of Aβ16–22 Peptide Aggregation by Small Molecules and Their Permeation through POPC Lipid Bilayer: Insight from Molecular Dynamics Simulation Study, J. Chem. Inf. Model., 2022, 62, 5193–5207.
19. M. Sugita, T. Fujie, K. Yanagisawa, M. Ohue and Y. Akiyama, Lipid composition is critical for accurate membrane permeability prediction of cyclic peptides by molecular dynamics simulations, J. Chem. Inf. Model., 2022, 62, 4549–4560.
20. R. M. Venable, A. Kramer and R. W. Pastor, Molecular dynamics simulations of membrane permeability, Chem. Rev., 2019, 119, 5954–5997.
21. D. Bemporad, J. W. Essex and C. Luttmann, Permeation of small molecules through a lipid bilayer: a computer simulation study, J. Phys. Chem. B, 2004, 108, 4875–4884.
22. M. Orsi, W. E. Sanderson and J. W. Essex, Permeability of small molecules through a lipid bilayer: a multiscale simulation study, J. Phys. Chem. B, 2009, 113, 12019–12029.
23. P. Kedia, Y. Badhe, R. Gupta, S. Kausley and B. Rai, Modeling the effect of pH on the permeability of dried chitosan film, J. Food Eng., 2023, 358, 111682.
24. P. Dutta, R. Gupta, D. S. Jain and B. Rai, Deep learning models for the estimation of free energy of permeation of small molecules across lipid membranes, Digital Discovery, 2023, 2, 189–201.
25. Y. Mitsuta, T. Asada and Y. Shigeta, Calculation of the permeability coefficients of small molecules through lipid bilayers by free-energy reaction network analysis following the explicit treatment of the internal conformation of the solute, Phys. Chem. Chem. Phys., 2022, 24, 26070–26082.
26. C. J. Knight and J. S. Hub, MemGen: a general web server for the setup of lipid membrane simulation systems, Bioinformatics, 2015, 31, 2897–2899.
27. L. Martínez, R. Andrade, E. G. Birgin and J. M. Martínez, PACKMOL: a package for building initial configurations for molecular dynamics simulations, J. Comput. Chem., 2009, 30, 2157–2164.
28. J. Huang and A. D. MacKerell Jr, CHARMM36 all-atom additive protein force field: validation based on comparison to NMR data, J. Comput. Chem., 2013, 34, 2135–2145.
29. K. Vanommeslaeghe, E. Hatcher, C. Acharya, S. Kundu, S. Zhong, J. Shim, E. Darian, O. Guvench, P. Lopes and I. Vorobyov, CHARMM general force field: a force field for drug-like molecules compatible with the CHARMM all-atom additive biological force fields, J. Comput. Chem., 2010, 31, 671–690.
30. E. Neria, S. Fischer and M. Karplus, Simulation of activation free energies in molecular systems, J. Chem. Phys., 1996, 105, 1902–1921.
31. S. Jo, T. Kim, V. G. Iyer and W. Im, CHARMM-GUI: a web-based graphical user interface for CHARMM, J. Comput. Chem., 2008, 29, 1859–1865.
32. B. R. Brooks, C. L. Brooks III, A. D. Mackerell Jr, L. Nilsson, R. J. Petrella, B. Roux, Y. Won, G. Archontis, C. Bartels and S. Boresch, CHARMM: the biomolecular simulation program, J. Comput. Chem., 2009, 30, 1545–1614.
33. J. Lee, X. Cheng, J. M. Swails, M. S. Yeom, P. K. Eastman, J. A. Lemkul, S. Wei, J. Buckner, J. C. Jeong and Y. Qi, CHARMM-GUI input generator for NAMD, GROMACS, AMBER, OpenMM, and CHARMM/OpenMM simulations using the CHARMM36 additive force field, J. Chem. Theory Comput., 2016, 12, 405–413.
34. K. Vanommeslaeghe and A. D. MacKerell Jr, Automation of the CHARMM General Force Field (CGenFF) I: bond perception and atom typing, J. Chem. Inf. Model., 2012, 52, 3144–3154.
35. K. Vanommeslaeghe, E. P. Raman and A. D. MacKerell Jr, Automation of the CHARMM General Force Field (CGenFF) II: assignment of bonded parameters and partial atomic charges, J. Chem. Inf. Model., 2012, 52, 3155–3168.
36. N. Schmid, A. P. Eichenberger, A. Choutko, S. Riniker, M. Winger, A. E. Mark and W. F. Van Gunsteren, Definition and testing of the GROMOS force-field versions 54A7 and 54B7, Eur. Biophys. J., 2011, 40, 843–856.
37. A. K. Malde, L. Zuo, M. Breeze, M. Stroet, D. Poger, P. C. Nair, C. Oostenbrink and A. E. Mark, An automated force field topology builder (ATB) and repository: version 1.0, J. Chem. Theory Comput., 2011, 7, 4026–4037.
38. K. B. Koziara, M. Stroet, A. K. Malde and A. E. Mark, Testing and validation of the Automated Topology Builder (ATB) version 2.0: prediction of hydration free enthalpies, J. Comput.-Aided Mol. Des., 2014, 28, 221–233.
39. G. Bussi, D. Donadio and M. Parrinello, Canonical sampling through velocity rescaling, J. Chem. Phys., 2007, 126, 14101.
40. M. Parrinello and A. Rahman, Polymorphic transitions in single crystals: a new molecular dynamics method, J. Appl. Phys., 1981, 52, 7182–7190.
41. B. Hess, H. Bekker, H. J. C. Berendsen and J. G. E. M. Fraaije, LINCS: a linear constraint solver for molecular simulations, J. Comput. Chem., 1997, 18, 1463–1472.
42. T. Darden, D. York and L. Pedersen, Particle mesh Ewald: An N·log(N) method for Ewald sums in large systems, J. Chem. Phys., 1993, 98, 10089–10092.
43. M. J. Abraham, T. Murtola, R. Schulz, S. Páll, J. C. Smith, B. Hess and E. Lindahl, GROMACS: high performance molecular simulations through multi-level parallelism from laptops to supercomputers, SoftwareX, 2015, 1, 19–25.
44. W. Humphrey, A. Dalke and K. Schulten, VMD: visual molecular dynamics, J. Mol. Graphics, 1996, 14, 33–38.
45. J. S. Hub, B. L. De Groot and D. van der Spoel, g_wham A Free Weighted Histogram Analysis Implementation Including Robust Error and Autocorrelation Estimates, J. Chem. Theory Comput., 2010, 6, 3713–3720.
46. G. Hummer, Position-dependent diffusion coefficients and free energies from Bayesian analysis of equilibrium and replica molecular dynamics simulations, New J. Phys., 2005, 7, 34.
47. C. T. Lee, J. Comer, C. Herndon, N. Leung, A. Pavlova, R. V. Swift, C. Tung, C. N. Rowley, R. E. Amaro and C. Chipot, Simulation-based approaches for determining membrane permeability of small compounds, J. Chem. Inf. Model., 2016, 56, 721–733.
48. J. P. M. Jambeck and A. P. Lyubartsev, Exploring the free energy landscape of solutes embedded in lipid bilayers, J. Phys. Chem. Lett., 2013, 4, 1781–1787.
49. J. P. Thomas, C. O. Baughn, R. G. Wilkinson and R. G. Shepherd, A new synthetic compound with antituberculous activity in mice: ethambutol (dextro-2,2′-(ethylenediimino)-di-1-butanol), Am. Rev. Respir. Dis., 1961, 83, 891–893.
50. R. Goude, A. G. Amin, D. Chatterjee and T. Parish, The arabinosyltransferase EmbC is inhibited by ethambutol in Mycobacterium tuberculosis, Antimicrob. Agents Chemother., 2009, 53, 4138–4146.
51. L. Zhang, Y. Zhao, Y. Gao, L. Wu, R. Gao, Q. Zhang, Y. Wang, C. Wu, F. Wu and S. S. Gurcha, Structures of cell wall arabinosyltransferases with the anti-tuberculosis drug ethambutol, Science, 2020, 368, 1211–1219.
52. J. Suarez, K. Ranguelova, A. A. Jarzecki, J. Manzerova, V. Krymov, X. Zhao, S. Yu, L. Metlitsky, G. J. Gerfen and R. S. Magliozzo, An oxyferrous heme/protein-based radical intermediate is catalytically competent in the catalase reaction of Mycobacterium tuberculosis catalase-peroxidase (KatG), J. Biol. Chem., 2009, 284, 7017–7029.
53. J. H. Hoofnagle, Drug-induced liver disease, Elsevier, 2013, pp. 725–732.
54. A. D. MacKerell Jr, J. Wiorkiewicz-Kuczera and M. Karplus, An all-atom empirical energy function for the simulation of nucleic acids, J. Am. Chem. Soc., 1995, 117, 11946–11975.
55. J. B. Klauda, B. R. Brooks, A. D. MacKerell, R. M. Venable and R. W. Pastor, An ab initio study on the torsional surface of alkanes and its effect on molecular simulations of alkanes and a DPPC bilayer, J. Phys. Chem. B, 2005, 109, 5300–5311.
56. J. B. Klauda, R. M. Venable, J. A. Freites, J. W. O’Connor, D. J. Tobias, C. Mondragon-Ramirez, I. Vorobyov, A. D. MacKerell Jr and R. W. Pastor, Update of the CHARMM all-atom additive force field for lipids: validation on six lipid types, J. Phys. Chem. B, 2010, 114, 7830–7843.
57. A. P. Lyubartsev and A. L. Rabinovich, Force field development for lipid membrane simulations, Biochim. Biophys. Acta, Biomembr., 2016, 1858, 2483–2497.
58. J. Hermans, H. J. C. Berendsen, W. F. Van Gunsteren and J. P. M. Postma, A consistent empirical potential for water–protein interactions, Biopolymers, 1984, 23, 1513–1518.
59. L. D. Schuler, X. Daura and W. F. Van Gunsteren, An improved GROMOS96 force field for aliphatic hydrocarbons in the condensed phase, J. Comput. Chem., 2001, 22, 1205–1218.
60. I. Chandrasekhar, M. Kastenholz, R. D. Lins, C. Oostenbrink, L. D. Schuler, D. P. Tieleman and W. F. van Gunsteren, A consistent potential energy parameter set for lipids: dipalmitoylphosphatidylcholine as a benchmark of the GROMOS96 45A3 force field, Eur. Biophys. J., 2003, 32, 67–77.
61. A. N. Leonard, E. Wang, V. Monje-Galvan and J. B. Klauda, Developing and testing of lipid force fields with applications to modeling cellular membranes, Chem. Rev., 2019, 119, 6227–6269.
62. S. Kanchi, M. Gosika, K. G. Ayappa and P. K. Maiti, Dendrimer Interactions with Lipid Bilayer: Comparison of Force Field and Effect of Implicit vs. Explicit Solvation, J. Chem. Theory Comput., 2018, 14, 3825–3839.
63. Z. Zhang, Y. Pen, R. G. Edyvean, S. A. Banwart, R. M. Dalgliesh and M. Geoghegan, Adhesive and conformational behaviour of mycolic acid monolayers, Biochim. Biophys. Acta, Biomembr., 2010, 1798, 1829–1839.
64. M. Villeneuve, M. Kawai, H. Kanashima, M. Watanabe, D. E. Minnikin and H. Nakahara, Temperature dependence of the Langmuir monolayer packing of mycolic acids from Mycobacterium tuberculosis, Biochim. Biophys. Acta, Biomembr., 2005, 1715, 71–80.
65. O. Sten-Knudsen, Passive transport processes, Springer, 1978.
66. P. Sharma, R. Vaiwala, A. K. Gopinath, R. Chockalingam and K. G. Ayappa, Structure of the Bacterial Cell Envelope and Interactions with Antimicrobials: Insights from Molecular Dynamics Simulations, Langmuir, 2024, 40(15), 7791–7811.
67. B. J. Bennion, N. A. Be, M. W. McNerney, V. Lao, E. M. Carlson, C. A. Valdez, M. A. Malfatti, H. A. Enright, T. H. Nguyen and F. C. Lightstone, Predicting a drug's membrane permeability: a computational model validated with in vitro permeability assay data, J. Phys. Chem. B, 2017, 121, 5228–5237.
68. M. Orsi and J. W. Essex, Permeability of drugs and hormones through a lipid bilayer: insights from dual-resolution molecular dynamics, Soft Matter, 2010, 6, 3797–3808.
69. M. S. Poka, M. Milne, A. Wessels and M. Aucamp, The Effect of Isoniazid–Maltitol Solid Dispersions on Aqueous Solubility and Permeability, Crystals, 2023, 13, 1568.
70. C. F. de Faria, T. Moreira, P. Lopes, H. Costa, J. R. Krewall, C. M. Barton, S. Santos, D. Goodwin, D. Machado, M. Viveiros, M. Machuqueiro and F. Martins, Designing new antitubercular isoniazid derivatives with improved reactivity and membrane trafficking abilities, Biomed. Pharmacother., 2021, 144, 112362.

---

Footnote

† Electronic supplementary information (ESI) available: Results of unbiased MD simulation together with the corresponding figures, additional figures related to snapshots of the systems' final conformations obtained from the unbiased MD, the number of hydrogen bonds, the radial distribution function of water oxygens, the area per lipid, electrostatic potential and electric field intensity across the MA monolayer, PMF convergence, umbrella sampling (US) histograms, the orientation of the drug molecules, and tables related to the physico-chemical properties of TB drug molecules. See DOI: https://doi.org/10.1039/d4cp02659d

---