{'fraction_non_alphanumeric': 0.03586747522555835, 'fraction_numerical': 0.017268155598284277, 'mean_word_length': 4.4057565460723564, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Precise Particle Tracking Against a Complicated Background: Polynomial Fitting with Gaussian Weight

We present a new particle tracking software algorithm designed to accurately track the motion of low-contrast particles against a background with large variations in light levels. The method is based on a polynomial fit of the intensity around each feature point, weighted by a Gaussian function of the distance from the centre, and is especially suitable for tracking endogeneous particles in the cell, imaged with bright field, phase contrast or fluorescence optical microscopy. Furthermore, the method can simultaneously track particles of all different sizes, and allows significant freedom in their shape. The algorithm is evaluated using the quantitative measures of accuracy and precision of previous authors, using simulated images at variable signal-to-noise ratios. To these we add a new test of the error due to a non-uniform background. Finally the tracking of particles in real cell images is demonstrated. The method is made freely available for non-commencial use as a software package with a graphical user-inferface, which can be run within the Matlab programming environment.


Introduction
Particle tracking-the use of computer analysis to measure the motion of individual particles in a sequence of images-is critical for a number of biophysical techniques, including particle tracking microrheology [1,2,3], magnetic [4] and optical tweezers [5]. It is also useful for studying the transport of fluorescentlylabelled single-molecules and organelles within the cell [6,7], and is widely used for particle imaging velocimetry [8,9].
A wide range of methods have been used to track particles in optical microscope images [10,11,12,13,14]. The most accurate for tracking small spherical particles (i.e. of diameter less than the wavelength of light) are based on the fitting of each intensity peak in the image with a two-dimensional Gaussian function. This method works well in images with a uniform background, but tends to fail if the background is complicated, or when the particles are larger or non-spherical [11]. For larger or non-spherical particles, tracking methods based on calculating the brightness-weighted centroid of a particle tend to give good results [11,12]. Here pixels are normally judged to be "inside" a particle if their brightness values are above a chosen threshold. However the latter method has an inherent drawback: that if a pixel on the edge of a particle is just above threshold in one frame, and just below threshold in the next, then a spurious jump in the centroid's position is recorded. This error becomes larger with increasing noise or decreasing diameter of the tracked particle. Methods based on calculating centroids also tend to fail if the image background is not uniform, since it becomes more difficult to exclude the background simply by subtraction or thresholding.
To date, the clearest analysis of errors in these tracking methods has been given by Cheezum et al. [11], who consider numerically simulated images of a single particle in a microscope image. The simulation is simple but realistic, and the authors find that the signal to noise ratio S/N is crucial in determining the tracking error. The error can be discussed in terms of accuracy and precision, which here respectively describe the determinate error-i.e. the bias B in the position with respect to the geometry of the system, and the indeterminate error-i.e. the standard deviation σ between measured and true positions. It is unfortunate that most scientists still rely on crude measures of their particle tracking accuracy or precision, such as measuring the apparent displacements of immobilised particles. This measure can only yield a lower bound to the error, as it ignores errors due to the displacement of the particle relative to the optical elements, CCD pixels, and visible features within the sample. Moreover, it is important to understand that tracking error is not analogous to the Rayleigh criterion, i.e. it is not a constant measured in nanometers for the particular microscope configuration. Instead it depends strongly on the signal intensity, that is, the contrast between the particle and its surroundings.
We present a new tracking method built on a different concept: a two-dimensional polynomial of 4th order is fitted to each particle image. The fit is weighted by a Gaussian function that is centred on the particle, with a decay length set to the particle's radius, which is measured from the image. The particle centre is itself determined by the polynomial's local maximum. This "polynomial-fit Gaussian weight" (PFGW) method is explained in detail below, but we can first see three obvious advantages of using a polynomial fit. Firstly, like Gaussian fitting, the fit is not sensitive to the edges of a particle. Secondly, like the centroid methods, the fit has significant flexibility towards the shape of a particle. Thirdly, unlike the previous methods, the procedure of determining the particle centre is inherently insensitive to the image background.
We evaluate the new fitting method using the analysis of Cheezum et al. [11], and describe a new test of the effect of a simple non-uniform background on particle tracking errors. Finally we demonstrate the new method on real cell images where its accuracy and precision can also be appreciated.
The PFGW method is made freely available for non-commencial use, in the form of a software package with a graphical user-inferface, named PolyParticleTracker. The package requires the Matlab programming environment, version 7.0 or higher.

Particle tracking method
The particle tracking method, as implemented in PolyParticleTracker consists of the following steps: 1. Image smoothing and noise reduction, 2. Identification of particles to track and estimation of particle coordinates, 3. Subpixel refinement of particle coordinates, 4. Particle discrimination, 5. Linking of particle postions between frames.
As with other tracking methods, it is the subpixel refinement algorithm (Step 3 above) that is the dominant factor which determines the accuracy and precision of the analysis [11]. This is the novel part in our method, so it is described in detail below. The other steps are similar to previous methods and are described in the Appendix for completeness.

Subpixel refinement of particle positions
The subpixel refinements are made as follows. For an identified particle with estimated position (x 0 , y 0 ), and estimated radius R, we find the closest fit (in the least-squares sense) of the intensity map of the image I(x, y), with a quartic polynomial function: where x and y denote pixel coordinates. The fit is given a weight at each pixel by the Gaussian function: where setting the decay constant α = 1 ensures that particle itself is weighted and the surroundings are excluded. The least-squares fit is sought in the domain x 0 − βR ≤ x ≤ x 0 + βR and y 0 − βR ≤ y ≤ y 0 + βR where setting the constant β = 2 ensures that the particle is adequately encompassed within the domain. We find that the subpixel correction itself is fairly insensitive to both α and β. Example quartic fits of particles in a real cell image are shown in figure 1.
Having fitted the particle image, we use the coefficients of the quadratic part of I f it (i.e. P 00 , P 10 , P 01 , P 20 , P 11 and P 02 ,) to evaluate the subpixel correction to the particle centre coordinates. The refined particle centre is defined by the extremum of the quadratic part: i.e. we take x 1 = x 0 +∆x and y 1 = y 0 +∆y, where ∆x = P 11 P 01 − 2P 02 P 10 4J , and J = P 20 P 02 − P 2 11 /4. (For a reference on quadratic surfaces, see [15].) However, we can make the particle centre coordinates more precise if the subpixel correction is made iteratively as follows. If either ∆x or ∆y are above a threshold fraction of a pixel Q, then the approximate centre is moved by Q, and the fitting procedure and calculation of (∆x, ∆y) are repeated: where ∆x, ∆y and the particle radius are recalculated in each iteration. Setting Q = 0.5 gives good precision, and the cycle is ended three iterations after the conditions ∆x, ∆y < Q are met. The cycle is also terminated if |x n − x 0 | > βR (i.e. the centre has wandered too far without finding the particle) or if J < 0 (i.e. the fit no longer resembles a particle with an extremum in intensity). Provided that images are acquired fast enough so that particles do no move further than their own radii, from one frame to the next, then the coordinates of a particle in frame T can be used as the initial estimate for its position in frame T + 1. This procedure provides a ready method of linking a particle's trajectory (Step 5 above). The coefficients P ij also provide the required information to calculate the particle's radius, eccentricity, angle of rotation and skewness (see Appendix).

Quantitative evaluation on simulated images
We first reproduce the error analysis of Cheezum et al. [11]. Second, we introduce a new test, to simulate the tracking of a particle against a simple non-uniform background.
The analysis of Cheezum et al. [11] is based on simulated images of a single fluorescent spherical particle, imaged in a typical microscope with a 100x oil-immersion lens. The microscope is equipped with a CCD camera, modelled by a square array of pixels, such that the final resolution of the image is 99 nm per pixel. In the model, the particle is first represented by a circular disk of radius R with uniform emission brightness on a high resolution grid of 9 nm per cell, such that 11 cells of this grid correspond to the width of a CCD pixel. The image is then broadened by the point-spread funtion: where a = 2πNA/λ, NA=1.3, and λ = 570 nm. The brightness measured by each pixel of the CCD array is the average brightness of the corresponding 11x11 array of high-resolution cells. The CCD image intensity is then scaled to a peak value S +10 and a background value 10, where S is the signal strength and intensity is measured in number of photoelectrons. CCD images display "shot noise" when the photons counted by each pixel is low enough for statistical fluctuations to appear in their number. This noise follows a Poisson distribution: if the mean photoelectron count is I, then the standard error is √ I. We used the Poisson noise function of Matlab (Natick, MA, USA) to add Poisson noise to the CCD images. The signal to noise ratio of the intensity peak is then given by S/N = S/ √ S + 10. Two particle radii are considered below: a "point particle" of 9 nm, and a larger particle of 1-µm diameter whose edge can be clearly resolved. A range of S/N values are considered; at each 1000 frames are generated, in which the particle is moved by shifting the high resolution image by a whole number of cells in each of the x and y directions. We simulate Brownian motion by shifting the image in each direction according to a Gaussian distribution of zero mean and standard deviation 27 nm (three high-resolution cells).
The tracks generated by our method are compared with the exact track as follows. The bias and standard deviation are measured as where the average is taken over all 1000 successive images. Here a represents each of the tracked particle coordinates x and y, andâ represents the corresponding exact coordinate.
The results for tracking the point particle with the PFGW method are as follows: figure 2 shows that both B and σ are almost uncorrelated with position relative to the pixel edges at both low and high signal strengths. Furthermore, B is approximately a factor of 10 smaller than σ at each S/N ; thus the position of a particle relative to the CCD pixels, as determined by the method, is effectively unbiased.
The mean precision for tracking both the point particle and the 1-µm particle is evaluated for a range of S/N and plotted in figure 3. Results for Gaussian-fit tracking of a point-particle are superposed for comparison (from Cheezum et al. [11]). PFGW obtains better precision for the point particle than the 1-µm particle, despite the fact that more pixels are included in the fit of the latter. This may be explained by considering that the brightness profile of the 1-µm particle is almost flat near the peak-thus changes in this region make a relatively large difference to the peak position of the polynomial fit. However, it is especially encouraging that PFGW is more precise than Gaussian-fitting for the point particle, and still works well for 1-µm particles. For the former, our method reaches 1-nm precision at S/N = 41.2.
In real images, especially images of cells in transmission light microscopy, there is usually a significant non-uniform background intensity around each particle. Even in fluorescence microscopy, there are normally more than one fluorescent particle, at various distances from each other, which can introduce new errors into their tracking.
We introduce a simple stylised background into the simulation: a series of "ridges"-bright thin parallel lines-are placed 200 cells apart (corresponding to 18.2 CCD pixels) in a high-resolution image of the same dimensions as above. This background image is scaled such that its maximum intensity is half that of the particle image, and the two images are added. The combined image is convolved with the PSF, before producing the low-resolution CCD image and generating the noise. In 1000 successive frames, the particle is shifted diagonally across the ridges, by 27 nm in both the x and y directions each time. The motion Figure 3: Tracking precision in pixels using the PFGW method, as a function of S/N . The best precision is obtained for the point-particle (+), where σ = 1 nm is reached at S/N = 41.2. 1-µm particles are tracked with slightly higher error (×). Previous results for Gaussian-fit tracking of a point (•) are plotted for comparison [11].
was tracked using the PFGW method (figure 4). We then applied the efficient and popular method of Crocker and Grier [12], which is based on a centroid calculation with background subtraction, to the same images. In their method, the background is approximated by a "boxcar average"-i.e. a convolution of the image with a filled square of width 2w + 1 pixels. The background is then subtracted from the image, replacing with zero any intensity that falls below zero. Choosing w is a compromise between two extremes: w must be larger than the particle diameter, otherwise the particle itself will be subtracted, but it must not be so large that it significantly influences a particle by subtracting from it other nearby objects. Their subtraction step is particularly useful because it yields a region of zero intensity surrounding each particle, provided that the background is smoothly varying. However, the background we consider is outside this provision. Figure 5 shows the accuracy and precision of the point particle, as a function of its distance from one ridge to the next. The PFGW method is compared to the Crocker-Grier method, taking w = 4 and w = 10 in the latter. (w = 4 just fulfils the requirement of being greater than the PSF diameter, while w = 10 is near the maximum particle-ridge separation.) The PFGW results are biased only where the images of the particle and ridge are overlapping (∼3 pixels either side of the ridge). The Crocker-Grier method results are biased within ∼ w pixels of the ridge, i.e. where the boxcar function overlaps the ridge. The precision of the PFGW method is unaffected by the background (0.02 pixels throughout), but the Crocker-Grier method is much less precise where the particle is within w pixels of the ridge-this is due to the noise of the background being subtracted from the particle image by the boxcar average. We note that the problems with background subtraction would be increased if the background were allowed to move, since the bias due to the background would be continually changing. This would be the case, for example, where two moving fluorescent particles approach each other.

Tracking of real particles
We demonstrate our method firstly on bright field images and secondly on images of fluorescent particles. The bright field images were taken in an Olympus IX71 inverted microscope (Olympus UK Ltd., London EC1Y 0TX, UK) with a 100x oil-immersion lens and a 1.6x magnifier, using a Photron FastCam PCI Figure 4: A simple test of particle tracking against a background: a point-particle (centre) diagonally traverses a series of thin ridges. The exact trajectory (· · ·) and PFGW-tracked course (-) are superposed. Figure 5: (a) Bias B and (b) standard deviation σ for a point-particle tracked against a background of sharp ridges, spaced periodically at 18.2 CCD pixels, where each ridge has a peak brightness of half the particle brightness, and S/N = 31.3. Results with the PFGW method (·) are compared with the method of Crocker and Grier [12] where the background has been subtracted using a boxcar average of width w = 4 (+) and w = 10 (•). For the results of PFGW, B is negligible when the particle and ridge images do not overlap, while σ is unaffected by the background, even where the particle overlaps the ridge. For the Crocker-Grier method, bias and σ are both large if the particle-ridge separation is less than w.
CCD camera (Photron (Europe) Ltd., Bucks SL7 1NX, UK). The fluorescent images were taken with a Leica SP2 AOBS confocal microscope (Leica Microsystems GmbH., 35578 Wetzlar, Germany). For live cell images, cells were incubated in DMEM supplemented with 10% fetal bovine serum and 15 µm HEPES, in chambers with glass coverslip bases, incubated at 37 • C. 0.5-µm-diameter latex beads (Sigma-Aldrich, Gillingham UK, cat. no. L2153) were immobilised on a glass coverslip by drying, and images captured at frame rates of 10,000 and 45,000 fps. Both sets of images were captured without adjusting the illumination or focus, then particles were tracked using the PFGW method. Figure 6(a) shows a sample frame, with the bead centre coordinates marked. Figure  6(b) shows the mean-square displacement function (MSD), which shows a static error at short timescales of √ MSD = 1.9 and 8.3 nm for 10,000 and 45,000 fps respectively [16]. (At larger timescales a 150-Hz vibration from the building is apparent.) The reason that the static error is higher for the faster image, is that the amount of light per frame has been reduced by a factor of 4.5, and hence the signal to noise ratio has also been reduced. However, in this case, the precision is not properly described by the shot-noise model described above for two reasons: the sensitivity of the CCD pixels is far below single-photon levels, and the pixel brightness is discretized with an 8-bit colormap. For the images at 10,000 fps, we measure a mean peak pixel brightness I peak = 143 and standard deviation ( (I peak − I peak ) 2 ) 1/2 = 1.15, where the averages are taken over all time frames. These values compare with a mean brightness of the whole image I = 82.3 and standard deviation ( (I − I ) 2 ) 1/2 = 0.92. Clearly the fluctuation amplitudes are approximately constant over the image, unlike shot-noise, with an amplitude of one brightness unit. Thompson et al. [17] have given theoretical expressions for the precision of Gaussian fitting in the two limits of shot noise and constant background noise. The latter limit (taking background noise b = 1) predicts the present error remarkably well. By substituting our peak signal intensity of S in their theory, we obtain a modified form of their equation 13 [17]: where s = 1.1 is the standard deviation of the PSF in pixels. Taking the signal to background noise ratio of the 10,000 fps data as S/b = (143 − 82.3)/1 gives σ = 2.0 nm. For the 45,000 fps data, equation 7 then predicts the standard deviation 4.5 times higher at 9.1 nm. Both values of σ are close to the measured static errors, thus it appears that the precision is limited, in this case, by discretization of pixel brightness. The PFGW method can be used to track low-contrast endogeneous particles in cells. Figure 7(a) shows an example of particle tracks in a well-spread L929 cell. Bright field images were captured at 1000 fps, with a final resolution of 106 nm/pixel. Some particles display only random motion, but many display active motion along tracks of various shapes (examples highlighted), some of which abruptly change direction. Although we cannot distinguish the particles or motors involved, the tracking method is of sufficient precision to see the discrete steps of many particles along the fibres of the cytoskeleton. An example run of 36-nm steps is shown in figure 7(b); in this case, the particle position is plotted as an arc length along the fibre, calculated by taking the particle's nearest position along a "smoothed" track. The smoothed track was produced by convolving the particle's x and y coordinates with a Gaussian function of large width. The MSD of each particle present for more than 10,000 frames is plotted in figure 7(c), and shows two different scaling regimes similar to those identified by previous authors [18]. Most particles show MSD∼ t 1/2 at short times t 0.1 s, indicative of tethered thermal motion on stiff filaments [19], and the exponent ≈ 2 at long times t 0.1 s, indicative of active motion. The fact that the MSD has not reached a plateau at the fastest timescale of 1 ms shows that the measured tracks are dominated by the real motion of the particles; thus the precision is better than MSD(1 ms) = 4 nm.
Another example is particle motion in streaming pseudopodia of Amoeba proteus, as shown in figure  8. The presence of a large number of endogeneous particles allows its internal motion to be measured at Finally we used fluorescence microscopy to track motion in HeLa cells of the same latex beads mentioned above, which fluoresced in green light. The beads were endocytosed by the cells: figure 9 (inset) shows a portion of a cell with bead tracks superposed. The brightness fluctuations at each particle peak were measured and show a reasonable correspondence to the shot noise model in figure 9; each peak brightness S is plotted against standard deviation of fluctuations N . The PFGW method successfully tracked particles with S/N 3.
We are currently applying the PFGW method to the microrheology and particle imaging velocimetry of cortical oscillations and amoeboid migration (article in preparation).

Conclusion
Our particle tracking method, based on polynomial fitting of feature points with a Gaussian weighting function, is distinguished from previous methods of particle tracking by its inherant ability to track particles against a complex background, while maintaining a similar precision to Gaussian fitting at a range of signal to noise ratios. It is ideal for particle tracking in situations where a complicated background or neighbouring particles produce errors in other methods, e.g.: intracellular microrheology and tracking in dense colloidal situations.

Instructions for obtaining PolyParticleTracker
The package PolyParticleTracker can be downloaded from our website: http://personalpages.manchester.ac.uk/staff/salman.rogers/polyparticletracker/. The package includes a graphical user interface, and requires Matlab version 7.0 or higher.

A Details of PolyParticleTracker
A.1 Image smoothing and noise reduction Following Crocker and Grier [12], we smooth each image to correct for discretization noise, by convolving the image with a Gaussian function: where the correlation length λ is set to 1, and setting w = 3λ sufficiently approximates the unbounded convolution. Note the background is not subtracted as in previous particle tracking methods.
A.2 Identification of particles to track and estimation of particle coordinates Particles are identified by the presence of either a local maximum or minimum in I(x, y). For each of these extrema (x 0 , y 0 ), the adjacent points of inflexion are found along the orthogonal lines defined by (x, y 0 ) and (x 0 , y). The mean of the distances to the four points of inflexion is used as the initial radius estimate R 0 . x 0 , y 0 , R 0 are then input into the subpixel correction routine, described above.

A.3 Particle discrimination
For each particle, the intensity map I(x, y) and the parameters of the polynomial fit P ij allow calculation of the particle's eccentricity, rotation, radius, average brightness and skewness, as follows. Eccentricity e and rotation θ, i.e. the angle between the major axis of the particle's brightness peak and the x axis, are both obtained by considering the quadratic terms which dominate the vicinity of the peak. These terms define a surface whose horizontal section is an ellipse: i.e. i+j=2 i=0,j=0 P ij x i y j = 0. Hence [20] θ = 1 2 cot −1 P 20 − P 02 P 11 , e = 1 + (P 20 − P 20 ) 2 + P 2 11 + (P 20 + P 20 ) (P 20 − P 20 ) 2 + P 2 11 − (P 20 + P 20 ) .
Calculation of the particle's radius R is more complicated, as the true shape of the particle is not known, and the background intensity is non-negligible. However, a good estimate of R is given by the geometric mean of the distances from the peak of I f it (x, y) to the four points of inflexion in the directions of the particle's major and minor axes, i.e. the directions θ, θ + π/2, θ + π, θ + 3π/2. The radius calculation is illustrated in figure 10. In terms of the polynomial coefficients: where A =P 20 cos 2 θ − P 11 cos θ sin θ + P 02 sin 2 θ , B =P 20 sin 2 θ − P 11 cos θ sin θ + P 02 cos 2 θ , C =P 40 cos 4 θ − P 31 cos 3 θ sin θ + P 22 cos 2 θ sin 2 θ − P 13 cos θ sin 3 θ + P 04 sin 4 θ , D =P 40 sin 4 θ + P 31 sin 3 θ cos θ + P 22 sin 2 θ cos 2 θ + P 31 sin θ cos 3 θ + P 04 cos 4 θ . is the transparent mesh. The peak of I f it (x, y) is marked by a triangle and its four points of inflexion in the directions θ, θ + π/2, θ + π, θ + 3π/2 are marked by circles. The estimated radius is the mean distance from the peak to the points of inflexion.
The cubic terms in I f it (x, y) determine how skewed or lopsided the particle peak is. We can define an approximate skewness of a particle by analogy to the skewness of a statistical distribution, as the ratio of the cubic terms to the quadric terms at x, y ≈ R:
The eccentricity, brightness, radius and skewness can all be used to discriminate whether a particle is acceptable for tracking or not.
A.4 Linking of particle postions between frames As mentioned above, the coordinates of a particle in frame T are used as the initial estimate for its position in frame T + 1, providing a reliable method of linking a particle's trajectory, provided the particle does not move further than its own radius, from T to T + 1. If, however, the particle has moved further, PolyParticleTracker will look for a nearby extremum in I(x, y), which is nearer than the next-nearest identified particle and has a similar brightness to the original particle in frame T . PolyParticleTracker can also be set to search for new particles every N frames.