Tracking Particles with Large Displacements using Energy Minimization

Tracking Particles with Large Displacements using Energy Minimization

Rostislav Boltyanskiy Department of Physics, Yale University, New Haven, Connecticut 06520, USA    Jason W. Merrill Desmos, Inc., San Francisco, California 94103, USA    Eric R. Dufresne Department of Materials, ETH Zürich, 8092 Zürich, CH.
July 13, 2019

We describe a method to track particles undergoing large displacements. Starting with a list of particle positions sampled at different time points, we assign particle identities by minimizing the sum across all particles of the trace of the square of the strain tensor. This method of tracking corresponds to minimizing the stored energy in an elastic solid or the dissipated energy in a viscous fluid. Our energy-minimizing approach extends the advantages of particle tracking to situations where particle imaging velocimetry and digital imaging correlation are typically required. This approach is much more reliable than the standard squared-displacement minimizing approach for spatially-correlated displacements that are larger than the typical interparticle spacing. Thus, it is suitable for particles embedded in a material undergoing large deformations. On the other hand, squared-displacement minimization is more effective for particles undergoing uncorrelated random motion. In the Supplement, we include a flexible MATLAB particle tracker that implements either approach. This implementation returns an estimation of the strain tensor for each particle, in addition to its identification.

Valid PACS appear here


In a wide range of pure and applied sciences, the motion of objects needs to be tracked over time. The central difficulty is that while an object’s trajectory is continuous through space and time, its position can only be sampled at discrete moments. When the objects of interest are far away from each other or otherwise distinguishable, tracking is simple. However, when many identical objects (which we refer to as particles) are near each other, tracking can be difficult.

Particle tracking is employed across physics, biology, and many other fields. In soft matter physics, it has shed light on mechanical properties like rheology Squires and Mason (2010), and revealed microscopic processes underlying phase transitions Peng et al. (2011), among many other applications. In biology, particle tracking is used across many length scales from the movement of single molecules, to the transport of organelles within cells, to movement of whole organisms such as flies, birds, fish, and humans Ardekani et al. (2013); Straw et al. (2011); Katz et al. (2011); Rodriguez et al. (2011).

In the absence of Brownian motion, particles embedded in a material reveal its deformation field. In fluid mechanics, particle tracking reveals flow fields in the Lagrangian frame of reference (Xu et al. (2006)). Traction force microscopy (TFM) quantifies the forces applied to a solid surface by tracking and analyzing motion of particles embedded within the solid. TFM was originally developed to quantify the forces exerted by adherent cells (Munevar et al. (2001)), but has recently been applied to a wide range of problems in biology and physics Style et al. (2014).

The standard method of particle tracking assigns identities by minimizing the sum of squared displacements across time points Crocker and Grier (1996). While this method is rigorously correct for objects undergoing Brownian motion, it works very well across a wide range of applications, provided that the particles move a distance that is small compared to the typical inter-particle separation. The improvement of this basic tracking approach has been an active area of research in recent years, driven by applications in the biomedical community. For a useful comparison of these particle tracking approaches, see Chenouard et al. (2014). When the displacements are large, it is helpful to employ a tracking algorithm that exploits knowledge of the system’s kinematics. For example, particle tracking in dense turbulent fluid flows has been greatly improved by employing a “predictive tracker,” which exploits the inertial character of high Reynolds number fluid flow Ouellette et al. (2006).

Particle Imaging Velocimetry (PIV) and Digital Imaging Correlation (DIC) are widely used to assess the flow of fluids and deformation of solids Adrian and Westerweel (2010); Pan et al. (2009); Bar-Kochba et al. (2015). These methods assign velocities or displacements to locations in a material by cross-correlating patches of an image across time points. This is a robust and successful approach, even for systems with large displacements. However, it is inappropriate in cases where individual particle identities need to be followed over time or where the displacements need to be known at the resolution of individual particles.

Here, we extend particle tracking to the large-displacement regime where PIV and DIC are usually employed. Our algorithm is optimized for tracking the motion of particles undergoing large spatially-correlated displacements, typical for tracers embedded in a deformed elastic solid or flowing viscous fluid. We start by posing the tracking problem generally, and review the maximum likelihood method employed to track Brownian particles. Then, we describe our energy minimization approach and directly compare the results of the two methods for simulated and real data.

Optimal Assignment

Suppose a set of particles, , where is the particle index, are found at locations at time and a set of particles are found at at time . We wish to connect particle positions across time to make trajectories. The particle identities are represented by a list, , where if particle at time becomes particle at time .

Figure 1: Schematic of strain-based particle tracking. (a) Particles at time . Dotted circle represents a region around the particle of interest, , with radius within which the particle neighbors are considered. (b) Particles at time with candidate particles circled as in (a). (c) Particle and its neighbors overlapped with and its neighbors. Arrows show tracked displacements between neighbors of the candidate particles. (d) Particle and its neighbors overlapped with and its neighbors. Arrows show tracked displacements between neighbors of the candidate particles.

We find an optimal assignment of identities by minimizing a cost function. We associate a cost, , to each possible particle pairing. The total cost, , is the sum of the costs of individual pairings across all of the particles, . We assign identities that minimize the cost using the Hungarian algorithm. Kuhn (2005); Buehren (2004).

Since assigning particles is a combinatoric process, it rapidly becomes computationally overwhelming as the system size increases. To avoid this difficulty, it is essential to limit the number of combinations by ruling out unphysical assignments, as described below.

Minimization of Squared-Displacements

The standard particle-tracking algorithm minimizes the sum of the squared displacements across all of the particles Crocker and Grier (1996). This yields the most likely assignment of identities when particles undergo Brownian motion. In that case, the probability that the particle will be displaced by in a time interval, , is:


Here, is the diffusion coefficient and is the number of spatial dimensions Berg (1993). Therefore, the probability of observing a specific set of displacements of identical particles moving independently is:


The most likely assignment of particle identities maximizes the total probability, , which is equivalent to minimizing the total squared-displacement . Furthermore, we can assign costs for potential pairings as , and minimize the cost function , as defined in the previous section.

To accelerate the data analysis, one needs to rule out unphysical assignments. Typically, one specifies a maximum displacement a particle could have between time points, and assigns an infinite cost to all of the elements of that correspond to displacements that would exceed this maximum. These elements are ignored in the combinatorial optimization.

In practice, minimizing the squared-distance works very well whenever particle displacements are small compared to the distance of a particle to its nearest neighbors.

Minimization of Energy

In the case of large displacements, square-distance minimizing trackers may not work effectively. In many cases, however, these displacements may have strong spatial correlations. For example, particles embedded in a rigid material undergoing translation and rotation are fixed relative to their neighbors. Similarly, displacements of particles embedded in a liquid or solid undergoing large deformations are typically strongly correlated to their neighbors. In both these cases, we do not wish to penalize displacements in the laboratory frame, but displacements relative to neighboring particles.

In continuum mechanics, the variation of displacements over space is characterized by the strain. In the case of a linear, isotropic, elastic solid of Young’s modulus, , the stored energy is


where and are the strain and displacement fields, respectively. In component form, these two are related as Landau and Liftshitz (1986). Analogously, for a fluid of dynamic viscosity sheared with a strain rate , the rate of energy dissipation is Kundu and Cohen (2008):


If one is considering only two time points, then experimental estimates of and only differ by a pre-factor.

Motivated by these two physical examples, we propose to assign particle identities by minimizing the stored/dissipated energy. The basic hypothesis is that tracking errors tend to exaggerate the strain/strain-rate, increasing the apparent stored/dissipated energy. Therefore, we implement a cost function . The main challenge of the proposed algorithm is to identify the strain associated with a given particle pairing. Calculating the strain requires positions of particles and their nearest neighbors at two times.

The nearest neighbors of particle , can be found in 2-D using the Delaunay triangulation. To ensure that the strain is uniform over the region where it is calculated, we only consider neighbors within a distance of . In Fig. 1(a,b), neighbors included in the strain calculation are within the dotted circles.

We use the method of Falk and Langer Falk and Langer (1998) to calculate the strain associated with candidate particle pairs and their neighbors. The strain about a particle at position from the position of neighbors can be estimated as follows:


Note that in order to calculate strain accurately, at least neighbors must be included for particle tracking in -dimensions.

To properly calculate strain around a particle, one needs to accurately track its neighbors across time. Consider a particle at time , (in the center of the dotted circle in Fig. 1a) and a candidate corresponding particle at time , (in the center a dotted circle in Fig. 1b). To connect the neighbors of with those of , we translate the candidate particles and their neighbors to the same location, as shown in Figure 1(c,d). We then track the neighbors of the candidate pair by minimizing the square of the residual displacements. To accelerate the calculation, we rule out relative displacements that would exceed an expected maximum value of the strain. If the candidate pair produces at least tracked neighbors, we calculate the strain for the candidate pairing according to Eq. (8) and assign the pairing a cost . Otherwise, the pairing is ruled out by setting the cost to infinity. Once all possible pairings have been assigned costs, the Hungarian algorithm is used to minimize the total cost across all of the particles.

Implementation of tracking algorithms

We implement the squared-distance-minimizing “diffusion tracker” and energy-minimizing “strain tracker” algorithms in MATLAB. In the supplement, we include the main particle-tracking function, Tracker.m. We also include a script, Example.m, that allows the user to explore the examples described below. Details regarding use of the code are found in the comments as well as the text of ReadMe.doc. A convenient by-product of our tracker is that it returns the symmetrized strain matrix for each tracked particle.

Comparison of tracking strategies with simulated data

Figure 2: Comparison of the diffusion and strain trackers for simulated data of particles moving by diffusion (first column), translation (second column), shear (third column), and stretch (fourth column). The particle positions at the first time point (black dots) and second time point (magenta dots) are shown in the first (top) row. The correct displacments (green arrows) are shown in the second row. The displacments calculated from identities returned by the diffusion and strain trackers are shown in the third and fourth rows, respectively. In all cases, arrows are drawn to scale. In the third and fourth rows, the numbers in gray boxes correspond to percentages of particle trajectories calculated correctly and incorrectly by each tracker. The percentages do not add to 100% because the trackers were not able to find partners for some particles.

We compare the performance of the diffusion and strain trackers with simulated data. As shown in Figure 2, we consider four types of particle motion: diffusion, translation, shear, and stretch. For each example, the first row of Figure 2 shows the positions of particles at the first (black) and second (magenta) time points. The second row shows the displacements for correctly tracked particles in green. The third and fourth rows show the displacements determined by the diffusion tracker (red) and the strain tracker (blue).

Figure 3: Histograms of individual particle displacements (a-d) and strain eigenvalues (e-h) from analysis of simulated data in Figure 2. (a) and (b) are histograms of horizontal and vertical displacements, respectively, from simulated diffusion data. (c) and (d) are histograms of horizontal and vertical displacements, respectively, from simulated translation data. (e) and (f) are histograms of strain eignevalues from a simulated shear. (g) and (h) are histograms of strain eigenvalues from a simulated stretch. In all panels, results from the diffusion tracker are in red, those of the strain tracker are in blue, and the overlap is in purple. In green is an outline of the histogram of correct displacements in (a),(b) and a vertical line corresponding to the correct values of displacements and strains in (c)-(h).

Not surprisingly, the diffusion tracker out-performs the strain tracker for the case of diffusion. For the example shown in Figure 2, the diffusion tracker returns results for all of the particles. About 96% of these tracks are correct. On the other hand, the strain tracker provides identifications for only 84% of the particles, and about 15% of these are incorrect. These tracking errors impact quantitative measures of the particle motion. Even though the distributions of the particle displacements for the two trackers look similar, Figure 3, tracking errors tend to introduce counts in the tails of the displacement distribution that have a significant impact on the mean-squared particle displacement (MSD). In this example, the MSD calculated from the diffusion tracker is within 2% of the correct value. However, the value calculated using the strain tracker is 30% greater than the correct value.

On the other hand, for the cases of large displacements due to translation, shear, and stretch, the strain tracker is a much more reliable choice. For the examples in Figure 2, the diffusion tracker returns tracks for greater than 84% of the particles, but less than 8% of them are correct. On the other hand, the strain tracker identifies greater than 73% of the particles and less than 3% of them are incorrect. Furthermore, for the case of the strain tracker, dropped particles and errors tend to be localized to the boundaries of the field of view, which are easily discarded when greater accuracy is required.

In these cases, the strain tracker accurately quanitfies particle displacements and strains. The expected displacements and strains for the cases of translation, shear, and stretch are shown as vertical green lines in Figure 3. The red histograms display the values for each particle returned by the diffusion tracker and the blue histograms report those from the strain tracker. The strain tracker consistently returns the correct values, while the incorrect tracks from the diffusion tracker return broadly distributed incorrect values.

Comparison of tracking strategies with experimental data

In this section, we compare the diffusion and strain tracker with some experimental data.

Figure 4: Particle tracking examples with experimental data.
(a)-(b) Trajectories of particles embedded in a silicone gel undergoing uniform compression, calculated with the diffusion tracker (a) and the strain tracker (b). Arrows are scaled by a factor of 2. (c)-(f) Histograms of strain eigenvalues as calculated with the diffusion tracker (c), (e) and the strain tracker (d),(f). (g)-(h) A contractile fibroblast cell adherent on a silicone gel. Arrows overlaid on top represent particle displacements from cell traction forces as calculated with the diffusion tracker (g), and the strain tracker (h). Arrows are scaled by a factor of 5. Scale bars are

First, we assess both trackers’ ability to identify particles during a relatively large and homogeneous strain. Here, we deform a silicone gel with embedded fluorescent tracers using a macroscopic deformation. We image the tracers over a small field of view, where we expect a reasonably uniform compression. We manually applied a uniform affine transformation to the particle locations to estimate the strain, resulting in eigenvalues of -0.06 and -0.04. Then,we calculated particle displacements with the strain and diffusion trackers, as shown in Figure 4(a,b). While both trackers return comparable results where the displacements are small (upper-left corner of the field of view), they disagree where the displacements are large. In this region, the diffusion tracker returns an uncorrelated displacement field, while the strain tracker returns the expected compression. Our strain tracker returns a strain tensor associated with the displacement of each tracked particle. Histograms of the strain eigenvalues of each particle are plotted in Figure 4 c-f. For both trackers, the peaks of the histogram agree well with the expected values from the manual affine transformation. However, the diffusion tracker reports a much broader distribution, including some unphysical positive values.

Next, we consider an example from Traction Force Microscopy (TFM) where the strains have a much stronger spatial heterogeneity. In TFM, forces exerted by small objects, such as cells, are quantified by measuring the deformation of an elastic material they are adhered to. To quantify the deformation, fluorescent particles are embedded in the elastic material. Displacements caused by a fibroblast cell adherent to a silicone gel are presented in Figure 4 g,h. The cell is fluorescently tagged and is displayed in inverted contrast. Overlaid on top of the cell image are displacements of the surface underneath the cell (scaled by a factor of 5), determined by the (g) the diffusion tracker and (h) strain tracker. While displacements measured by the strain tracker and diffusion tracker mostly agree, the diffusion tracker identifies some unexpected large strains near the center of the cell.


We have introduced a particle tracking algorithm based on the minimization of energy. This approach out-performs the conventional squared-displacement minimizing particle tracker when the displacements are larger than the interparticle spacing. On the other hand, our strain tracker may have difficulty when the strain changes significantly over the typical interparticle spacing. This can occur in the vicinity of strain singularities, such as cracks, and when particles are moving randomly, as in Brownian motion. Although, the code included in the supplement is designed for two time points and two spatial dimensions, expanding to three dimensions and many time points is a straightforward modification.


We thank Madhusudhan Venkadesan, Katharine Jensen, and Robert Style for helpful discussions and support from the Army Research Office Multidisciplinary University Research Initiative (ARO MURI) W911NF-14-1-0403.


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