PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
Med Image Comput Comput Assist Interv. Author manuscript; available in PMC 2010 June 21.
Published in final edited form as:
Med Image Comput Comput Assist Interv. 2009; 12(Pt 2): 673–681.
PMCID: PMC2888536
NIHMSID: NIHMS128408

Actin Filament Tracking Based on Particle Filters and Stretching Open Active Contour Models

Abstract

We introduce a novel algorithm for actin filament tracking and elongation measurement. Particle Filters (PF) and Stretching Open Active Contours (SOAC) work cooperatively to simplify the modeling of PF in a one-dimensional state space while naturally integrating filament body constraints to tip estimation. Existing microtubule (MT) tracking methods track either MT tips or entire bodies in high-dimensional state spaces. In contrast, our algorithm reduces the PF state spaces to one-dimensional spaces by tracking filament bodies using SOAC and probabilistically estimating tip locations along the curve length of SOACs. Experimental evaluation on TIRFM image sequences with very low SNRs demonstrates the accuracy and robustness of the proposed approach.

1 Introduction

Actin proteins are present in all eukaryotic cells. Their ability to polymerize into long filaments (over 10 µm in length and ~7 nm in thickness) underlies basic processes of cell life such as cell motility, cytokinesis during cell division, and endocytosis. The kinetics of polymerization of individual actin filaments in vitro have been studied extensively using Total Internal Reflection Fluorescence Microscopy (TIRFM) [1], [2], [3] (Figure 1(a) and 1(c)). In these experiments, actin filaments are attached to a glass slide by surface tethers that act as pivot points that can be used as fiducial markers to help distinguish the elongation of each end [1], [2], [3]. The two ends of an actin filament, the “barbed” (B-end) and “pointed” (P-end) end, respectively, grow at distinctively different rates. Two basic features of actin kinetics that can be extracted from TIRFM are (i) the average rate of filament elongation at each end, and (ii) the fluctuations in the average rate. Both of these two numbers depend in a unique way on the details of the microscopic mechanism of monomer addition to the ends of the filament [1], [2], [4].

Fig. 1
An example of actin filament tracking using the proposed method.

In [5], we presented the stretching open active contours (SOAC) for filament segmentation and tracking. This automated method allowed simultaneous measurements of multiple filaments, thus enabling the extraction of statistics on filament elongation. In earlier studies, such measurements were typically performed using manual or semi-automated methods [3]. Related methods have been applied to tracking microtubule (MT) filaments. Hadjidemetriou et al. [6] minimized an image-based energy function to segment MTs at each time step using consecutive level sets method. Saban et al. [7] automatically detected tips in the first frame and then proceeded to track each tip separately by searching for the closest match in subsequent frames. However, the above automated methods did not utilize temporal coherency of the motion or the growth of filaments.

Particle Filter (PF) based methods have also been proposed to track MTs. In [8], PFs were utilized to track the tips of polymerizing microtubules. In these images, MTs were labeled by plus-end tracking proteins that associate with the tips but not with the body of growing MTs. Without using supporting information from the MT body, Markov Random Fields were employed to model the joint transition model of multiple tips; this required that the posterior pdfs of different tips to be sampled jointly. A limiting factor in this work was that its high-dimensional state space added complexity to modeling and tracking computation. Kong et al. [9] employed a PF-based method similar to [8] to track MT tips. In this work, MT tip locations were estimated recursively using PF, and then MT bodies were segmented using active contours based on those estimations.

Although in [5] filament bodies were tracked accurately, we reported errors on tip location estimation because of the low SNR of filament tips. In this paper, we present a novel actin filament tracking algorithm which combines particle filters (PF) [10] with the stretching open active contours (SOAC) to address this problem. An example of our tracking results is shown in Figure 1. By construction, SOACs stretch along bright ridges in an image. At each time step t, SOACs stretch along filament bodies. But unlike SOACs in [5], they are forced to grow over distances that exceed the tip locations. Subsequently, particles are spread along each SOAC according to transition models, which predict the tip “length” (location) at time t + 1, given the states and likelihoods of the particles at time t. Using the image information at time t + 1, each particle is associated with a likelihood that the particle represents the correct tip. The particles’ states and likelihoods are then used to estimate the posterior pdf describing the probability distribution of tip length at time t + 1. Before tracking, this PF-based method requires the transition and likelihood models to be properly defined by providing an initial estimate of the value of the filament elongation rate, to within ~20% of the actual value.

The main contribution of this paper is that we reduce the state space of PF to a one dimensional space by implicitly modeling tip “length” as the state vector of our PF. By construction, tip locations are distributed at SOACs that grow along filament bodies; thus we only need to track the “length” of tips. This novel framework naturally integrates filament body constraints to tip estimation. This is a computationally more efficient method as compared with the 4-dimensional space in [9] and the even higher dimensional space in [8], and different from other filament tracking algorithms which either track tips [7], [8], [9] or entire bodies [6] of MTs. Experimental evaluation and comparison showed accurate and robust tracking performance of our algorithm.

2 Stretching Open Active Contour Models

In [5], open active contour models were presented to segment actin filaments. Let r(s) = (x(s), y(s)), s [set membership] [0, 1] parametrically represent an open curve. Treating this curve as an active contour, we seek to minimize its overall energy E, which consists of internal energy Eint and external energy Eext, i.e., E = Eint + k · Eext, where k controls the relative contributions of the internal and external energy. The internal energy, Eint=01(α(s)|rs(s)|2+β(s)|rss(s)|2)ds, controls the smoothness of the curve. The external energy, Eext, represents external image constraints and consists of an image term Eimg and an intensity-adaptive stretching term Estr, i.e., Eext = Eimg + kstr · Estr, where kstr balances the contributions of Eimg and Estr. The image term Eimg is defined by the Gaussian filtered image, i.e., Eimg = Gσ * I, to make the active contour converge to filament locations, which correspond to bright ridges in the images. The stretching term stretches the open active contour along a filament and stops stretching when its ends meet filament tips. Tip locations are determined by using an intensity-based criterion.

The main problem of [5] is that, in some very noisy TIRFM image sequences, the intensity-adaptive stretching term Estr may lead to large errors on tip estimation. Low contrast near filament tips may result in active contours “over-growing” (see Figure 2(a)), while intensity gaps on filament bodies may lead to active contours “undergrowing” (see Figure 2(b)). It was difficult to define an external energy term that copes with both scenarios for all filaments.

Fig. 2
Problems observed in our previous method [5]. (a) An “over-grown” active contour because of the low SNR near the tip, and (b) an “under-grown” active contour because of the intensity gap on filament body. (c) An over-grown ...

Although tips were difficult to be identified by intensity or contrast alone in noisy images, we observed that the over-grown active contours followed filament bodies accurately and were able to cover tip locations even when they over grew (See Figure 2(c)). This means that we can search for tips along an over-grown active contour. Therefore, we propose a new stretching open active contour (SOAC) model similar to the one proposed in [5] but with a non-intensity-adaptive stretching term Estr, which always makes a SOAC grow over distances that exceed the tip locations:

Estr(r(s))={ rs(s)|rs(s)|s=0,rs(s)|rs(s)|s=1,0otherwise.
(1)

Therefore, the overall energy of a SOAC model is defined by

E=Eint+k01Eimg(r(s))+kstrEstr(r(s))ds,
(2)

where kstr is a constant that controls the stretching force that makes the active contour always grow over tips. The above SOAC model enables us to reduce the search space of filament tips to a one-dimensional space (along the SOAC’s curve length) and naturally adds a continuous body constraint: a tip must be connected to a filament body.

To address the filament intersection problem, we propose a strategy different from [5] since the new SOAC models do not identify intersections during their deformation. When computing each SOAC’s external energy, we simply set pixels covered by other SOACs to have background intensity. This new strategy has been demonstrated to be more effective and resulted in fewer tracking failures (Section 5.2).

3 Actin Filament Tracking using Particle Filters

A SOAC model provides an estimation of filament body and therefore simplifies the problem of tip tracking to searching for and tracking tip patterns in a one-dimensional space along the SOAC curve’s length. This is one of the major advantages when compared with approaches that work in a 4-dimensional space in [9] and the even higher dimensional space in [8]. To systematically search along an over-grown SOAC for the optimal tip location, we employ the widely used Sequential Importance Resampling (SIR) PF [11] to estimate the locations of both B-end and P-end of a filament.

3.1 Bayesian Tracking based on Particle Filters

In the PF framework, the state vector of a target at time t, Xt [set membership] Rn, is given by Xt = ft(Xt−1, vt−1), where vt−1 [set membership] Rn is an i.i.d. process noise vector, and ft : Rn × RnRnis a possibly nonlinear function that is modeled by a known transition model Pt (Xt|Xt−1). The measurements at time t, Zt, relate to the state vector by Zt = ht(Xt, nt), where nt is an i.i.d. measurement noise vector, and ht : Rn × RnRn is a possibly nonlinear function that is modeled by a known likelihood model Pl(Zt|Xt). It is assumed that the initial posteriori pdf P(X0|Z0) [equivalent] P(X0) is available as a priori.

The objective of tracking is then to recursively estimate Xt from measurements. Suppose that P(Xt−1|Z1:t−1) at time t − 1 is available. In principle, the pdf P(Xt|Z1:t) can be obtained, recursively, in two stages: prediction and update.

Prediction:  P(Xt|Z1:t1)=Pt(Xt|Xt1)P(Xt1|Z1:t1)dXt1,
(3)

Update:  P(Xt|Z1:t)Pl(Zt|Xt)P(Xt|Z1:t1).
(4)

In the SIR PF algorithm, at each time t, a group of N particles {Xt1(i),wt1(i)}i=1N is propagated from t − 1 to characterize P(Xt|Z1:t) using the Monte Carlo principle: P(Xt|Z1:t)i=1Nwt(i)δ(XtXt(i)). It is usually assumed that P(Xt|Zt) = P(Xt|Z1:t).

Clearly, three aspects of the PF framework need to be specified:

  • The state vector Xt, which models the system in a n-dimensional space;
  • The transition model Pt (Xt|Xt−1), which models the transition function ft ;
  • The likelihood model Pl (Zt|Xt), which models the measurement function ht.

For the state vector in our filament tip tracking problem, as we have simplified tip location estimation to a 1D space using SOAC, we propose to use the “length” of the tip at time t, St, as the state of our PF, i.e., Xt = S t. The “length” of a tip is defined as the filament length from the tip point to a reference point along the SOAC representing the filament. Reference points can be randomly chosen on SOACs during initialization. Tracking the length of a tip is equivalent to tracking its location on a SOAC.

For the transition model, because we choose tip length as the state vector, the transition model should describe the change of one tip’s length over time. Therefore, it can also be interpreted as the tip elongation model. As mentioned in Section 1, the B-end and P-end of a filament grow at different rates. We model the transition probabilities of B-end and P-end tips separately. For the transition model of B-ends, we used a normal density, PBt(Xt|Xt1)~𝒩(Xt1+μb,σb2), where µ b is the average elongation rate of B-ends. Obviously, the more accurate the estimation of µ b is, the more robust the tracking results are. For P-ends, because they grow much slower than B-ends, we set PPt(Xt|Xt1)~𝒩(Xt1,σp2) .

When a SOAC has just been initialized, the B-end and the P-end of the filament it represents cannot be distinguished. Therefore in the first few frames, we dispatch half of the particles following the B-end transition model and the other half following the P-end model. After several tracking steps that take into account image measurement information, the B-end and P-end of the filament can be distinguished by comparing the posterior pdfs estimated using B-end and P-end particles respectively.

For the likelihood model, we use an appearance template based approach. In particular, we use a 10 × 4 pixel rectangle template containing a mean appearance image µ T and a standard deviation image σT, both computed from manually selected tip image patches. When a filament A is intersecting another filament B, the tip of A would be occluded by B’s body. This naturally implies that any pixel covered by filament B should have low confidence or certainty when being used to compute the likelihood of the tip of A. Therefore, at each time t, a confidence map Mt for each filament is created, in which pixels covered by other filaments are given the value 0.5, and all other pixels are given the value 1. The likelihood model is then defined by

Pl(Zt|Xt)exp{1ni=1nMt(si)2σT(si)|F(Xt)(si)μT(si)|},
(5)

where n is the total number of pixels in the template tip patch, s i is its ith pixel’s coordinates, and F(Xt) and Mt are the image patch and the confidence map patch of a tip with state Xt, respectively, after being translated and rotated to the template tip’s coordinate.

3.2 SOAC Registration and Length Measurement

To perform measurements on the length and elongation rate of a tip, a reference point on the corresponding filament needs to be specified. However, a TIRFM image sequence shows drift such as translation between contiguous frames [3]. To recover the same reference point on a filament over time, the converged SOACs in consecutive frames representing the filament are registered simultaneously using the Iterative Closest Points (ICP) algorithm [12]. Since the drift is small, the ICP has achieved accurate registration in all our experiments.

4 The Algorithm

We summarize the proposed filament tracking algorithm as follows. Let r t−1 represent the SOAC active contour tracking the filament at time t − 1. From the particle sample set {Xt1(i),wt1(i)}i=1N at time step t −1 characterizing the pdf P(Xt−1|Zt−1), we can construct a “new” sample-set {Xt(i),wt(i)}i=1N approximating P(Xt|Zt) at time t:

  • Initialize [r with macron]t using rt−1 (Figure 3(a)). Deform it by minimizing Eqn. (2) and make sure both ends of [r with macron]t grow over the tips of its corresponding filament (Figure 3(b)).
    Fig. 3
    Illustration of the algorithm. (a) A SOAC rt−1 (red) at time t−1, (b) initialize [r with macron]t (blue) using rt−1 and deform it by minimizing Eqn. (2), (c) before registration of rt−1 and [r with macron]t (Red ‘*’ ...
  • Register [r with macron]t to rt−1 using the ICP method and recover the reference point on [r with macron] t (Figure 3(c–d)).
  • Select a sample X¯t(i)=Xt1(j) with probability wt1(j).
  • Predict Xt(i) to approximate P(Xt|Zt−1) according to Eqn. (3) by sampling from
    PBt(Xt|Xt1=X¯t(i))   or   PPt(Xt|Xt1=X¯t(i)),
    (6)
    depending on the tip type (B-end or P-end) the particle X¯t(i) represents (Figure 3(e)).
  • Measure and weight the new particle using the measurement Zt according to Eqn. (4) and (5). Then normalize all N particles to estimate P(Xt|Zt):
    wt(i)=Pl(Zt|Xt=Xt(i)),   then   wt(i)=wt(i)j=1Nwt(j).
    (7)
  • Estimate the mean of P(Xt|Zt) at time t by
    (Xt)i=1Nwt(i)Xt(i).
    (8)
  • Cut [r with macron]t according to x2130(Xt) to generate rt (Figure 3(f)).
  • Go back to the initialization step for t + 1.

5 Application To Experimental Data

5.1 Experimental Image Data

We used two TIRFM image sequences from [2]. In these experiments, polymerization of muscle Mg-ADP-actin was monitored in the presence of varying concentrations of inorganic phosphate (Pi) and actin monomers. 30% of the actin was labeled on lysine side chains with Alexa green. The pixel size was 0.17 µm. There were 20 images in sequence I and and 34 images in sequence II. The time interval between frames was 30 sec in sequence I and 10 sec in sequence II.

5.2 Evaluation and Comparison with The Previous Method [5]

For both sequences, we set α = 0.05, β = 0.1, k = 0.6, and kstr=0.55. µb for sequences I and II were set according to average elongation rates of B-ends measured by a manual method [3]. We set large values to σp and σb to avoid imposing any strong prior on tip length estimation. For sequence I, µb = 11.2524 mon/sec, σb = 4.5 mon/sec, and σp = 1 mon/sec. For sequence II, µb = 11.5760 mon/sec, σb = 5 mon/sec, and σp = 2.5 mon/sec. We used 50 particles for each tip and initialized each SOAC using the segmentation method proposed in [5] to obtain P(X0).

Figure 4 illustrates 3 examples of our tracking algorithm. Taking advantage of temporal coherence and filament body constraints, our algorithm tracked filaments accurately and showed robust performance against filament intersection.

Fig. 4
Three examples of tracking filaments. (1–2) Tracking filaments in sequence I, and (3) tracking a filament in sequence II.

We observed that filament bodies were always tracked accurately by our algorithm. Therefore, we evaluated the algorithm by measuring errors on tip location estimation. We selected 10 and 5 actin filaments from image sequence I and II respectively to measure tracking errors. For all selected filaments, we manually labeled their two tips in each frame as ground truth and calculated L2 distances between the ground truth and our algorithm’s results. We also compared with tip location errors obtained by a previous method [5], which did not utilize temporal coherence and used an intensity criterion to determine tip locations. During the tracking process, when we observed a failure, which means a SOAC stretched onto a different filament, we reinitialized the SOAC in the next frame by hand and resumed tracking. Table 1 shows tip tracking error statistics of our algorithm and of the previous method; our new one-dimensional PF-based algorithm clearly outperforms the previous method. The tracking errors for sequence II were higher because of its very low SNR.

Table 1
Tip tracking error statistics of selected filaments in both image sequences I (Figure 4(1–2)) and II (Figure 4(3)). (Unit: pixel)

6 Conclusion

In this paper, we introduce a novel actin filament tracking algorithm based on Stretching Open Active Contour (SOAC) models and one-dimensional Particle Filters (PF). Taking advantage of filament body estimated by the SOAC models, our method is able to find the tip using PF with a one-dimensional state vector. Such simplification naturally integrates filament body constraints to guarantee continuity between the estimated tip and body. A template based likelihood model and the stochastic nature of the PF framework also make our algorithm robust to noise and filament intersections. Experimental evaluation on TIRFM image sequences with low SNRs and comparison with a previous method demonstrate the accuracy and robustness of the proposed approach.

Acknowledgments

We would like to thank Ikuko Fujiwara (NHLBI/NIH) and Thomas Pollard (Yale) for providing the data, and Matthew Smith (Lehigh) for his suggestions. This work was supported by NIH grant R21GM083928.

References

1. Fujiwara I, Takahashi S, Tadakuma H, Funatsu T, Ishiwata S. Microscopic analysis of polymerization dynamics with individual actin filaments. Nat. Cell Biol. 2002;4:666–673. [PubMed]
2. Fujiwara I, Vavylonis D, Pollard TD. Polymerization kinetics of ADP- and ADP-Pi-actin determined by fluorescence microscopy. Proc. Natl. Acad. Sci. USA. 2007;104:8827–8832. [PubMed]
3. Kuhn JR, Pollard TD. Real-time measurements of actin filament polymerization by total internal reflection fluorescence microscopy. Biophys. J. 2005;88:1387–1402. [PubMed]
4. Vavylonis D, Yang Q, Pollard TD. Actin polymerization kinetics, cap structure, and fluctuations. Proc. Natl. Acad. Sci. USA. 2005;102:8543–8548. [PubMed]
5. Li H, Shen T, Smith M, Fujiwara I, Vavylonis D, Huang X. Automated actin filament segmentation, tracking and tip elongation measurements based on open active contour models. ISBI. 2009 [PMC free article] [PubMed]
6. Hadjidemetriou S, Toomre D, Duncan J. Motion tracking of the outer tips of micro-tubules. Medical Image Analysis. 2008;12:689–702. [PubMed]
7. Saban M, Altinok A, Peck A, Kenney C, Feinstein S, Wilson L, Rose K, Manjunath B. Automated tracking and modeling of microtubule dynamics. ISBI. 2006;1:1032–1035.
8. Smal I, Draegestein K, Galjart N, Niessen W, Meijering E. Particle filtering for multiple object tracking in dynamic fluorescence microscopy images: Application to microtubule growth analysis. IEEE Trans. on Medical Imaging. 2008;27:789–804. [PubMed]
9. Kong K, Marcus A, Giannakakou P, Wang M. Using particle filter to track and model microtubule dynamics. ICIP. 2007;5:517–520.
10. Gordon N, Salmond D, Smith A. Novel approach to nonlinear/nongaussian Bayesian state estimation. IEE Proceedings-F (Radar and Signal Processing) 1993;140:107–113.
11. Isard M, Blake A. Condensation–conditional density propagation for visual tracking. IJCV. 1998;29:5–28.
12. Besl P, McKay H. A method for registration of 3-D shapes. TPAMI. 1992;14:239–256.