Molecular Origins of the Ultra-Low Friction Exhibited by Biocompatible Zwitterionic Polymer Brushes
Surface grafted polymer brushes comprised of poly(2-(methacryloyloxy)ethyl phosphorylcholine) (pMPC) have attracted the attention of researchers due to their unique lubricating properties. These biocompatible materials exhibit friction coefficients orders of magnitude lower than those used in current artificial joints yet also readily withstand wear. This property is believed to arise from a mechanism called hydration lubrication where the polymer brushes tightly bind surrounding water molecules yet still offer a fluid-like response to shear. However, the molecular origins of this behavior, specifically where and how the water molecules bind, remain unknown. Here, we use molecular dynamics simulation to study a pMPC brush in aqueous solution. Our results show that the phosphate group is the predominant contributor to the hydration lubrication mechanism but that the choline group may also serve a stabilizing role, which potentially contributes to the brush’s remarkable durability.
The present study provides insight into the highly efficient lubrication mechanism in the human synovial joint and synthetic polymer systems that mimic it.7
Experimentally, zwitterionic polymer brushes have been found to exhibit friction coefficients as low as 0.0004.6 These low friction coefficients are even lower than those exhibited by the mammalian synovial joint which can yield friction coefficients of 0.001 to 0.1; moreover, surfaces functionalized with pMPC brushes have been found to yield the lowest friction coefficients recorded in surfaces functionalized with zwitterionic polymer brushes while under pressure.5, 8, 7 Additionally, pMPC is particularly attractive for bio-lubrication applications due to its high biocompatibility.9
Therefore, pMPC can be considered as an optimal candidate for implementation in artificial synovial joints. Its biocompatibility yields the potential for grafting to prosthetic materials without health complications, where it would facilitate a more lubricated joint with friction coefficients greatly reduced compared to those found in a normal prosthetic implant.
Experimental research has supported the hypothesis that the mechanism of lubrication of pMPC brushes occurs via hydration lubrication but is limited to studying the macroscopic brush behavior. While it is understood that an aqueous environment facilitates low friction coefficients the behavior of individual polymers and the water around them is not known. There is no definitive understanding of what functional groups in the individual monomers of the polymer support hydration lubrication. Additionally, the role of intra- and inter-monomer interactions is not well understood, specifically, how they influence the hydration lubrication exhibited by these systems.5, 8, 1
Molecular dynamics simulation is a way to create in silico models of extremely small structures such as pMPC polymers. Simulation allows researchers to investigate the polymers at the molecular level. Specifically, simulating an individual pMPC polymer, using parameters that are known to accurately mimic the experimental structures of pMPC, allows molecular insight into the behavior of the polymer and the solvent with a specific analysis on the groups actively contributing to hydration lubrication. This insight will ultimately further the understanding of hydration lubrication with respect to zwitterionic polymers.
MATERIALS AND METHODS.
The pMPC 50-mer was solvated in a rhombic dodecahedron, which reduces the number of required solvent molecules, of volume 40 nm3 leaving at least 2 nm of space between the polymer and its periodic images. Simulations were performed in GROMACS11 using the OPLS-aa forcefield3,4 and the newly refined TIP4P-FB water model.10 The TIP4P-FB water molecule is most applicable to this study in that it better reproduces diffusivity than previous models. A short NVT equilibration was run after energy minimization at 300 K using the Nose-Hoover thermostat. Lennard-Jones interactions had a cutoff radius of 0.8 nm with a switching function between 0.8-0.9 nm. The system was run with a time step of 2 fs. Bond constraints were placed on all hydrogen bonds. Electrostatic attractions were modeled using the fast smooth particle-mesh Ewald algorithm with grid spacing of 0.1 nm and a cutoff of 1 nm.
One 83 ns simulation was performed using the aforementioned parameters and preliminary results were obtained. Analysis was performed with snapshots taken every 100ps. GROMACS and MDTraj10 tools were used to conduct the analysis addressed in the present study. In order to record the configuration of pMPC50 throughout the simulation, the GROMACS analysis tools were used to calculate the radius of gyrationnd radial distribution function (RDF), and MDTraj analysis tools were used to calculate the root mean square deviation (RMSD). To study the active roles of pMPC50 charged groups, GROMACS analysis tools were used to calculate the RDF and the solvent accessible surface area (SASA).11,10
Configuration of pMPC50.
Within Figure 1A the intersection between the plots of the y and z axis supports that the polymer has assumed a folded state at around 60 ns. Radius of gyration calculations show that pMPC50 folds on itself after 60 ns of simulation. This is supported qualitatively in Figure 2C where the polymer at 0 ns is shown with no vertical or horizontal symmetry. After 60 ns in Figure 1A the x, y, and z axis intersect and after about 78 ns the lines have nearly converged. This is indicative of a more spherical configuration, or a more condensed polymer. Additionally, this is supported in Figure 2B where the two ends of the polymer are nearly touching. A root mean square deviation (RMSD) calculation of the polymer for the entire simulation (data not shown) supplements the radius of gyration data by showing that the polymer at the end of the simulation deviated from its starting configuration by 3 nm.
Figure 1A) The radius of gyration of solvated pMPC50 for the entire 83 ns simulation along the x, y, and z axis of the simulation box, as well as a total radius of gyration expressed as Rg, colored as blue, red, green, and black respectively. B) Intra-monomer interactions that influence the configuration of pMPC50 using radial distribution functions of atoms from given reference points. The reference point to atom RDF calculated, nitrogen to phosphorous, and choline to oxygen.
Finally, we find that the nitrogen atoms of MPC have an affinity for phosphorous atoms of neighboring monomers. Figure 1B shows two radial distribution functions (RDF) between 1) choline and phosphate oxygens and 2) between the choline nitrogen and phosphate phosphorous. The large peaks for both sets at 0.5 nm shows the intra-monomer interactions – the distance between these groups within a single monomer. However, the RDF calculation showed that the nitrogen and phosphorous atoms of neighboring monomers may have a high affinity for each other as demonstrated by the second set of peaks on the nitrogen to phosphorous data, which, due to the larger radius, may be between atoms of neighboring monomers. This supports separate monomers forming contacts between the phosphorous and nitrogen atoms.
Figure 2: A) A single monomer of pMPC shown with surrounding water molecules. Dashed lines indicate proximity to the phosphate oxygens in the middle of the monomer. B) Simulated pMPC50 at simulation time = 80 ns. C) Simulated pMPC50 at simulation time = 0 ns
Active roles of charged groups.
Figure 3 shows the radial distribution function between phosphate oxygens. This shows that the phosphate oxygens may be supportive of hydration shells at distances that correspond to the peaks on the graph. That is, at the distance where a peak occurs, the amount of water hydrogen or oxygen is high and a hydration shell exists.
The choline group, also analyzed in Figure 3, shows minimal binding to water at 0.2 nm. The inset of Figure 3 shows the solvent accessible surface area (SASA) for the choline groups and the phosphate oxygens of pMPC50 over the entire duration of the simulation. Measurements show that choline has a higher solvent accessible surface area than do the phosphate oxygens. Although the choline groups and the phosphate oxygens exist in pMPC at a 3:2 ratio, when the SASA measurements are reduced to a 1:1 ratio choline still exhibits a higher SASA than the phosphate oxygens.
Figure 3: Displays the roles of charged groups in terms of water attraction using a radial distribution function (RDF). Reference groups to atoms calculated are choline to oxygen and hydrogen atoms of water (Owater -Choline and Hwater-Choline respectively), and oxygen atoms of MPC to hydrogen and oxygen atoms of water (Hwater-O2 and Owater-O2 respectively). The inset displays solvent accessible surface area for two different reference groups of MPC. Solvent accessible surface area is shown in nanometers squared for the O2 groups and choline groups of MPC over the entire simulation.
After 60 ns of simulation time pMPC50 begins to adopt a folded configuration. The folding can be seen starting at 60 ns in Figure 1A and is qualitatively supported by the folded shape of pMPC50 in Figure 1B. Knowledge of this folded configuration may have implications to the biomechanics of pMPC brush systems where the configuration of the polymer may be important in the function of the system. Not only does the overall shape of the polymer seem to fold, but it is likely that monomers bend to some degree as well as bind together. However, the RDF calculation showed that the nitrogen and phosphorous atoms of neighboring monomers have a high affinity for each other. Shown in Figure 1B, the first large set of peaks show interactions within the same monomers due to their small radius of 0.5 nm on the x axis. The smaller set of peaks existing between 0.5 and 1 nm show interaction between phosphorous and nitrogen atoms of neighboring monomers due to their larger radius. The measurements from the RDF support that the phosphorous, at the center of the phosphate oxygens, and the nitrogen, at the center of the choline groups at the head of the monomer, have a strong affinity for each other. Since inter-monomer attractions are present, the data shows that there may be attraction between neighboring monomers which would facilitate a durable and tightly packed polymer system, this could contribute to decreased wear during shearing of MPC brush layers. However, investigations of the full polymer brush layers will be required to confirm this hypothesis. Visually and based on the RDF peaks, the phosphate oxygens are oriented outwards while the choline groups fold onto neighboring phosphate groups. An outwards orientation of the phosphate oxygens contributes to their ability to form hydration shells.
Phosphate oxygens of MPC monomers are hydrogen bonded to the hydrogen atoms of the surrounding water, thus forming the tightest layer of the hydration shells. The data shown in Figure 3 of the RDF function for water-O2 displays large peaks, with the graph of Hwater-O2 having the highest peaks, showing that there is strong hydrogen bonding between water and phosphate oxygens. This suggests hydration shells existing with respect to the distances corresponding to the peaks on the graph around the phosphate oxygens and the phosphate group, thus hydration shells may exist around individual monomers as water is preferentially localizing itself within proximity of the phosphate oxygens. This indicates that the phosphate group of pMPC50 which houses the two active phosphate oxygen atoms could play a large role in the experimentally observed hydration lubrication of MPC brush layers. Conversely, the choline groups, which are displayed on the same graph show barely any peaks and is indicative of there being much less powerful water-choline bonding in pMPC50. Furthermore, the choline groups have a higher SASA, indicated by the inset of Figure 3. This gives some insight to the extent to which the phosphate group has affinity for the surrounding solvent as the RDF indicates that although water has less SASA; the phosphate oxygens have more bonds to water than the choline does.
It appears that the phosphate groups of pMPC50 have the greatest affinity to the solvent and the highest influence on solvent interaction due to their likelihood to form strong hydrogen bonds and their outwards orientation towards the solvent in qualitative observation. It is more likely that the choline groups play a greater role in the mechanical stability of the polymer during folding due to their central nitrogen atom interacting with the phosphate groups within the polymer. Future directions include analyzing the structure to see if the same folding occurs among multiple polymers in a brush system. An understanding of the structural dynamics and hydration structure of a brush system could be used in the application of these polymers in order to utilize their frictional properties. Atomistic knowledge of the dynamics and hydration structure of pMPC could lead to more targeted applications and the discovery of other zwitterionic polymers with enhanced performance. Ideally, the enhanced performance of zwitterionic polymers can translate to biomedical advancements. To date, pMPC is the zwitterionic polymer primarily studied for such applications. By utilizing similar simulation techniques to analyze different zwitterionic molecules we hope to accelerate the discovery and design of new materials for biolubrication.
I would like to thank the members of the Multiscale Molecular Simulation Laboratory at Vanderbilt University for their mentorship and help with developing this study. Additionally, thanks go to Dr. Mulligan and Dr. Howdyshell for the opportunity to conduct research through the REHSS program.
Figure S1: The above graph shows the mean square deviation of the bulk water in the simulation. The diffusivity constant (m2 s-1) for the bulk water in the simulation is shown alongside plotted data.
Figure S2: Shows the root mean square deviation of the pMPC50 from the beginning of the simulation to the end of the simulation.
1Banquy, X., Burdyńska, J., Lee, D. W., Matyjaszewski, K., & Israelachvili, J. (2014). Bioinspired bottle- brush polymer exhibits low friction and Amontons-like behavior. Journal of the American Chemical Society, 136(17), 6199-202.
2Holz, M., Heil, S. R., & Sacco, A. (2000). Temperature-dependent self-diffusion coefficients of water and six selected molecular liquids for calibration in accurate 1H NMR PFG measurements. Physical Chemistry Chemical Physics, 2(20), 4740-4742.
3Jorgensen, W. L., Jorgensen, W. L., Maxwell, D. S., Maxwell, D. S., Tirado-Rives, J., & Tirado-Rives, J. (1996). Development and Testing of the OLPS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc., 118(15), 11225-11236.
4Kaminski, G. a, Friesner, R. a, Tirado-rives, J., & Jorgensen, W. L. (2001). Comparison with Accurate Quantum Chemical Calculations on Peptides †. Journal of Physical Chemistry B, 105(28), 6474–6487.
5Klein, C., Iacovella, C., McCabe, C., & Cummings, P. T. (2015). Tunable Transition from Hydration to Monomer-Supported Lubrication in Zwitterionic Monolayers Revealed by Molecular Dynamics Simulation. Soft Matter, 11, 3340-3346.
6Klein, J. (2013). Hydration lubrication. Friction, 1(1), 1-23.
7Klein, J. (2006). Molecular mechanisms of synovial joint lubrication. Proceedings of the Institution of Mechanical Engineers Part J Journal of Engineering Tribology, 220(8), 691-710.
8Kobayashi, M., Terayama, Y., Hosaka, N., Kaido, M., Suzuki, A., Yamada, N., … Takahara, A. (2007). Friction behavior of high-density poly(2-methacryloyloxyethyl phosphorylcholine) brush in aqueous media. Soft Matter, 3(6), 740.
9Kyomoto, M., Moro, T., Takatori, Y., Kawaguchi, H., & Ishihara, K. (2011). Cartilage-mimicking, high- density brush structure improves wear resistance of crosslinked polyethylene: A pilot study. Clinical Orthopaedics and Related Research, 469(8), 2327-2336.
10Mcgibbon, R. T., Beauchamp, K. A., Schwantes, C. R., Wang, L., Hern, C. X., Herrigan, M. P., … Pande, V. S. (2014). MDTraj : a modern , open library for the analysis of molecular dynamics trajectories, 9–10.
11Pronk, S., Pańll, S., Schulz, R., Larsson, P., Bjelkmar, P., Apostolov, R., Lindahl, E. (2013). GROMACS 4.5: A high-throughput and highly parallel open source molecular simulation toolkit. Bioinformatics, 29(7), 845-854.
12Wang, L. P., Martinez, T. J., & Pande, V. S. (2014). Building force fields: An automatic, systematic, and reproducible approach. Journal of Physical Chemistry Letters, 5(11), 1885-1891.