Brownian dynamics simulation of analytical ultracentrifugation experiments
© Diez et al; licensee BioMed Central Ltd. 2011
Received: 11 February 2011
Accepted: 2 March 2011
Published: 2 March 2011
We have devised a protocol for the Brownian dynamics simulation of an analytical ultracentrifugation experiment that allows for an accurate and efficient prediction of the time-dependent concentration profiles, c(r, t) in the ultracentrifuge cell. The procedure accounts for the back-diffusion, described as a Brownian motion that superimposes to the centrifugal drift, and considers the sector-shaped geometry of the cell and the boundaries imposed by the meniscus and bottom.
Simulations are carried out for four molecules covering a wide range of the ratio of sedimentation and diffusion coefficients. The evaluation is done by extracting the molecular parameters that were initially employed in the simulation by analyzing the profiles with an independent tool, the well-proved SEDFIT software. The code of simulation algorithm has been parallelized in order to take advantage of current multi-core computers.
Our Brownian dynamics simulation procedure may be considered as an alternative to other predictors based in numerical solutions of the Lamm equation, and its efficiency could make it useful in the most relevant, inverse problem, which is that of extracting the molecular parameters from experimentally determined concentration profiles.
Since the invention of the analytical ultracentrifuge by Svedberg , the technique of analytical ultracentrifugation (AUC) has been a classical - and, thanks to advances in instrumentation and analysis software, it is still a most modern - technique for characterization of macromolecules and nanoparticles in solution. The reader may grasp the recent importance of this field in monographs [2–4] and thematic issues of other journals [5–7].
where k B is Boltzmann constant, T is the absolute temperature and R = k B N A is the constant of perfect gases.
Eq. 6 is written in cylindrical coordinates, because the AUC geometry is radial, and sector-shaped cells are used (devised so that the radial trajectories would not collide with the lateral walls).
In spite of the basic nature of the concepts involved, the solution of the Lamm partial differential equation, with the geometry and boundary conditions of the AUC experiment is extremely difficult, and requires numerical methods [10, 11] based on time discretization and a finite-elements description of the AUC cell. . Nonetheless, in modern AUC analysis programs, like SEDFIT , a Lamm-equation solver is embodied, enabling the prediction of computed c(r, t) for estimations of s, D and M, which are to be optimized as to fit the experimental c(r, t). In these procedures (apart from the fitting or optimization algorithms, a central piece is the prediction of c(r, t).
In this paper we investigate an alternative predictor of AUC concentration profiles. Instead of starting from the balance of the macroscopic flows established by the Lamm equation, we consider a microscopic description of the motion of particles under the simultaneous effect of a deterministic force, and the random forces characteristic of Brownian motion, so that the latter replaces the Fick's law description of macroscopic diffusion. In a simple (and somewhat naïve) approach, we formulate a Brownian dynamics algorithm to simulate the trajectories of particles. With our microscopic perspective, our method also discretizes time (Brownian simulation steps) and instead of finite elements use discrete particles. Carrying out such simulation for a sufficiently large number of particles, we can determine the time-dependent concentration profile in the AUC cell. In the next section we describe the procedure and demonstrate adequacy of its results, and finally we shall discuss on its performance, eventual advantages and possible extensions for further applications.
The simulated system
Note that, apart from several constants relative to the instrument or the simulation the concentration is not determined by n(i, j), but by n(i, j)/r i . This takes into account what is called in the AUC terminology, the radial dilution effect.
What we simulate corresponds to a highly diluted solution, in which there are no particle interactions. Thus we can generate trajectories of individual particles, independent of each other.
Brownian dynamics simulation
where δr sed = ω2srδt is the deterministic sedimentation drift of the particle with instantaneous, position-dependent velocity ω2sr, while the random Brownian displacement has zero mean and variance .
Thus the algorithm based on eqs. 11, 12 and 13 could be applicable to arbitrarily large time steps (even as large as the time interval τ between registers). This is essentially true if there were no end effects, i.e., in infinite, unbound AUC cell. For the sake of simplicity, we still adopt the simple criteria that particles stop at the meniscus and bounce at the bottom. Thus the only defect introduced by this procedure would be an inaccurate prediction of the concentration near the meniscus and bottom. In this regard, we note that the end-effects also affect other prediction procedures, like those based in Lamm-equation solvers, and influence the experiment itself, so that it is a common practice to discard the two terminal regions in the analysis of AUC experiments.
Summarizing from the previous description, Brownian dynamics trajectories are simulated for a large number of particles, N part . The trajectory of one particle is monitored, determining at successive times t j the interval of radial position r j . Then the counter for those interval and position is increased n(i, j) → n(i, j) + 1.
If t run is the total time of the experiment (usually, several hours), the simulation for each particle consists of a number of steps N steps = t run /Δt. For recording purposes, a number of registers (scans in the AUC experiment), equal to N scans T = t run /τ is made, observing the radial distance, and therefrom the index i, of the slice at which the particle is at that instant, j = t run /τ. Then one unit is added to the counter, n(i, j). Trajectories are simulated for a sufficiently large number of particles N part . At the end of the simulations, concentration profiles c(r i , t j ) are obtained from the n(i, j) counter using eq. 9.
The data employed in the simulation reflect those of real AUC experiments. The geometry of the cell is given by r m = 5.80 cm, r b = 7.20 cm and sector angle ϕ = 3 degrees (0.05 radians). The rotor speed ω and the duration of the experiment t run was varied depending on the molecule being simulated, and the mode (velocity or equilibrium) of the experiment. In some cases also the meniscus position was varied for equilibrium experiments. We found that N part = 105 particles suffices to obtain a rather low level of noise as it will be shown below.
t run (secs.)
10−3M calc (g/mol)
23 278 (-10)
s calc (S)
107D calc (cm3/g)
t run (secs.)
10−3M calc (g/mol)
Results and discussion
In order to evaluate the accuracy of the results for the prediction of sedimentation profiles in the velocity mode, we have compared them with calculations with the well-known SEDFIT software . This tool, which includes a sophisticated Lamm-equation solver has an utility for predicting the concentration profiles. Comparison of our simulated profiles with the SEDFIT predictions resulted in a good agreement, as shown in Figure 2. The determination of concentration profiles has the final purpose of determining the molecular parameters of interest, so the most relevant evaluation of their prediction is the confirmation of whether their analysis provides correct values for those parameters. In order to do so with our simulation results, we employed the analysis tool of SEDFIT as an independent and robust criterion. Indeed, SEDFIT uses a reliable Lamm-equation solver to determine the parameters by means of powerful optimization algorithms. In the single-species mode of SEDFIT, the program provides the values of the molecular weight, M, and the sedimentation and diffusion coefficients, s and D of the sample. We made the SEDFIT analysis with typically 50 scans (t values) each covering 1400 radial positions r values.
In Table 1 we report the s, M, D and values for the four samples considered in our study. We note the high accuracy of the recovered values of these three parameters, reflected in their very small deviations (particularly in the case of s) from the values employed as input in the simulation. The agreement is good for the four cases, which, as indicated above, cover a wide range of the sedimentation-to-diffusion ratio and molecular mass, including in the velocity experiments a molecule as small as cyclodextrin.
which, as usual, is linearized in the form of a plot of ln c(r) vs. whose slope is , from which the molecular weight M is extracted. The results are listed in Table 2. Here we observe that the recovered values of M are not as good as those obtained in the velocity mode, with deviations of about 15%. Among other reasons (like the smaller amount of information resulting from equilibrium experiments, and the errors inherent to the determination of the slope in the fit to the linearized equation) this may be because the simple eq. 17 neglects radial dilution.
The extreme simplicity of the algorithm that simulates the trajectory of one molecule makes the simulation scheme very well adapted for parallelization, thus taking benefit of present multi-core platforms, because each trajectory can be generated in a separate thread/core. We have implemented OpenMP directives in our Fortran 90 code and tested the performance in a DELL T5500 workstation with two Intel Xeon X5660 processors. 1000 simultaneous simulation with N part = 105, N r = 100 radial positions and N t = 50 recorded scans took 46.5 CPU seconds. Broadly speaking, our algorithm, which is easily parallelized, is able to run one thousand c(r, t) calculations in less than one minute. The computing speed of the algorithm is crucial in its main use, namely, in the analysis of c(r, t) experimental profiles by any kind of fitting to computed profiles, and the speed of our procedure seems suitable for that purpose.
A nice feature of the Brownian simulation of ultracentrifugation is that it allows to visualize the simulated trajectories using computer graphics. This may be of utility for demonstrative purposes (e.g. in teaching AUC principles). We have produced two videos, showing the evolution of the solute particles in the AUC as the sedimentation proceeds, one at high rotor speed, in the velocity mode, and another at low speed, that reaches the sedimentation-diffusion equilibrium. These videos accompany this paper as additional files 1 and 2.
In this work we have presented a simple procedure, based on a Brownian dynamics, microscopic simulation, for predicting the time/position-dependence of concentration (concentration profiles c(r, t)) during the AUC experiment, which can be regarded as an alternative to the numerical solution of the macroscopic Lamm equation. The correctness of the procedure has been tested comparing the profiles with those computed by the Lamm-equation solution, and the concordance of the molecular parameters recovered from them with those used in the simulation.
Having presented this proof-of-concept, some advantages of our scheme can be hinted - although they all remain to be evaluated in future work. An important aspect to be considered is computational efficiency. As commented above, our computational procedure has the feature of being perfectly parallelizable. A comparison with the numerical solution of Lamm solution requires further the labor, implementing codes for those procedures and making a side-to-side analysis of computing speed and requirements. Further work is planned in this regard.
Apart from computing efficiency, the BD scheme has the potential advantage that it can treat easily cases of arbitrary complexity. In our proof-of-concept we have restricted ourselves to the simplest case of identical, non interacting particles (extension to different but still non-interacting particles is trivial). An extraordinary utility of AUC is the characterization of macromolecular interactions. For problems with interacting particles, the BD scheme can be easily adapted. In BD, detailed interactions between the particles can be modeled whereas the continuum method only allows for averaged density dependent potentials. Even if computing efficiency would suffer in such, more complex problems, still the BD approach, which is based on first principles and allows explicit description of interaction between particles (or the effect of special conditions in the AUC experiments) may be a valuable tool for testing other approaches.
Fortran 90 source code with OpenMP directives will be freely downloaded from our web site, http://leonardo.inf.um.es/macromol/
This work was performed within a Grupo de Excelencia de la Región de Murcia (grant 04531/GERM/06). Support also provided by grant CTQ-2009-08030 from Ministerio de Educación y Ciencia (MEC), including FEDER funds. A.O. acknowledges a postdoctoral fellowship from Universidad de Murcia, and A.I.D. is recipient of a predoctoral fellowship from MICINN.
- Svedberg T, Pedersen K: The Ultracentrifuge. 1940, Oxford: Clarendon PressGoogle Scholar
- Mächtle W, Börger L: Analytical Ultracentrifugation of Polymers and nanoparticles. 2006, Berlin, Heidelberg: Springer-VerlagGoogle Scholar
- Harding S, Rowe A, Horton J: Analytical Ultracentrifugation: Techniques and Methods. 1992, Cambridge, Chester: Royal Society of ChemistryGoogle Scholar
- Scott D, Harding S, Rowe A: Analytical Ultracentrifugation: Techniques and Methods. 2005, Cambridge: Royal Society of ChemistryView ArticleGoogle Scholar
- Solovyova A: 17th international symposium on analytical ultracentrifugation and hydrodynamics. Eur Biophys J. 2010, 39: 345-346. 10.1007/s00249-009-0406-4.View ArticleGoogle Scholar
- Cölfen H: Analytical Ultracentrifugation. Macromol Biosci. 2010, 10: 687-688.View ArticleGoogle Scholar
- Zhao J, Schuck P: Methods.
- van Holde KE, Johnson WC, Ho PS: Physical Biochemistry. 1998, New Jersey: Prentice HallGoogle Scholar
- Lamm O: Die Differentialgleichung der Ultrazentrifugierung. Ark Mat Astr Fys. 1929, 21B: 1-4.Google Scholar
- Claverie JM, Dreux H, Cohen R: Sedimentation of generalized systems of interacting particles. I. Solution of systems of complete Lamm equations. Biopolymers. 1975, 14: 1685-1700. 10.1002/bip.1975.360140811.View ArticleGoogle Scholar
- Dishon M, Weiss G, Yphantis D: Numerical solutions of the Lamm equation. I. Numerical procedure. Biopolymers. 1966, 4: 449-455. 10.1002/bip.1966.360040406.View ArticleGoogle Scholar
- Schuck P: Size distribution analysis of macromolecules by sedimentation velocity ultracentrifugation and Lamm equation modeling. Biophys J. 2000, 78: 1606-1619. 10.1016/S0006-3495(00)76713-0.View ArticleGoogle Scholar
- Pavlov G, Korneeva E, Smolina N, Shubert U: Hydrodynamic properties of cyclodextrin molecules in dilute solutions. Eur Biophys J. 2010, 39: 371-379. 10.1007/s00249-008-0394-9.View ArticleGoogle Scholar
- Müller J: Prediction of the rotational diffusion behavior of biopolymers on the basis of their solution or crystal structure. Biopolymers. 1991, 31: 149-160.View ArticleGoogle Scholar
- Luo Z, Zhang G: Scaling for Sedimentation and Diffusion of Poly(ethylene glycol) in Water. J Phys Chem. 2009, 113: 12462-12465.View ArticleGoogle Scholar
- Crothers D, Zimm B: Viscosity and Sedimentation of DNA from Bacteriophages T2 and T7 and Relation to Molecular Weight. J Mol Biol. 1965, 12: 525-536. 10.1016/S0022-2836(65)80310-2.View ArticleGoogle Scholar
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.