Kinetics of microribbon formation in a simplified model of amelogenin biomacromolecules

Kinetics of microribbon formation in a simplified model of amelogenin biomacromolecules

Wei Li Department of Physics, Lehigh University, Bethlehem, PA 18015    Ya Liu Department of Physics, Lehigh University, Bethlehem, PA 18015    Toni Perez Department of Physics, Lehigh University, Bethlehem, PA 18015    J. D. Gunton Department of Physics, Lehigh University, Bethlehem, PA 18015    C. M. Sorensen Department of Physics, Kansas State University, Manhattan, KS 66503    A. Chakrabarti Department of Physics, Kansas State University, Manhattan, KS 66503
July 15, 2019

We show that the kinetics of microribbon formation of amelogenin molecules is well described by a combination of translational and rotational diffusion of a simplified anisotropic bipolar model consisting of hydrophobic spherical colloid particles and a point charge located on each particle surface. The colloid particles interact via a standard depletion attraction while the point charges interact through a screened Coulomb repulsion. We study the kinetics via a Brownian dynamics simulation of both translational and rotational motions and show that the anisotropy brought in by the charge dramatically affects the kinetic pathway of cluster formation and our simple model captures the main features of the experimental observations.

To be submitted to PRL

The self-assembly of colloidal particles and globular proteins with inhomogenous properties into various desired structures is an important area of current research, in which the inhomogeneity (such as an amphiphilic group and a residual group) plays a crucial role in the kinetic pathway of cluster formation and eventually impacts on the stabilized crystal structures. Examples include Janus particles and DNA-coated colloidal particles that lead to the next generation of building blocks of new materials and have potential applications in fabricating phontonic crystals, targeted drug delivery, and electronic equipment Qian11 (); Flavio11 (); Hong08 (); Mirkin10 (). Another important example is amelogenin, the chief hydrophobic protein in enamel matrix with a hydrophilic 25-amino acid C terminus. Amelogenin is involved in the mineral deposition and is responsible for the major structural process during enamel biomineralization. A recent experiment shows that the hydrophilic C-terminus is essential for the self-assembly of amelogenin into microribbons, which may also be relevant to the elongated and oriented growth of apatite crystals during biomineralization Du05 (); Ame07 (). The experiment also revealed that the microribbon is formed by a hierarchical organization of amelogenin molecules into a chain of nanospheres.

In this article we present a coarse-grained model of amelogenin and explain how the microribbon forms through self-assembly. The amelogenin molecule is hydrophobic with a charged hydrophilic tail. Our model describes this in a simplified way, representing the monomer as a spherical molecule, with the charged hydrophilic tail replaced by a single tethered point charge located on the surface of the molecule. Using Brownian dynamics simulations, we investigate the static and dynamic properties of the self-assembly process. We show that the anisotropy brought by the charge dramatically affects the kinetic pathway of cluster formation and our simple model captures the main features of the experimental observations in the early stages.

In the experimental studies of amelogenin assembly process, one typically adds salt and a precipitant such as polyethylene glycol (PEG). Although PEG-protein interactions are quite complicated prl06 (), we model this as a depletion interaction that we take to be given by the Asakura-Oosawa potential plus a repulsive hard-core-like interaction, depending on the center-center distance between spherical particles, , where


The cut-off range , where is the size-ratio between a PEG coil and a colloidal amelogenin particle that controls the range of the depletion interaction, and is the value of the PEG volume fraction that controls the strength of the interaction described by the absolute value of the minimum potential depth AO54 (); Vrij76 (). The hard core potential is given by


We set , since the values of have been reported to lead to anomalies when a mimic of the hard-core potential is required in the potential A36 (). The point charges interact with each other through a screened Coulomb potential,


in which magnitude controlled by and range controlled by Debye screening length , are varied in the simulations. We find that the early morphology of the self-assembly is more sensitive to the size of rather than . These Coulomb charges exert torques on adjacent molecules and hence produce a rotational motion of the molecules that is included in our Brownian equations of motion, read as:


where , , , , are the mass, moment of inertia, position vector, angular velocity, and torque, respectively, of the ith colloidal particle. The mass of the point charge is ignored in this model. () is the translational (rotational) friction coefficient. and are the random forces and torques acting on the ith colloidal particle respectively, which satisfy a fluctuation-dissipation relation , VG82 (). We choose , and the time step in reduced time units of with Li11 (). The other parameters are chosen according to the range of values of the experimental data displayed in Du et al’s paper Du05 (). In particular, here the results shown are for the values =0.1, , (all length scales are measured in units of monomer diameter and energies are scaled by ). We consider a small volume fraction of the molecules, typically of the order of or . Periodic boundary conditions are enforced to minimize wall effects. All simulations start from a random initial monomer conformation and the results for the kinetics are averaged over several (5-10) runs.

We show typical configurations in the early stage of microribbon formation in Fig. 1. We see that initially the monomers form oligomers and the oligomers self-assemble into larger aggregates contained (approximately) within nanospheres, then nanospheres associate together into a chain structure, quite similar to the process of amelogenin self-assembly (see Fig. 3 in reference Du05 ()).

Figure 1: (color online). (a) The model of amelogenin molecule consists of a spherical colloid particle and a point charge located on its surface that preserves a bipolar nature. (b) and (c) Oligomerization of amelogenin molecules occurs by means of hydrophobic interactions and modification due to Coulomb repulsions. (d) Nanosphere structures are formed through aggregation of monomers and oligomers. (e) Further association of nanospheres results in larger assemblies among which chains are formed.

The Coulomb forces cause the majority of the charges to point outwards inside the nanospheres, as shown in Fig. 2. This again is reminiscent of the situation with the hydrophilic tails of amelogenin. To characterize the orientational distribution of point charges within a nanosphere, we introduce a parameter , defined as , where is the position vector of a point charge referenced to the center of its host colloid and the position vector of a colloid referenced to the center of mass of this nanosphere. Therefore if the charge points outwards and if it points inwards. Due to fluctuations and our choice of short repulsion interaction range, the distribution of is spread around a small nonzero peak position, as shown in Fig. 2, for the configuration in inset. Subsequent aggregation forms via necks between nanospheres, leading to larger structures that eventually form (flexible) microribbons.

Figure 2: (color online). Distribution of ; shows whether the point charge orients outwards () or inwards (). Inset: Nanosphere structure, large (grey) spheres indicate amelogenin molecules while small (red) ones for point charges.

We characterize the morphology of the clusters in terms of their fractal dimension . The q-dependence of shown in Fig. 3(a) is given by a slope of -2 on a log-log plot over a reasonable range of q. This value of is larger than the typical diffusion-limited-cluster-cluster aggregation (DLCA) value of 1.8 as the ribbon-like clusters observed here have significant short range ordering Porod (). For this reason we speculate that the repulsive Coulomb force causes the surface particles to reorganize, which reduces the surface roughness of our model, in contrast to the case without the Coulomb interactions (see reference Li11 ()). We find that a peak develops in the structure factor whose magnitude increases with time. Correspondingly, the peak position decreases as a function of wave number with increasing time; these features are typical of spinodal decomposition Gunton (). However, we note that for deep quenches that lead to cluster growth considered in this model, the system is controlled by two characteristic lengths that evolve differently in time and an apparent scaling for the structure factor can only be observed over some period of time when these two characteristic length scales become comparable to each other Amit04 (); PRE98 ().

Figure 3: (color online). Log-log plot of structure factor at several times. The dashed line indicates fractal clusters with , where , while the dotted line indicates the Porod regime , .

The kinetics of the cluster growth process in this simple model for amelogenin self-assembly is consistent with a cluster-cluster aggregation mechanism, as shown in Fig. 4. The number of clusters decreases inversely with time (i.e. the kinetic exponent ) and the radius of gyration increases as a power law with an exponent of KF (); Amit (). The relation between and involves the fractal dimension in the following way: . Thus the kinetic exponents are consistent with a fractal dimension of as well.

Figure 4: (color online). Plot of (a) number of clusters, and (b) radius of gyration, as a function of time in a log-log plot. The dotted line in (a) has a slope of ; the dashed line in (b) has a slope of .

In order to check the stability of the clusters produced in straight simulation from a quench () into the two phase gas-solid region, we “heat” the system formed after steps from the initial quench to and run the simulation for another time steps. This heating corresponds to the change in PEG concentration as it is done in the experiment. A typical run is shown in Fig. 5. As shown in the figure, we find that the microribbons are quite robust under such a change of conditions, as found experimentally for a wide range of conditions (including a change in the PEG concentration). We have also checked the structure factor behavior for the heating process and found that the structure factor curve has the same exponent of -2 as at the cooler temperature (see Fig. 3). This is consistent with the robustness of the major cluster morphology.

Figure 5: (color online). (a) Morphology of entire system (, ) from another time steps shallow quench of after first quench into , time steps. (periodic boundary condition are enforced). (b) One of the clusters inside (a).

We characterize these microribbons by computing three eigenvalues of the gyration tensor, denoted as , , and in descending order gt85 (). We calculate the ratios among these three principal values, which characterize the relative colloid (or charge) particle distribution. A typical result is given in Table 1, for the colloid particles and point charges separately, based on the configuration in Fig. 5(b). These ratios of eigenvalues demonstrate the anisotropy of these clusters, suggesting that the cluster is ribbon-like. The difference of the ratios between colloidal particles and point charges may indicate that under electrical repulsion, the rotation of monomers inside the cluster leads to the small reduction of anisotropy in the point charges’ spatial distribution.

Colloid particles 9.57 9.03
Point charges 9.25 8.74
Table 1: Typical ratios of eigenvalues of gyration tensor

From the gyration tensor we can determine the typical length scales of these clusters by taking square roots of the eigenvalues. To illustrate this, we use the colloid configuration shown in Fig. 5(b), which gives , , and for the three length scales. Since the chain structure is the association of several nanospheres (as shown in Fig. 1), the typical diameter of the nanosphere in this configuration is at the most . We can also use this estimate of the diameter to obtain an estimate of the minimum wavenumber that corresponds in Fig. 3 to the crossover to the Porod regime (corresponding to compact clusters on a local scale). This yields a value of , which is consistent with the behavior shown in Fig. 3(b). Since the hydrodynamic radius of the amelogenin macromolecule is about nm, the of the nanosphere for the configuration in Fig. 5(b) is nm; this lies in the range of values of the radii of typical nanospheres, nm, found in the experiment (cf. Fig. 3, reference Du05 ()).

In summary, we have developed an anisotropic, bipolar model for the hierarchical self-assembly of amelogenin molecules and have carried out Brownian dynamics simulations of the self-assembly process. Simulations show a hierarchical self-assembly process where the molecules aggregate to form dimers, hexamers to nanospheres, and the assembly of the nanospheres then lead to the formation of microribbons in agreement with experimental findings. The relative strengths of the interaction parameters can lead to a phase diagram where the nanospheres are stable against microribbon formation. Thus a control over the morphology of the clusters can be achieved. Such a study is currently underway.

This work was supported by grants from the Mathers Foundation and the National Science Foundation (Grant DMR-0702890). Work at Kansas State University was supported by NSF NIRT grant CTS0609318.


  • (1) Qian Chen, Sung Chul Bae, and Steve Granick, Nature 469, 381, (2011).
  • (2) Flavio Romano and Francesco Sciortino, Nature material 10, 171, (2011).
  • (3) L. Hong, A. Cacciuto, E. Luijten, and S. Granick, Langmuir 24, 621-625 (2008).
  • (4) Chad A. Mirkin, MRS Bulletin 35, 532-539 (2010).
  • (5) S. Gajjeraman, K. Narayanan, J. Hao, C. Qin and A. George, J. Biological chemistry 282, 1193 (2007).
  • (6) C. Du, G. Falini, S. Fermani, C. Abbott, J. Moradian-Oldak, Science 307, 1450 (2005).
  • (7) J. Bloustine, T. Virmani, G. M. Thurston, and S. Fraden, Phys. Rev. Lett. 96, 087803 (2006).
  • (8) S. Asakura and F. Oosawa, J. Chem. Phys. 22, 1255 (1954).
  • (9) A. Vrij, Pure Appl. Chem. 48, 471 (1976).
  • (10) K. G. Soga, J. R. Melrose, and R. C. Ball, J. Chem. Phys. 108, 6026 (1998); K. G. Soga, J. R. Melrose, and R. C. Ball, J. Chem. Phys. 110, 2280 (1999).
  • (11) W. F. Van Gunsteren, and H. J. C. Berendsen, Mol. Phys. 45, 637 (1982).
  • (12) Wei Li, J. D. Gunton, Siddique J. Khan, J. K. Schoelz, and A. Chakrabarti, J. Chem. Phys. 134, 024902 (2011).
  • (13) G. Porod, Zeit. Kolloid, 124, 53 (1951). See also, O. Glatter and O. Kratky. Small-angle X-ray Scattering, (Academic Press, London 1982).
  • (14) J. D. Gunton, M. San Miguel, and P. S. Sahni, in Phase Transitions and Critical Phenomena, edited by C. Domb and J. L. Lebowitz (Academic, London, 1983), Vol. 8.
  • (15) J. J. Cerda, T. Sintes, C. M. Sorensen, and A. Chakrabarti, Phys. Rev. E 70, 051405 (2004).
  • (16) H. Huang, C. Oh, and, C.M. Sorensen, Phys. Rev. E 57, 875 (1998).
  • (17) S. K. Friedlander. Smoke, Dust, and Haze: Fundamentals of Aerosol Behavior (Wiley, New York, 1977).
  • (18) S. J. Khan, C.M. Sorensen, and A. Chakrabarti, J. Chem. Phys. 131, 194908 (2009).
  • (19) Doros N. Theodorou and Ulrich W. Suter, Macromolecules, 18, 1206 (1985).
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description