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

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 MEEECL20-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. This journal is © The Royal Society of Chemistry.

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.[1] 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.[16] 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,[19] while an instantaneous release can be obtained if both blocks becomes soluble and the self-assembled structures fall apart.[20]

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.[21] This polymer is soluble at room temperature but becomes insoluble above the so-called coil-to-globule transition temperature Tcg (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.[9] Differently from the widely used poly(N-isopropylacrylamide) (PNIPAM),[13] PMEEECL is fully biocompatible.[22] Cheng et al.[22] copolymerized PMEEECL with poly(γ-octyloxy-ε-caprolactone) (PMEEECL-b-POCTCL) by ring-opening polymerization. Although the obtained block copolymers showed a low cytotoxicity and a Tcg 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.[23] 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] MEEECL20 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,[28] have been used to determine the equilibrium morphology of MEEECL20–AE self-assembled structures as a function of the hydrophobic block length M. Based on the SCF results, three PMEEECL20–PAE (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 Tcg 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 PMEEECL20–PAE5 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.

Methods

MD simulation

Both all-atom and MARTINI CG MD simulations were performed by LAMMPS[35] 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 PMEEECL20–PAE 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[36] was used.

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

All-atom simulations protocol

For the all-atom simulations, the OPLS-AA,[37] GROMOS 54A7[38] and TIP3P[39] force fields have been used for PMEEECL–PAE, DOX and water molecules, respectively. Firstly, a dispersion with a concentration of 4 wt% of PMEEECL20–PAE 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[40] (accuracy of 10−5) 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 PMEEECL20–PAE 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 Tcg, the radius of gyration of the blocks (Rg), 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 Rg of PMEEECL was fitted by the Boltzmann sigmoid function:[41]where ΔT is the transition width[42] which describes the steepness of the curve and it has been chosen 5 K for PMEEECL.

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 Pmicelle (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.[45] Between 400 ns and 600 ns the value of Pmicelle 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 Rmicelle and core size Rcore were calculated using the following equations:[46]where Np is the number of micelle beads, rcom is the position of centre-of-mass of micelle, r and rPAE, 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 mencDOX and mtotDOX are the weight of encapsulated DOX and weight of total DOX, respectively. Moreover, mmicellepolymer is the weight of polymer in the micelle.

The diffusion coefficient of DOX beads (DDOX) was calculated from the mean square displacement (MSD) as follows:[47,48]where r(0) and r(t) are the position of bead i at initial time and time t, respectively.

SCF theory computations

Numerical self-consistent field theory computations were performed using the method developed by Scheutjens and Fleer.[26] The computations were executed using the program SFbox.[49] Spherical, cylindrical and a planar lattice geometries are assumed to study the formation of spherical, cylindrical micelles and vesicles, respectively.[50] 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 PMEEECL20–PAE copolymer is modelled as A20–(BC), 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 pKSCFa. In the SCF procedure these constants depend on the lattice constant, and can be determined from the experimental pKa values:

In the computations we used pKSCFaW = 15.35 and pKSCFaC = 9.5, the latter estimated from literature values.[52] Since two potentially charged groups are present on each PAE monomer, the protonation of C is modelled according to the equilibrium:

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 χ (Table S2†). For the interactions between A, B, C, and W the χ values where estimated from the Hansen solubility parameters[53] obtained via the group contribution method.[54] All other interaction parameters were set to 0.

Results and discussion

Preferred morphology of PMEEECL20–PAE

SCF computations have been performed first to determine the preferred morphology of PMEEECL20–PAE 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 ΔGmic is related with ϕbulkvia the relation ΔGmic = kBT ln ϕbulk.[56] The results indicate that at neutral pH and at room temperature PMEEECL20–PAE 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.

Effect of hydrophobic block on Tcg

The typical signature of the Tcg is a large drop in the polymer size, expressed via the value of Rg, which is due to a collapse of the chains. By fitting the Rg values obtained from the simulation results with a Boltzmann sigmoid function (see eqn (1))[41] (see Fig. 2a and Table S3†) Tcg values = 312, 302 and 298 K with 2.5 K error bars estimated from the figure have been found for PMEEECL block in PMEEECL20–PAE5, PMEEECL20–PAE10, and PMEEECL20–PAE15, respectively. The value of Tcg decreases with increasing hydrophobic PAE block length. Tcg as a function of (Fig. 2b and Table S4†) of the PMEEECL mass fraction f can be reasonably described using Tcg = 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 (Rg) for the hydrophilic (PMEEECL) part of PMEEECL20–PAE 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 (Tcg) of PMEEECL in the block copolymer as a function of PMEEECL mass fraction. For f = 1 the Tcg has been derived from ref. 24. The dot-dashed curve follows Tcg = 26.89 ln(f) + 318.96.

To gain better understanding on the coil-to-globule transition behaviour of PMEEECL20–PAE 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 PMEEECL20–PAE 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 Rg, 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 Tcg.

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 Tcg (Fig. S5†), testifying that the hydrophilic–hydrophobic transition is caused by a temperature-induced dehydration, in agreement with previous literature on LCST polymers.[57] 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.

Self-assembly of PMEEECL20–PAE

Once the relation between Tcg and the PAE block length was established, the self-assembly behaviour of the three selected copolymers, PMEEECL20–PAE5, PMEEECL20–PAE10, and PMEEECL20–PAE15, 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 PMEEECL20–PAE 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 PMEEECL20–PAE5 forms small spherical micelles, PMEEECL20–PAE10 forms cylindrical micelles and PMEEECL20–PAE15 forms bilayer structures.

Fig. 3

Snapshot series from micellization of PMEEECL20–PAE 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 PMEEECL20–PAE 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.

Response to stimuli

Response to pH by SCF computations

The self-assembly of PMEEECL20–PAE5, PMEEECL20–PAE10 and PMEEECL20–PAE15 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 PMEEECL20–PAE block copolymer with M = 5, 10, and 15 in water.

DOX encapsulation

As mentioned in the introduction section, 310–315 K would be an appropriate range of the coil-to-globule transition temperature Tcg for the amphiphilic block copolymer for drug delivery applications.[13] By using all-atom MD simulations, Tcg of PMEEECL20–PAE with M = 5, 10, and 15 was calculated and the appropriate amphiphilic block copolymer for drug delivery application was found to be PMEEECL20–PAE5 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, PMEEECL20–PAE5 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 PMEEECL20–PAE5) 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.[58] 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 PMEEECL20–PAE5 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 PMEEECL20–PAE5) resulting from MARTINI CG force field simulations.

The self-diffusion coefficient of DOX (DDOX) 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, DDOX lies in range ∼4800–6500 μm2 s−1 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 μm2 s−1. 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 DDOX at three different concentrations in a solution of spherical micelle composed of PMEEECL20–PAE5 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 PMEEECL20–PAE5 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.[58] 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.[58]

Fig. 9

Simulation box snapshots from the last simulation step including PMEEECL20–PAE5 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.

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 MEEECL20–AE5 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.