Dual responsive PMEEECL–PAE block copolymers: a computational self-assembly and doxorubicin uptake study

Abstract

The self-assembly behaviour of dual-responsive block copolymers and their ability to solubilize the anticancer drug doxorubicin (DOX) has been investigated using all-atom molecular dynamics (MD) simulations, MARTINI coarse-grained (CG) force field simulation and Scheutjens–Fleer self-consistent field (SCF) computations. These diblock copolymers, composed of poly{γ-2-[2-(2-methoxyethoxy)ethoxy]ethoxy-ε-caprolactone} (PMEEECL) and poly(β-amino ester) (PAE) are dual-responsive: the PMEEECL block is thermoresponsive (becomes insoluble above a certain temperature), while the PAE block is pH-responsive (becomes soluble below a certain pH). Three MEEECL₂₀–AE_M compositions with M = 5, 10, and 15, have been studied. All-atom MD simulations have been performed to calculate the coil-to-globule transition temperature (T_cg) of these copolymers and finding appropriate CG mapping for both PMEEECL–PAE and DOX. The output of the MARTINI CG simulations is in agreement with SCF predictions. The results show that DOX is solubilized with high efficiency (75–80%) at different concentrations inside the PMEEECL–PAE micelles, although, interestingly, the loading efficiency is reduced by increasing the drug concentration. The non-bonded interaction energy and the RDF between DOX and water beads confirm this result. Finally, MD simulations and SCF computations reveal that the responsive behaviour of PMEEECL–PAE self-assembled structures take place at temperature and pH ranges appropriate for drug delivery.

---

1. Introduction

The self-assembly of amphiphilic block copolymers in water results in different types of supramolecular structures, such as spherical micelles, cylindrical micelles and vesicles, depending on the chemical composition and relative mass of the blocks.¹ These self-assembled structures are composed of an insoluble core made of the hydrophobic blocks, surrounded by a hydrophilic corona which provide stability in solution.^2,3 Bio-compatible and bio-degradable amphiphilic block copolymers have attracted the attention of researchers for their potential use in bio-applications such as drug delivery.^4,5 In fact, the core of the assemblies can be used to encapsulate drug while the corona can provide protection against immune recognition (opsonization).^6–8 Recently, it was demonstrated that block copolymers which can respond to external stimuli such as temperature and pH are efficacious for drug delivery applications,^9–12 as tumor tissues are warmer^9,13 and more acidic^14,15 than normal ones.

Thermoresponsive polymers show a change in solubility upon changing temperature. Differently, pH-responsive polymers undergo a hydrophilic to hydrophobic transition upon changing pH due to charge effect.¹⁶ These two types of responsive polymers can be combined into block copolymers where the solubility of each block can be controlled independently. This gives rise to a dual-responsive self-assembly behaviour which can be exploited to release drugs in a controlled fashion.^17,18 For instance, a slow drug release can be achieved when both blocks are in a hydrophobic state resulting in a compact self-assembled structure,¹⁹ while an instantaneous release can be obtained if both blocks becomes soluble and the self-assembled structures fall apart.²⁰

In this in silico study we designed a new class of dual-responsive block copolymers and we computationally investigated their self-assembly behaviour and their ability to solubilize the anticancer drug doxorubicin (DOX). As corona-forming block the thermoresponsive poly{γ-2-[2-(2-methoxyethoxy)ethoxy]ethoxy-3-caprolactone} (PMEEECL) has been selected.²¹ This polymer is soluble at room temperature but becomes insoluble above the so-called coil-to-globule transition temperature T_cg (i.e. the temperature that marks the transition from swollen to collapsed state of the polymer). This kind of phase behaviour is named lower critical solution temperature (LCST) and is usually observed in solutions containing hydrophilic polymers which interact with water via hydrogen bonds.⁹ Differently from the widely used poly(N-isopropylacrylamide) (PNIPAM),¹³ PMEEECL is fully biocompatible.²² Cheng et al.²² copolymerized PMEEECL with poly(γ-octyloxy-ε-caprolactone) (PMEEECL-b-POCTCL) by ring-opening polymerization. Although the obtained block copolymers showed a low cytotoxicity and a T_cg around 38 °C, the loading efficiency of doxorubicin (DOX) was not sufficiently high due to the characteristics of the core-forming blocks. Min et al.²³ prepared micelles from block copolymers using the pH-responsive biodegradable poly(β-amino ester) (PAE) as core-forming block for the delivery of camptothecin (CPT) and observed a significant increase in therapeutic efficacy with minimum side effects using these micelles containing CPT. Moreover, the appropriate loading efficiency around 96% and 82% was observed with this block copolymer when CPT was fed as the ratio of 5 wt% and 10 wt% to polymeric micelles, respectively. Therefore, PAE has been selected as core-forming block in this study. On the basis of these results one could expect that block copolymers made of PMEEECL and PAE (PMEEECL–PAE) exhibit the desired dual-responsive behaviour combined with a high drug loading ability.

A step-by-step approach has been followed to study the self-assembly behaviour of PMEEECL–PAE, their dual-responsivity and their uptake of DOX. Based on previous investigations^24,25 MEEECL₂₀ has been chosen as hydrophilic block. Scheutjens–Fleer self-consistent field (SCF) computations,^26,27 which are suitable to perform systematic investigations due to low computational cost,²⁸ have been used to determine the equilibrium morphology of MEEECL₂₀–AE_M self-assembled structures as a function of the hydrophobic block length M. Based on the SCF results, three PMEEECL₂₀–PAE_M (M = 5, 10, and 15) compositions, forming, respectively, spherical micelles, cylindrical micelles and vesicles (or flat bilayers) have been selected. All atom molecular dynamics (MD) simulations have been performed on these copolymers to compute the value of T_cg as a function of the PAE block length. Subsequently, MARTINI coarse-grained (CG) force field simulations^29,30 have been used to study the self-assembly behaviour of the copolymers, particularly the structure and the equilibrium morphology of the resulting assemblies. These computational methods have been previously used to study block copolymer self-assembly providing good agreement with experimental results.^31–34 It is best suited to study macromolecular aggregates because the coarse-graining allows describing the system with less degrees of freedom, making the simulation much faster than all-atom MD. The encapsulation of DOX within PMEEECL₂₀–PAE₅ spherical micelles at different DOX concentrations has been investigated by MARTINI CG force field simulations. Finally, the effect of temperature and pH on the self-assembled structure was studied by MARTINI CG force field and SCF theory, respectively. The results obtained with the different methods were in good agreement and showed the great potential of PMEEECL–PAE block copolymer for drug delivery applications. To our knowledge, such a combination of computational techniques has not yet been used to study dual-response block copolymers.

2. Methods

2.1. MD simulation

Both all-atom and MARTINI CG MD simulations were performed by LAMMPS³⁵ under isothermal–isobaric NPT (constant pressure P, temperature T and number of particles in the box N) conditions. Details of the simulations are provided in Table S1.† Moreover, by using the Nose–Hoover algorithm the temperature and pressure were controlled with a coupling time of 0.1 ps and a coupling time of 1 ps, respectively. To calculate the solvent access surface area (SASA), the radial distribution function (RDF) between the PMEEECL₂₀–PAE_M amphiphilic block copolymer and hydrogen atoms of water molecules as well as the snapshots of simulation systems for both all-atom MD and MARTINI CG force field simulations, Visual Molecular Dynamic (VMD) package³⁶ was used.

Below, a detailed description of the applied simulation protocols and the performed data analysis is reported.

2.1.1. All-atom simulations protocol. For the all-atom simulations, the OPLS-AA,³⁷ GROMOS 54A7³⁸ and TIP3P³⁹ force fields have been used for PMEEECL–PAE, DOX and water molecules, respectively. Firstly, a dispersion with a concentration of 4 wt% of PMEEECL₂₀–PAE_M chains with M = 5, 10, or 15 was equilibrated for 10 ns and subsequent simulations were performed for 50 ns. The last 35 ns of the simulation runs were used for statistical analysis. The simulations were performed using a 2 fs time step. In order to compare the bond and angle distribution for validating the CG mapping (see below) the simulation run was performed for 10 ns. The particle/particle particle/mesh (PPPM) method⁴⁰ (accuracy of 10⁻⁵) was used to evaluate the long-range electrostatic interactions. A 12 Å cut-off distance was used to estimate the non-bonded dispersive interactions. A single PMEEECL₂₀–PAE_M chain was always placed in the simulation box at an appropriate distance from the periodic images.

To study the thermoresponsive behaviour of the copolymers and identify the T_cg, the radius of gyration of the blocks (R_g), the SASA, the non-bonded interaction energy between block copolymer and solvent molecules and the RDF between the oxygen of PMEEECL block and hydrogen of water molecules were calculated from the simulation results.

The temperature dependence of R_g of PMEEECL was fitted by the Boltzmann sigmoid function:⁴¹

where ΔT is the transition width⁴² which describes the steepness of the curve and it has been chosen 5 K for PMEEECL.

2.1.2. Coarse-grained simulations protocol. The MARTINI CG force field^29,30 and a 20 fs time step have been used in the CG simulations to study the self-assembly of the block copolymers and the DOX encapsulation processes. The initial simulation configuration was constructed by randomly filling the simulation box with the block copolymer chains, water and DOX molecules. Moreover, the shift function starting at 9 Å with a cut-off of 12 Å was used for the dispersion interactions. The CG mapping of PMEEECL–PAE block copolymers and DOX, shown in Fig. 1, have been chosen (after testing several mapping for CG methods) such that the specific bond and angle distributions obtained with all-atom MD and MARTINI CG force field simulations were in agreement (Fig. S1 and S2 in ESI†). Although there is some difference between the angular distributions between the all-atom and CG models for small and large angles, see Fig. S2,† the most probable angles of the distributions are close for the all-atom and CG model simulations. Inspection of the bending coefficient resulting from using the MARTINI and GROMOS force fields reveals that it is lower for the MARTINI model. This explains a wider range of probable angles of the distribution. It is noted that such different angle distributions of all-atom and CG models were discussed previously.^43,44

Fig. 1 MARTINI coarse-grained (CG) mapping for (a) a PMEEECL–PAE block copolymer, (b) DOX used for MARTINI CG force field simulations and (c) SCF coarse-graining of PMEEECL–PAE block copolymers.

The evaluation of the order parameter P_micelle (Fig. S3†), was monitored to set an appropriate simulation timed. It has been defined as where η_micelle is the micelle dimensionless density (volume fraction) and has been calculated by the volume average of the difference between local and overall density squared.⁴⁵ Between 400 ns and 600 ns the value of P_micelle changed hardly varied, so a simulation time of 600 ns was set. The self-assembled structural properties of block copolymer such as total radius of the micelle R_micelle and core size R_core were calculated using the following equations:⁴⁶

where N_p is the number of micelle beads, r_com is the position of centre-of-mass of micelle, r_i and r_PAE,i are the position of outmost CG beads of micelle and core, respectively.

To quantify the loading efficiency (LE) and capacity (LC) of DOX into the PMEEECL–PAE self-assembled structures the following equations have been used:

where m^enc_DOX and m^tot_DOX are the weight of encapsulated DOX and weight of total DOX, respectively. Moreover, m^micelle_polymer is the weight of polymer in the micelle.

The diffusion coefficient of DOX beads (D_DOX) was calculated from the mean square displacement (MSD) as follows:^47,48

where r_i(0) and r_i(t) are the position of bead i at initial time and time t, respectively.

2.1.3. SCF theory computations. Numerical self-consistent field theory computations were performed using the method developed by Scheutjens and Fleer.²⁶ The computations were executed using the program SFbox.⁴⁹ Spherical, cylindrical and a planar lattice geometries are assumed to study the formation of spherical, cylindrical micelles and vesicles, respectively.⁵⁰ Within these lattice approaches concentration gradients were accounted for in one direction only. In the spherical and cylindrical lattice, the gradient is set along the radial direction. In the flat lattice the gradient is perpendicular to the lattice plane. Each lattice was composed of L = 200 lattice layers in the gradient direction; for larger L values the results are unaffected by L for the systems studied. As coordination number of the lattice z = 3 was used, in agreement with previous investigations.^50,51 Mirroring conditions are assumed at the outermost (for all) and innermost (for flat geometry only) lattice layers. Each lattice site was considered to have a size λ = 1.5 nm, estimated from the size of a PMEEECL monomer calculated from equilibrated state all-atom MD simulation results. Each PMEEECL monomer (A) and water molecule (W) are considered to occupy a single lattice site, while each monomer of the PAE block is modelled as composed of two units, as shown in Fig. 1c (B and C), each occupying one lattice site. Hence a PMEEECL₂₀–PAE_M copolymer is modelled as A₂₀–(BC)_M, where the C unit is the pH-responsive unit.

All the computations were performed in the presence of a fixed amount of monovalent salt in the lattice to adjust the ionic strength. To this end the volume fraction of the cations (I⁺) in the lattice is set to ϕ⁺ = 0.03, while the anions (I⁻) are used as charge neutralizers, meaning that their volume fraction (ϕ⁻) is automatically adjusted to achieve the overall charge neutrality of the lattice. Both cations and anions are set to occupy single lattice sites.

The acid–base reactions of water and the pH responsive C group are modelled as chemical equilibria by specifying the reaction acid dissociation constants pK^SCF_a. In the SCF procedure these constants depend on the lattice constant, and can be determined from the experimental pK_a values:

pK^SCF_a = pK_a − log₁₀[10⁻²⁴λ³N_av] (8)

In the computations we used pK^SCF_aW = 15.35 and pK^SCF_aC = 9.5, the latter estimated from literature values.⁵² Since two potentially charged groups are present on each PAE monomer, the protonation of C is modelled according to the equilibrium:

C + 2W⁺ ⇌ C²⁺ + 2W (9)

It should be noted that this simple description was interpreted to align with the coarse-grained model description. The short-range interactions are modelled via a set of Flory–Huggins interaction parameters χ_ij (Table S2†). For the interactions between A, B, C, and W the χ_ij values where estimated from the Hansen solubility parameters⁵³ obtained via the group contribution method.⁵⁴ All other interaction parameters were set to 0.

3. Results and discussion

3.1. Preferred morphology of PMEEECL₂₀–PAE_M

SCF computations have been performed first to determine the preferred morphology of PMEEECL₂₀–PAE_M self-assembled structures at pH = 7 and room temperature as a function of the PAE block length 2 < M < 20. For each copolymer the equilibrium value of the critical aggregation concentration, ϕ^bulk, have been computed in each of the three lattice geometries. For details about this procedure we refer to previous literature.^26,28,55 The preferred geometry is the one associated with the lowest ϕ^bulk value, since the free energy of micelle formation ΔG^mic is related with ϕ^bulk via the relation ΔG^mic = k_BT ln ϕ^bulk.⁵⁶ The results indicate that at neutral pH and at room temperature PMEEECL₂₀–PAE_M block copolymers in water assemble into spherical micelles for M ≤ 6, into cylindrical micelles for 7 ≤ M ≤ 12 and into vesicles or flat bilayers for M ≥ 13. Based on these results three copolymer compositions (M = 5, 10, 15) assembling in three different geometries, have been selected for in-depth studies.

3.2. Effect of hydrophobic block on T_cg

The typical signature of the T_cg is a large drop in the polymer size, expressed via the value of R_g, which is due to a collapse of the chains. By fitting the R_g values obtained from the simulation results with a Boltzmann sigmoid function (see eqn (1))⁴¹ (see Fig. 2a and Table S3†) T_cg values = 312, 302 and 298 K with 2.5 K error bars estimated from the figure have been found for PMEEECL block in PMEEECL₂₀–PAE₅, PMEEECL₂₀–PAE₁₀, and PMEEECL₂₀–PAE₁₅, respectively. The value of T_cg decreases with increasing hydrophobic PAE block length. T_cg as a function of (Fig. 2b and Table S4†) of the PMEEECL mass fraction f can be reasonably described using T_cg = 26.89 ln(f) + 318.96. Hence, the thermoresponsive behaviour of the copolymer can be controlled by adjusting the PAE block length. It should be noted that no observable size change was detected for the PAE blocks because it is not a thermosensitive polymer (see Table S3†).

Fig. 2 (a) All-atom MD simulations results for the radius of gyration (R_g) for the hydrophilic (PMEEECL) part of PMEEECL₂₀–PAE_M block copolymer with M = 5, 10, and 15 as a function of temperature. The dotted curve is a fit using eqn (1). (b) Coil-to-globule transition temperature (T_cg) of PMEEECL in the block copolymer as a function of PMEEECL mass fraction. For f = 1 the T_cg has been derived from ref. 24. The dot-dashed curve follows T_cg = 26.89 ln(f) + 318.96.

To gain better understanding on the coil-to-globule transition behaviour of PMEEECL₂₀–PAE_M block copolymers, the non-bonded interaction energies (including van der Waals and electrostatic energies) between the copolymer blocks and the solvent as well as the solvent accessible surface area (SASA) have been calculated (Fig. S4†). These parameters describe the affinity between blocks and solvent molecules which enable to quantify the transition behaviour of the blocks. For instance, SASA is related to the surface area of PMEEECL₂₀–PAE_M block copolymers which is accessible to the water molecules. It can be calculated by extending the van der Waals radius for each atom by 1.4 Å. A reduction of the SASA value implies a reduction of the access surface area for water molecules along the PMEEECL chain. Therefore, in agreement with the temperature-dependence of R_g, the non-bonded interaction energies and SASA values of PMEEECL decreased significantly, indicating a reduction in the affinity between the water molecules and PMEEECL. This corresponds to the observed trends in T_cg.

Coil-to-globule transition behaviour can also be observed in the radial distribution function (RDF) of the hydrogen atom of water with respect to the PMEEECL oxygen atoms (Fig. S5†), which have been computed at two different temperatures, 285 K (hydrophilic state) and 335 K (hydrophobic state). These RDFs describe the probability to find water hydrogens at a certain distance from the PMEEECL oxygen atoms, hence they are related with the PMEEECL-water hydrogen bonding probability. A substantial in the RDF in proximity to the oxygen atoms can be observed above T_cg (Fig. S5†), testifying that the hydrophilic–hydrophobic transition is caused by a temperature-induced dehydration, in agreement with previous literature on LCST polymers.⁵⁷ Unlike the PMEEECL block, the RDF of centre of mass (COM) of PAE block and hydrogen of water molecules (Fig. S6†) shows no reduction upon increasing the temperature which is in agreement with previous results.

3.3. Self-assembly of PMEEECL₂₀–PAE_M

Once the relation between T_cg and the PAE block length was established, the self-assembly behaviour of the three selected copolymers, PMEEECL₂₀–PAE₅, PMEEECL₂₀–PAE₁₀, and PMEEECL₂₀–PAE₁₅, at 300 K has been investigated by MARTINI CG force field simulations (Fig. 3). In Fig. 3 the self-assembly process is depicted with a series of snapshot after different simulation times where the water beads are removed for clarity. Initially, 200 PMEEECL₂₀–PAE_M chains were randomly placed in the simulation box with 250 000 water beads with an approximate concentration of 8 wt%, 10 wt%, and 12 wt% for M = 5, 10, and 15, respectively. After 200 ns of simulation time the onset of self-assembly could be observed and by further increasing the simulation time the self-assembled structure where formed. These remained stable and did not coalesce for the last 100 ns of simulation. In addition, by increasing the PAE block length, the self-assembly process speeds up, as shown by the evolution of the order parameter (Fig. S3†). Details of the spherical and cylindrical and bilayer morphologies are presented in the fifth column of Fig. 3 for more clarification. The quantitative description of the self-assembly behaviour, i.e. the number of cluster, as a function of simulation time for spherical morphology is plotted in Fig. S7† which reveals two spherical micelles remain after a simulation time of 600 ns. It follows in agreement with the SCF predictions, that PMEEECL₂₀–PAE₅ forms small spherical micelles, PMEEECL₂₀–PAE₁₀ forms cylindrical micelles and PMEEECL₂₀–PAE₁₅ forms bilayer structures.

Fig. 3 Snapshot series from micellization of PMEEECL₂₀–PAE_M with M = 5, 10, and 15. Water beads are removed for clarity and blue and red beads show PMEEECL and PAE beads, respectively.

The self-assembly of the three copolymers have also been studied with SCF theory which, similarly to the CG MARTINI method, enables to generate radial concentration profiles (Fig. 4). These profiles describe the spatial distribution of the components as a function of the distance from the centre of the self-assembled structures. For both methods at room temperature and neutral pH the concentration profiles predict the self-assembled structures to be composed of a dehydrated core (only 5–10% water present) rich in PAE, surrounded by a hydrated corona rich in PMEEECL.

Fig. 4 Radial concentration profiles of self-assembled PMEEECL₂₀–PAE_M block copolymer with M = 5, 10, and 15 resulting from MARTINI CG force field simulations (a, c and e) and SCF computations (b, d and f). Red, blue and gray curves are PMEEECL, PAE and water, respectively.

These regions are separated by an interface which appears more diffuse in the CG MARTINI simulations (Fig. 4). This difference is probably caused by three main factors. Firstly, the interaction potentials used in the CG MARTINI method are softer. Secondly, the interaction parameters used in the two methods are different. Thirdly, in the SCF computations the concentration profiles are directly calculated optimizing the free-energy and assuming perfect geometries while in the CG MARTINI simulations they are calculated from the beads distribution in space. The CG MARTINI self-assembled structures appear quite irregular in shape and exhibit areas where the lyophobic core is directly in contact with the solvent. When performing a radial average such irregularities result in an apparently larger interface between the lyophilic and lyophobic domains, which explains the jumps observed in the CG MARTINI water concentration profile at the interface between the two blocks. The SCF computations predict bigger sizes for the assemblies (as can be observed in the concentration profiles in Fig. 4) but this difference is simply due to the different coarse-graining procedures applied, which does not fully reflect the size, the steric hindrance and the flexibility of the copolymer chains. Overall, it may be concluded that the SCF predictions and computer simulation results are in qualitative and semi-quantitative agreement.

3.4. Response to stimuli

3.4.1. Response to pH by SCF computations. The self-assembly of PMEEECL₂₀–PAE₅, PMEEECL₂₀–PAE₁₀ and PMEEECL₂₀–PAE₁₅ have been studied as a function of pH in the range 8.5 < pH < 5. At each pH the ϕ^bulk value was determined together with the equilibrium aggregation number g. In the spherical geometry the aggregation number quantifies the number of copolymers per micelle. In the cylindrical geometry it indicates the number of copolymer per unit length of the cylindrical micelle, while in the flat geometry it corresponds to the number of copolymers per unit area. For all three block copolymers dissociation into unimers was observed at pH < 5.4 due to protonation of the PAE block. The transition between the assembled and unassembled states appears to be characterized by a sudden increase of ϕ^bulk (Fig. 5a) and a drop of g (Fig. 5b–d) values. This is probably an artifact originating from the simple way we model the protonation equilibrium of the C unit (to be either neutral or doubly charged, see eqn (8)). One could expect to observe a (more) gradual transition experimentally, where the protonation of the two PAE basic groups might occur sequentially.

Fig. 5 SCF results (data) for the pH dependency of (a) the critical aggregation concentration ϕ^bulk and (b–d) the aggregation number g of the self-assembled structures of PMEEECL₂₀–PAE_M block copolymer with M = 5, 10, and 15 in water.

3.5. DOX encapsulation

As mentioned in the introduction section, 310–315 K would be an appropriate range of the coil-to-globule transition temperature T_cg for the amphiphilic block copolymer for drug delivery applications.¹³ By using all-atom MD simulations, T_cg of PMEEECL₂₀–PAE_M with M = 5, 10, and 15 was calculated and the appropriate amphiphilic block copolymer for drug delivery application was found to be PMEEECL₂₀–PAE₅ with 312 ± 2.5 K for. These block copolymers self-assemble into a spherical morphology based on both SCF computations and CG simulations. The spherical geometry is an appropriate morphology for particles moving through a body environment. Transport of for instance cylindrical micelles for instance is easily hampered due to their significant lengths, and also vesicles are typically big. Therefore, PMEEECL₂₀–PAE₅ block copolymers with a radius R close to 7 nm have been selected to study the solubilisation of doxorubicin at 300 K at three DOX concentrations (1 wt%, 5 wt%, and 10 wt%; here the wt% is the amount with respect to the total amount of PMEEECL₂₀–PAE₅) with CG MARTINI force field simulations. We note that, since only two spherical micelles appear in the simulation box (see Fig. 3), the aggregation number and size computed from the simulations are approximate. The CG MARTINI concentration profiles (Fig. 6) enable to visualize the spatial distribution of DOX within the micelles and how this is affected by the DOX concentration. At 1 wt%, DOX is mainly located at the interface between the micelle's core and corona (Fig. 6a). With increasing DOX concentration the location of DOX beads shifts from the interface of core and corona towards the corona section (Fig. 6b). Still, for 10 wt% (Fig. 6c) there are some DOX beads at the external layer of micelle, although the core encapsulation dominates in comparison to 1 wt% and 5 wt% which can be occurred at guest, i.e. DOX, concentrations higher than the solubility.⁵⁸ It should be also noted that by increasing the DOX concentration the size of micelle is increased.

Fig. 6 Radial concentration profile of a spherical PMEEECL₂₀–PAE₅ micelle (a) without presence of DOX beads, and mixed DOX beads for three DOX concentrations: (b) 1 wt%, (c) 5 wt%, and (d) 10 wt% (wt% of the total amount of PMEEECL₂₀–PAE₅) resulting from MARTINI CG force field simulations.

The self-diffusion coefficient of DOX (D_DOX) was computed from the mean-square displacements (MSD) determined using the CG simulations. The MSDs were converted to the D values by using eqn (6) and (7) and plotted in Fig. 7 for three DOX concentrations. Prior to self-assembly, D_DOX lies in range ∼4800–6500 μm² s⁻¹ but after 600 ns of simulation, when the micelles were formed and DOX is solubilized, the drug diffusivity decreases about one order of magnitude to ∼500 μm² s⁻¹. This can be explained by the fact that the DOX beads first diffuse as individual beads, whereas they diffuse together with the block copolymer micelle after encapsulation. A series of snapshots from the first and last part of the simulations are provided in Fig. 8 (for full movies refer to ESI Videos S1–S3†). It can be observed that the DOX beads move from the solution to the micelles as the micelles form.

Fig. 7 Simulation time dependence of the self-diffusion coefficient of DOX beads D_DOX at three different concentrations in a solution of spherical micelle composed of PMEEECL₂₀–PAE₅ block copolymer and water beads.

Fig. 8 Snapshots from the start (upper graphs) and last parts (lower graphs) of MARTINI CG force field simulations including 100 chains of PMEEECL₂₀–PAE₅ and DOX 1 wt% (a and d), 5 wt% (b and e), and 10 wt% (c and f) in a spherical micelle. Blue, red, and green beads show PMEEECL, PAE and DOX beads, respectively.

The DOX loading efficiency, defined as the ratio between loaded in the core and at the interface of core and corona of the micelle with respect to the total amount of DOX (see Methods) was approximately 80%, 78%, and 75% when DOX was fed as the ratio of 1 wt%, 5 wt%, and 10 wt% to block copolymers, respectively, after 600 ns of simulation time. This small decrease in the solubilization might be related to the reduction of micelle capacity by increasing the concentration of DOX. As can be observed, see Fig. 8, all DOX beads were captured by micelles at the three different concentrations studied. This implies that the non-loaded DOX molecules reside in the corona of the micelles. The presence of DOX resulted in an increase in the final number of micelles. Hence the average aggregation number, i.e. the number of block copolymer chains added to the simulation box divided by the number of micelles, reduced from 26.25 to 15 by increasing the DOX concentration from 1 wt% to 10 wt%, respectively. This is caused by the interfacial adsorption of DOX, which reduces the interfacial tension of the core. This effect was predicted theoretically and verified experimentally previously.⁵⁸ Its effect can be seen for example in Fig. 9: in absence of DOX only two spherical micelles were found, but by adding DOX at three different concentrations more stable micelles with the same simulation time were observed, which corresponds to theoretical predictions and experimental observations.⁵⁸

Fig. 9 Simulation box snapshots from the last simulation step including PMEEECL₂₀–PAE₅ block copolymers (a) without DOX and with (b) 1 wt%, (c) 5 wt%, and (d) 10 wt%. CG beads are selected as periodic image in x-direction to show better micelles. Water beads are removed for clarity and blue, red, and green beads show PMEEECL, PAE, and DOX beads, respectively.

4. Conclusions

In this work, a dual-responsive PMEEECL–PAE block copolymer was designed in silico and the self-assembly behaviour of these polymers was studied with a combination of all-atom molecular dynamics (MD) simulations, MARTINI coarse-grained (CG) force field simulation and Scheutjens–Fleer self-consistent field (SCF) computations. The results reveal that the PMEEECL blocks are thermoresponsive, and become insoluble above a certain temperature, of which the value can be tuned by varying the copolymer composition via the PMEEECL mass fraction f. Upon decreasing f, spherical micelles, cylindrical micelles and bilayer self-assembled structures can be obtained. These self-assembled structures are also pH sensitive, and disassemble at pH < 5.4 due to PAE protonation. Responsivity at this pH range is ideal for drug delivery applications. Finally, it is demonstrated that spherical MEEECL₂₀–AE₅ micelles can solubilize the anticancer drug doxorubicin (DOX) with high efficiency (75–80%) at different concentrations. The solubilisation occurs via interfacial adsorption of DOX at the micelle core–corona interface at low DOX concentration. Upon increasing the DOX concentration the preferred DOX accumulation location shifts from the interface towards the corona area. The results provided by this unique combination of computational techniques show that PMEEECL–PAE block copolymers are great candidates for stimuli controlled drug delivery applications.

Conflicts of interest

There are no conflicts to declare.

References

1. Y. Mai and A. Eisenberg, Self-assembly of block copolymers, Chem. Soc. Rev., 2012, 41, 5969–5985,  10.1039/C2CS35115C.
2. Y. Tang, S. Y. Liu, S. P. Armes and N. C. Billingham, Solubilization and controlled release of a hydrophobic drug using novel micelle-forming ABC triblock copolymers, Biomacromolecules, 2003, 4, 1636–1645,  DOI:10.1021/bm030026t.
3. R. D. Wesley, C. A. Dreiss, T. Cosgrove, S. P. Armes, W. Bank, L. Thompson, F. L. Baines and N. C. Billingham, Structure of a hydrophilic–hydrophobic block copolymer and its interactions with salt and an anionic surfactant, Langmuir, 2005, 21, 4856–4861,  DOI:10.1021/la046830y.
4. H. Oliveira, E. Pérez-Andrés, J. Thevenot, O. Sandre, E. Berra and S. Lecommandoux, Magnetic field triggered drug release from polymersomes for cancer therapeutics, J. Controlled Release, 2013, 169, 165–170,  DOI:10.1016/j.jconrel.2013.01.013.
5. G. Liu, S. Ma, S. Li, R. Cheng, F. Meng, H. Liu and Z. Zhong, The highly efficient delivery of exogenous proteins into cells mediated by biodegradable chimaeric polymersomes, Biomaterials, 2010, 31, 7575–7585,  DOI:10.1016/j.biomaterials.2010.06.021.
6. Y. H. Bae and K. Park, Targeted drug delivery to tumors: myths, reality and possibility, J. Controlled Release, 2011, 153, 198–205,  DOI:10.1016/j.jconrel.2011.06.001.
7. G. T. Chang, T. Y. Ci, L. Yu and J. D. Ding, Enhancement of the fraction of the active form of an antitumor drug topotecan via an injectable hydrogel, J. Controlled Release, 2011, 156, 21–27,  DOI:10.1016/j.jconrel.2011.07.008.
8. Z. P. Zhang and S. S. Feng, The drug encapsulation efficiency, in vitro drug release, cellular uptake and cytotoxicity of paclitaxel-loaded poly(lactide)-tocopheryl polyethylene glycol succinate nanoparticles, Biomaterials, 2006, 27, 4025–4033,  DOI:10.1016/j.biomaterials.2006.03.006.
9. A. Gandhi, A. Paul, S. O. Sen and K. K. Sen, Studies on thermoresponsive polymers: phase behaviour, drug delivery and biomedical applications, Asian J. Pharm. Sci., 2015, 10, 99–107,  DOI:10.1016/j.ajps.2014.08.010.
10. E. S. Lee, K. T. Oh, D. Kim, Y. S. Youn and Y. H. Bae, Tumor pH-responsive flower-like micelles of poly(L-lactic acid)-b-poly (ethylene glycol)-b-poly(L-histidine), J. Controlled Release, 2007, 123, 19–26,  DOI:10.1016/j.jconrel.2007.08.006.
11. V. S. Miguel, A. J. Limer, D. M. Haddleton, F. Catalina and C. Peinado, Biodegradable and thermoresponsive micelles of triblock copolymers based on 2-(N,N-dimethylamino)ethyl methacrylate and e-caprolactone for controlled drug delivery, Eur. Polym. J., 2008, 44, 3853–3863,  DOI:10.1016/j.eurpolymj.2008.07.056.
12. S. Q. Liu, N. Wiradharma, S. J. Gao, Y. W. Tong and Y. Y. Yang, Bio-functional micelles self-assembled from a folate-conjugated block copolymer for targeted intracellular delivery of anticancer drugs, Biomaterials, 2007, 28, 1423–1433,  DOI:10.1016/j.biomaterials.2006.11.013.
13. N. Kamaly, B. Yameen, J. Wu and O. C. Farokhzad, Degradable controlled-release polymers and polymeric nanoparticles: mechanisms of controlling drug release, Chem. Rev., 2016, 116, 2602–2663,  DOI:10.1021/acs.chemrev.5b00346.
14. D. J. Adams, M. L. Wahl, J. L. Flowers, B. Sen, M. Colvin, M. W. Dewhirst, G. Manikumar and M. C. Wani, Camptothecin analogs with enhanced activity against human breast cancer cells. II. Impact of the tumor pH gradient, Cancer Chemother. Pharmacol., 2006, 57, 145–154,  DOI:10.1007/s00280-005-0008-5.
15. D. J. Adams and L. R. Morgan, Tumor physiology and charge dynamics of anticancer drugs: implications for camptothecin-based drug development, Curr. Med. Chem., 2011, 18, 1367–1372,  DOI:10.2174/092986711795029609.
16. Z. Luo, Y. Li, B. Wang and J. Jiang, pH-sensitive vesicles formed by amphiphilic grafted copolymers with tunable membrane permeability for drug loading/release: a multiscale simulation study, Macromolecules, 2016, 49, 6084–6094,  DOI:10.1021/acs.macromol.6b01211.
17. H. Wei, X. Z. Zhang, H. Cheng, W. Q. Chen, S. X. Cheng and R. X. Zhuo, Self-assembled thermo- and pH-responsive micelles of poly(10-undecenoic acid-b-N-isopropylacrylamide) for drug delivery, J. Controlled Release, 2006, 116, 266–274,  DOI:10.1016/j.jconrel.2006.08.018.
18. R. S. Kalhapure and J. Renukuntla, Thermo- and pH dual responsive polymeric micelles and nanoparticles, Chem.-Biol. Interact., 2018, 295, 20–37,  DOI:10.1016/j.cbi.2018.07.016.
19. Á. González García, A. Ianiro and R. Tuinier, On the colloidal stability of spherical copolymeric micelles, ACS Omega, 2018, 3, 17976–17985,  DOI:10.1021/acsomega.8b02548.
20. J. Sankaranarayanan, E. A. Mahmoud, G. Kim, J. M. Morachis and A. Almutairi, Multiresponse strategies to modulate burst degradation and release from nanoparticles, ACS Nano, 2010, 4, 5930–5936,  DOI:10.1021/nn100968e.
21. J. Hao, J. Servello, P. Sista, M. C. Biewer and M. C. Stefan, Temperature-sensitive aliphatic polyesters: synthesis and characterization of γ-substituted caprolactone monomers and polymers, J. Mater. Chem., 2011, 21, 10623–10628,  10.1039/C1JM11288K.
22. Y. Cheng, J. Hao, L. A. Lee, M. C. Biewer, Q. Wang and M. C. Stefan, Thermally controlled release of anticancer drug from self-assembled γ-substituted amphiphilic poly(ε-caprolactone) micellar nanoparticles, Biomacromolecules, 2012, 13, 2163–2173,  DOI:10.1021/bm300823y.
23. K. H. Min, J. H. Kim, S. M. Bae, H. Shin, M. S. Kim, S. Park, H. Lee, R. W. Park, I. S. Kim, K. Kim, I. C. Kwon, S. Y. Jeong and D. S. Lee, Tumoral acidic pH-responsive MPEG-poly(β-amino ester) polymeric micelles for cancer targeting therapy, J. Controlled Release, 2010, 144, 259–266,  DOI:10.1016/j.jconrel.2010.02.024.
24. A. Koochaki, M. R. Moghbeli and S. Javan Nikkhah, Coil-to-globule transition of thermo-responsive γ-substituted poly (ε-caprolactone) in water: a molecular dynamics simulation study, Curr. Appl. Phys., 2018, 18, 1313–1319,  DOI:10.1016/j.cap.2018.07.011.
25. A. Koochaki, M. R. Moghbeli and S. Javan Nikkhah, Effect of γ-substituted poly(ε-caprolactone) chain length on its coil-to-globule transition temperature in water: a molecular dynamics simulation study, Chem. Phys., 2019, 527, 110506,  DOI:10.1016/j.chemphys.2019.110506.
26. G. J. Fleer, M. A. Cohen Stuart, J. M. H. M. Scheutjens, T. Cosgrove and B. Vincent, Polymers at Interfaces, Springer, Netherlands, 1998.
27. P. N. Hurter, J. M. H. M. Scheutjens and T. A. Hatton, Molecular modeling of micelle formation and solubilization in block copolymer micelles. 1. A self-consistent mean-field lattice theory, Macromolecules, 1993, 26, 5592–5601,  DOI:10.1021/ma00073a010.
28. A. Ianiro, J. Patterson, Á. González García, M. M. J. van Rijt, M. M. R. M. Hendrix, N. A. J. M. Sommerdijk, I. K. Voets, A. C. C. Esteves and R. Tuinier, A roadmap for poly(ethylene oxide)-block-poly-ε-caprolactone self-assembly in water: prediction, synthesis, and characterization, J. Polym. Sci., Part B: Polym. Phys., 2018, 56, 330–339,  DOI:10.1002/polb.24545.
29. S. J. Marrink, H. J. Risselada, S. Yefimov, D. P. Tieleman and A. H. De Vries, The MARTINI force field: coarse grained model for biomolecular simulations, J. Phys. Chem. B, 2007, 111, 7812–7824,  DOI:10.1021/jp071097f.
30. S. J. Marrink, A. H. de Vries and A. E. Mark, Coarse grained model for semiquantitative lipid simulations, J. Phys. Chem. B, 2004, 108, 750–760,  DOI:10.1021/jp036508g.
31. S. M. Loverde, M. L. Klein and D. E. Discher, Nanoparticle shape improves delivery: rational coarse grain molecular dynamics (rCG-MD) of taxol in worm-like PEG-PCL micelles, Adv. Mater., 2012, 24, 3823–3830,  DOI:10.1002/adma.201103192.
32. D. J. Beltran-Villegas, M. G. Wessels, J. Y. Lee, Y. Song, K. L. Wooley, D. J. Pochan and A. Jayaraman, Computational Reverse-Engineering Analysis for Scattering Experiments on Amphiphilic Block Polymer Solutions, J. Am. Chem. Soc., 2019, 141, 14916–14930,  DOI:10.1021/jacs.9b08028.
33. Z. Li and E. E. Dormidontova, Kinetics of Diblock Copolymer Micellization by Dissipative Particle Dynamics, Macromolecules, 2010, 43, 3521–3531,  DOI:10.1021/ma902860j.
34. A. Prhashanna and S. Bor Chen, Chain exchange kinetics between linear ABA-type triblock copolymer micelles, Polymer, 2017, 118, 22–29,  DOI:10.1016/j.polymer.2017.04.049.
35. S. Plimpton, Fast Parallel Algorithms for Short-Range Molecular Dynamics, J. Comput. Phys., 1995, 117, 1–19,  DOI:10.1006/jcph.1995.1039.
36. W. Humphrey, A. Dalke and K. Schulten, VMD: visual molecular dynamics, J. Mol. Graphics, 1996, 14, 33–38,  DOI:10.1016/0263-7855(96)00018-5.
37. W. L. Jorgensen and J. Tirado-Rives, The OPLS [optimized potentials for liquid simulations] potential functions for proteins, energy minimizations for crystals of cyclic peptides and crambin, J. Am. Chem. Soc., 1988, 110, 1657–1666,  DOI:10.1021/ja00214a001.
38. 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,  DOI:10.1007/s00249-011-0700-9.
39. W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, Comparison of simple potential functions for simulating liquid water, J. Chem. Phys., 1983, 79, 926,  DOI:10.1063/1.445869.
40. R. W. Hockney and J. W. Eastwood, Computer Simulation Using Particles, Taylor & Francis, Inc., 1988.
41. D. J. Finney, Bioassay and the practice of statistical inference, International Statistical Review/Revue Internationale de Statistique, 1979, 47, 1–12,  DOI:10.2307/1403201.
42. T. E. de Oliveira, C. M. Marques and P. A. Netz, Molecular dynamics study of the LCST transition in aqueous poly(N-n-propylacrylamide), Phys. Chem. Chem. Phys., 2018, 20, 10100–10107,  10.1039/C8CP00481A.
43. S. Nawaz and P. Carbone, Coarse-Graining Poly(ethylene oxide)–Poly(propylene oxide)–Poly(ethylene oxide) (PEO–PPO–PEO) Block Copolymers Using the MARTINI Force Field, J. Phys. Chem. B, 2014, 118, 1648–1659,  DOI:10.1021/jp4092249.
44. L. Qiu, J. Liu, R. Alessandri, X. Qiu, M. Koopmans, R. W. A. Havenith, S. J. Marrink, R. C. Chiechi, L. J. A. Koster and J. C. Hummelen, Enhancing doping efficiency by improving host-dopant miscibility for fullerene-based n-type thermoelectrics, J. Mater. Chem. A, 2017, 5, 21234,  10.1039/c7ta06609k.
45. D. Mu, X.-R. Huang, Z.-Y. Lu and C.-C. Sun, Computer simulation study on the compatibility of poly(ethylene oxide)/poly(methyl methacrylate) blends, Chem. Phys., 2008, 348, 122–129,  DOI:10.1016/j.chemphys.2008.03.015.
46. R. Mackenzie, J. Booth, C. Alexander, M. C. Garnett and C. A. Laughton, Multiscale modeling of drug-polymer nanoparticle assembly identifies parameters influencing drug encapsulation efficiency, J. Chem. Theory Comput., 2015, 11, 2705–2713,  DOI:10.1021/ct501152a.
47. Y. Chen, Q. L. Liu, A. M. Zhu, Q. G. Zhang and J. Y. Wu, Molecular simulation of CO₂/CH₄ permeabilities in polyamide–imide isomers, J. Membr. Sci., 2010, 348, 204–212,  DOI:10.1016/j.memsci.2009.11.002.
48. S. Rasouli, M. R. Moghbeli and S. J. Nikkhah, A comprehensive molecular dynamics study of a single polystyrene chain in a good solvent, Curr. Appl. Phys., 2018, 18, 68–78,  DOI:10.1016/j.cap.2017.10.010.
49. The SCF code was provided by F. A. M. Leermakers and developed by F. A. M. Leermakers and J. van Male from Wageningen University, The Netherlands.
50. F. A. M. Leermakers, J. C. Eriksson and J. Lyklema, Association colloids and their equilibrium modelling, in Fundamentals of Interface and Colloid Science, ed. J. Lyklema, 2005, vol. 5, pp. 4.1–4.123,  DOI:10.1016/S1874-5679(05)80008-X.
51. Y. Lauw, F. A. M. Leermakers and M. A. Cohen Stuart, Self-consistent-field analysis of the micellization of carboxy-modified poly(ethylene oxide)–poly(propylene oxide)–poly(ethylene oxide) triblock copolymers, J. Phys. Chem. B, 2006, 110, 465–477,  DOI:10.1021/jp053795a.
52. A. V. Rayer, K. Z. Sumon, L. Jaffari and A. Henni, Dissociation constants (pK_a) of tertiary and cyclic amines: structural and temperature dependences, J. Chem. Eng. Data, 2014, 59, 3805–3813,  DOI:10.1021/je500680q.
53. E. Stefanis and C. Panayiotou, Prediction of hansen solubility parameters with a new group-contribution method, Int. J. Thermophys., 2008, 29, 568–585,  DOI:10.1007/s10765-008-0415-z.
54. T. Lindvig, M. L. Michelsen and G. M. Kontogeorgis, A Flory–Huggins model based on the Hansen solubility parameters, Fluid Phase Equilib., 2002, 203, 247–260,  DOI:10.1016/S0378-3812(02)00184-X.
55. J. G. J. L. Lebouille, F. A. M. Leermakers, M. A. Cohen Stuart and R. Tuinier, Design of block-copolymer-based micelles for active and passive targeting, Phys. Rev. E, 2016, 94, 042503,  DOI:10.1103/PhysRevE.94.042503.
56. H. J. Butt, K. Graf and M. Kappl, Physics and Chemistry of Interfaces, Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim, FRG, 2003.
57. T. E. De Oliveira, D. Mukherji, K. Kremer and P. A. Netz, Effects of stereochemistry and copolymerization on the LCST of PNIPAm, J. Chem. Phys., 2017, 146, 034904,  DOI:10.1063/1.4974165.
58. A. Ianiro, Á. González García, S. Wijker, J. P. Patterson, A. C. C. Esteves and R. Tuinier, Controlling the spatial distribution of solubilized compounds within copolymer micelles, Langmuir, 2019, 35, 4776–4786,  DOI:10.1021/acs.langmuir.9b00180.

---

Footnote

† Electronic supplementary information (ESI) available: Computer simulations (Table S1), Flory–Huggins interaction parameters (Table S2), radii of gyration of PMEEECL and PAE parts of block copolymer (Table S3), the hydrophilic mass fraction dependency of coil-to-globule transition temperature of PMEEECL in block copolymers (Table S4), the bond and angle distribution of PMEEECL, PAE (Fig. S1), and DOX (Fig. S2), the simulation time-dependency of micelle order parameter (Fig. S3), solvent access surface area as a function of interaction energy for both PMEEECL and PAE blocks (Fig. S4), the RDF of PMEEECL block and water molecules (Fig. S5), the RDF of PAE block and water molecules (Fig. S6), the number of cluster for spherical micelles as a function of simulation time (Fig. S7), and the RDF of PMEEECL beads and water beads as a function of temperature (Fig. S8). Additionally, see Videos S1–S3 which show the encapsulation process of DOX into spherical micelle composed of PMEEECL₂₀–PAE₅ when it was fed as the ratio of 1 wt%, 5 wt%, and 10 wt% to block copolymers. See DOI: 10.1039/c9ra09066e

---