Model based dynamics analysis in live cell microtubule images
© Altınok et al; licensee BioMed Central Ltd. 2007
Published: 10 July 2007
Skip to main content
Volume 8 Supplement 1
© Altınok et al; licensee BioMed Central Ltd. 2007
Published: 10 July 2007
The dynamic growing and shortening behaviors of microtubules are central to the fundamental roles played by microtubules in essentially all eukaryotic cells. Traditionally, microtubule behavior is quantified by manually tracking individual microtubules in time-lapse images under various experimental conditions. Manual analysis is laborious, approximate, and often offers limited analytical capability in extracting potentially valuable information from the data.
In this work, we present computer vision and machine-learning based methods for extracting novel dynamics information from time-lapse images. Using actual microtubule data, we estimate statistical models of microtubule behavior that are highly effective in identifying common and distinct characteristics of microtubule dynamic behavior.
Computational methods provide powerful analytical capabilities in addition to traditional analysis methods for studying microtubule dynamic behavior. Novel capabilities, such as building and querying microtubule image databases, are introduced to quantify and analyze microtubule dynamic behavior.
Microtubules (MTs) are filamentous cytoskeletal structures composed of tubulin protein subunits. These subunits can associate with, or dissociate from, existing tubulin polymers rapidly, making MTs highly dynamic. Through these dynamic behaviors, MTs are critically involved in many essential cellular functions. MT dynamics are finely regulated in the cell, . It has been hypothesized that inadequate regulation of neuronal MT dynamics may underlie neuronal cell death in Alzheimer's and related dementias, . Additionally, drug induced modulation of MT dynamics underlies the effectiveness of various anti-cancer drugs, such as Taxol, . For these and a host of basic biology issues, the regulation of MT dynamics is a very active area of research in modern cell biology.
In this work, we propose a powerful approach for analyzing MT dynamic behavior. Briefly, we use an automated tracking method for measuring MT dynamics, which are then modeled as MT behavior patterns by Hidden Markov Models. The proposed methods go beyond the traditional analysis capabilities and offer new insights in investigating MT dynamic behavior.
The cytoskeleton of a eukaryotic cell consists of a network of fibers. MTs are one of the three principal types of cytoskeletal fibers. They are hollow cylindrical structures, 25 nm in diameter and up to several μ m in length, consisting of non-covalently bound tubulin protein subunits. MTs are constantly assembled and disassembled, making the cytoskeleton a dynamic system. MTs are critically involved in a number of essential cellular functions, such as chromosome segregation at mitosis and intracellular cargo transport. Additional background information on MT structure and function can be found in .
The growing and shortening dynamics of MTs are finely regulated, for example, by the action of MT-associated proteins (MAP) and MT-targeted drugs (MTD). A large body of evidence, reviewed by Feinstein and Wilson , suggest that cell viability requires that MT dynamics be properly regulated within a narrow range. Common conjecture is that certain diseases such as Alzheimers and cancer are at least correlated with the regulatory abnormalities in MT dynamics, [6–8]. Consequently, gaining a detailed mechanistic understanding of the regulatory activities of MAPs, [5, 6, 9], and MTDs, [3, 10], is a major focus of current research. A major challenge is assessing the activities of the large number of MAPs and their many isoforms, as well as the large number of MTDs and their many derivatives. For instance, the MAP tau consists of 441 amino acids, more than 25 of which can be phosphorylated in various combinatorial patterns. Whereas phosphorylation normally serves to regulate tau activity, excessive and abnormal phosphorylation correlates with cell death and dementia. Thus, to fully understand normal and pathological tau action, the regulatory effects of the many different combinational phosphorylation patterns of tau must be understood.
MTs are polar structures, possessing biochemically distinct minus and plus ends. Conventionally, the minus end of a MT is assumed to be fixed at the MT organizing center near the nucleus, and the other end -the plus end or the tip- is the dynamic end that is observed in most MT dynamics studies. Typically, in live cell studies, minus ends of the MTs are not visible because of the high density of MTs converging on the organizing center. Thus, in calculating the MT length, a point on the MT body is selected as a reference point, origin, after an initial observation of all frames in the time-lapse images, (Fig. 2).
Traditionally, time-lapse images of MT populations are collected following treatment with MTDs or MAPs. Dynamics parameters are then manually calculated from image sequences as follows. The positions of MT plus ends are manually tracked individually across all frames, (Fig. 1). MT lengths are approximated as the (Euclidean) distance between tracked tip positions and the origin, producing MT life histories or tracks, (Fig. 2). The change in MT length is computed between consecutive frames, and growth and shortening statistics are tabulated. Length changes below a threshold are marked as attenuation or pause, signifying undetectable change. Other biologically significant events are the conversion of a MT from a growing state to rapid disassembly, designated as a catastrophe, and a subsequent potential recovery from shortening to attenuation or growth, called rescue. To estimate the effects of a regulatory agent upon MT dynamics, these statistics are aggregated over a number of MTs from the same experimental condition. Resulting statistics of each condition are compared with the control behavior to quantify the effects of the examined agent on dynamics parameters.
In this fashion, regulatory effects of each individual agent are studied through a laborious set of tasks. Quantifying sufficient image data to achieve statistical significance and limited comparative capabilities in the presence of innumerous possible agents pose an enormous challenge to researchers. Example studies are [3, 5, 9–12].
Statistics obtained from the growth and shortening events treat these events independently, rather than as being part of a behavior pattern. For instance, a certain growth measurement is counted as the same event regardless of where it occurs in relation to preceding or subsequent events. Furthermore, studying event correlations between neighboring MTs are generally infeasible, despite potential biological significance.
There are no established non-manual methods for examining the similarities and differences in particular dynamic behaviors imposed by various agents. Furthermore, studying combined effects of multiple regulatory agents is difficult, due to the limitations imposed by the pairwise comparisons between experimental conditions. For example, consider a hypothetical MTD AB, derived from MTDs A and B. In order to understand the contributions of A and B, multiple individual experiments must be conducted. Therefore, quantifying behavioral similarities across experimental conditions may provide essential guidance in constructing hypotheses.
In this work, we propose an automated tracking and analysis method to address the limitations mentioned above. The tracking component provides behavioral features for subsequent analysis. We define the MT dynamic behavior as a sequence of changes in MT length over equal time intervals. Experimental conditions may exhibit a number of behavioral patterns, which are estimated in parametric form by a mixture of Hidden Markov Models. By using a model-based clustering technique, we propose to analyze the constituent parts of MT behavior in each experimental condition. Thus, each experimental condition can be described as a mixture of behaviors exhibited by different MT populations. Through estimated average behavior patterns, we introduce a probabilistic behavioral distance measure between experimental conditions. Furthermore, parameters of individual models may present significant information about the properties of corresponding behavioral patterns. We describe how model-based analysis can be effective in addressing the above limitations (see Discussion).
We present statistical models of MT behavior that are estimated using automatically tracked MT dynamics data. As a comparison, we provide models of manually collected MT tracks. We describe the results of automated tracking using visual samples and associated errors.
Quantitative results of MT tracking were given in . Evaluations against manually tracked data shows that the mean and the standard deviation of tracking error are 2.85 and 4.36 pixels, respectively. This error level is acceptable for biological studies. Recall that the MT width is 25 nm (see Background), which appears as curves that are 3 pixels in width. Thus, a growth or shortening event that is less than 3 pixels would correspond to an event that is too small to quantify reliably, and is considered as attenuation.
We note that the tracking performance is sensitive to the accuracy of initial tip detection step. Furthermore, the proposed approach requires multiple tips to be detected for reliable extraction of MT tracks by design. In other words, the tracking performance may be adversely affected in tracking MTs individually, which may limit the ability to track a particular MT in a cell. Finally, intersecting MTs may steal the body trace, as the geodesic distance will favor higher intensity levels, (Fig. 3d – 3f). This issue could be addressed with further constraints on the MT orientation and curvature. However, in this work, we limit the behavioral features to the observed change of length in the MT plus end, which only requires consistent estimation of the MT body.
While tracking performance may be improved as a consequence of higher image quality and suitable algorithms targeting frequent intersections, deformations, and intensity variations. In its current state, automated tracking can track and quantify 10 times more MT tracks per image sequence than manual methods. With this increase in analyzable data volume, we are able to estimate behavior models for different experimental conditions. Estimated statistical models of MT dynamic behaviors are presented in the next section.
In this work, we used MT time-lapse live cell images from . The authors of  investigate the hypothesis that resistance to Taxol may involve altered sensitivity to different tubulin isotypes. Chinese hamster ovary (CHO) cells were microinjected with rhodamine-labeled tubulin. A total of 111 sequences were acquired using fluorescence microscopy with a 100× objective lens (1000× magnification). 25 frames per sequence were captured at 4 second intervals, from five different conditions.
Growth and shortening rates were computed as the differences of a MT lengths between consecutive frames, measured in pixels. Thus, each track consists of an observation sequence composed of 25 points in time. Resulting observation sequences were in the range [-13.03, 11.22] pixels, where (-) and (+) denoting shortening and growth rates, respectively.
Experimental conditions and number of tracks collected, automatically (AT) and manually (MT).
β III-tubulin expressed, no Taxol
β III-tubulin expressed, plus Taxol
β III-tubulin not expressed, plus Taxol
β I-tubulin expressed, no Taxol
β I-tubulin expressed, plus Taxol
Correct classification rates for EX:A, B, C.
Correct AT (%)
Correct MT (%)
The same set of experiments were repeated with manually tracked MT data. Separation results are shown in Table 2. Similar classification rates with the automatically tracked experiments confirm the automated tracking as well as the applicability of model based analysis.
Change in correct classification rates vs. the number of models from EX:B. Separation peaks at W = 3.
Example emission distributions of λ1 from EC4, and λ2 from EC5.
Example emission distributions of λ1 from EC4, and λ2 from EC5.
Estimated models can provide more descriptive information about the behavior patterns than what is available through manual methods: (i) typical growth and shortening states of the modified behavior, and (ii) the transition probabilities between these states. For example, as a direct comparison with manual methods, besides the traditional catastrophe and rescue frequencies, transitions from small to larger events of the same type can be quantified. In essence, characteristics of behavior patterns are parametrically encoded in models, which can then be used in generating these behaviors. We describe further model-based analysis capabilities in the next section.
By quantifying behavioral comparisons between regulatory agents, studying combined effects of multiple regulatory agents may be guided with enhanced predictions. We envision a repository of MT dynamics data that can be probabilistically queried for behavioral similarities for a new regulatory agent, an isotype, or a combination. This can be done by evaluating p(O|EC) for an experimental condition EC, or evaluating p(O|w) for behavioral pattern w. Assuming that the tracking and modeling tasks were undertaken, a MT image database would contain a collection of individual MT tracks and model parameters representing w, in addition to original image sequences. Model based content retrieval provides additional advantages in query design. Hypothesized behaviors can be created and queried by manually selecting model parameters. Alternatively, query models can be estimated from a subset of MT tracks in the database.
This analysis provides the researcher with visual cues about regional dynamics within a cell. This may be especially important in studies of polarized cell types, such as neurons, where specific regional regulation of dynamics is critical to processes such as outgrowth and transport. Behavioral comparisons in adjacent populations may provide insight to the inner workings of flux between the soluble and polymeric tubulin fractions within the cytoplasm. The ratio between these two functionally distinct, but co-dependent phases may indicate cell-autonomous or drug-influenced regulation.
MT dynamics research seeks to understand the complex mechanisms that underlie cytoskeletal responses to changes in environmental conditions. A clear understanding of the regulation of MT dynamic behavior may elucidate causal factors in various diseases and may reveal new therapeutic targets and strategies. In this work, we introduce novel data collection and analysis capabilities based on computer vision and machine learning tools. With the proposed methods, researchers can study MT dynamics with improved spatial and temporal quantification.
The most notable contribution of the proposed method is the novel analysis capabilities that are beyond the current state-of-the-art. Other contributions are the improvements over the manual data collection methods, such as higher accuracy (length along the MT vs Euclidean estimate), increased number of analyzable MT tracks, and objective consideration of all MT tracks at a fraction of the normally required time. Our preliminary results support manually established findings, and show that automated analysis of spatial and temporal patterns offers previously unattainable insights. Most notably, the standardization of data collection and analysis facilitates a comparative platform for future biological research.
As the volume and number of dynamics datasets has increased in recent years, similarities between the behavioral influence of MAPs and MTDs upon dynamics have emerged, leading to speculation of similar mechanisms. Dynamics models may facilitate the union of previously isolated MAP and MTD datasets, furthering our understanding of regulatory mechanisms of MTs.
Despite the difficulties inherent in fluorescence imaging, the proposed approach confirms manual findings in both track computation and in analysis. For example, due to photobleaching, observation durations were generally limited to only a few minutes with very low signal-to-noise ratios in images. With emerging techniques in microscopy and probes, such as the tip-binding proteins (EB1), much longer acquisition times will be possible with superior image quality. Our goal is to track all MTs in live cell images at longer durations. In this direction, the tracking method can be improved by reliably identifying all MTs individually. The nature of live cell MT images requires that frequent intersections, abrupt intensity variations on a single MT body, and focusing issues must be addressed adequately.
The proposed analysis system evaluates MT dynamic behavior as a function of entire MT life histories through estimating statistical models from observations. A number of MT tracks per experimental condition is necessary for reliable estimates of model parameters. Thus, an automated tracking procedure was used in data collection.
Previous work on automated MT detection and tracking include [16–18]. In [13, 19], we described our tracking approach for live cell images and introduced the idea of model based analysis. In , the authors extract MT plus ends using a MT body and a tip model in a multi-scale operation. In  and , MTs are traced in segments from initially selected points and subsequently tracked. In , MTs are searched in a constrained space for tracking in subsequent frames.
An automated MT tracking method should address the following: (i) highly variable tubule shapes, (ii) accurate estimation of the MT length considering the nonlinear shape, (iii) frequent occlusions and intersections from surrounding MTs, and (iv) low signal-to-noise ratios with spatial and temporal variations in illumination.
To address these issues, we consider MTs as flexible open curves in the image plane, with a fixed minus end near the nucleus and a dynamic plus end. Formally, a single MT is modeled by the open curve C(s), where s ∈ [0, 1] is the curve parameter. The goal of the MT tracking task is to estimate the MT length by locating the tip and tracing the deformation of the MT body, in every frame.
where the derivative of the Gaussian is taken along orientations θ at position (x, y), and σ is chosen as the average MT width. The maximum filter response, I f (x, y), is then binarized to generate a mask showing MT polymer mass. The binary mask is used for determining tip candidates in each frame. Example tip detection results from consecutive frames are shown in (Fig. 9).
Once the tip candidates are located in each frame, correspondences are established between frames by using a multi-frame graph matching algorithm. The reasoning behind formulating the correspondence as a graph optimization problem is that by matching multiple tips at once, occasional spurious tips are removed. Furthermore, the graph matching algorithm provides the flexibility of skipping frames, which handles missing tips between frames.
where and u, v are two vertices in G for which the similarity is computed as in Eq.(4). Note that between two frames the best tracks can be computed as the maximum match of a bipartite graph. However, for multiple frames, the problem becomes NP-hard. Here, we adopt the approximation proposed in .
The described method is sufficient to track MT tips between different frames. However, without tracing the MT body, the best estimates of MT growth and shortening would be limited to Euclidean approximations between tip positions, (see Current analysis method). Since in live cell images, the MT body is typically non-linear, this approximation is a rough one in practice. Instead, we determine the MT body length in all frames.
In essence, we compute the MT body length along the body in each frame and determine the growth and shortening as consecutive length differences. Given the tip positions in each frame, we estimate the deformable curve constituting the MT body between these tips and a fixed point along the MT body. Note that the fixed point does not have to lie on the body of a specific MT for the purposes of computing the growth and shortening. In cases where the fixed point lies on another MT rather than the MT being measured, the resulting change in length is still a better estimate than the Euclidean case, so long as the fixed point taken consistently across frames. Details of fixing this point can be found in . Due to the constant deformations, the fixed point location may exhibit small variations, (Fig. 5d – 5f). This is the major contributor of errors in length estimation between frames. Finally, based on the estimated plus and minus ends of the MT, the MT body is extracted using active contours with ridge features.
A number of studies examined physical models for MT structure and dynamics. We refer the interested reader to [21–23], and the references therein, for a review of previous models of MT dynamic instability. For example, in , the authors use a simulation model to investigate the fluctuations in tubulin concentration in relation to MT dynamics. In contrast to previous dynamics models, we propose using machine learning methods for modeling various MT behavior patterns occurring in different experimental conditions.
Automated tracking is sufficient to quantify traditional dynamics parameters. We propose an analysis approach targeting behavioral information beyond what is provided by the traditional parameters. We begin with including contextual information in time. In other words, as opposed to parsing the growth and shortening events out of MT tracks (life histories), we keep the MT tracks intact. Therefore, each MT track is treated as an observation from some behavior pattern. For example, the tracks in (Fig. 12, top row and middle row) are observation instances from different behavior patterns. Thus, if g denotes a small, and G denotes a large growth events, then the observed tracks, ggggGGGG and ggGGggGG should be treated as different behaviors even if the average growth rates may be equal. This definition of a MT behavior pattern leads to new analysis capabilities. Each behavior pattern can be described by a model. Subsequently, estimated models are used in analyzing MT dynamic behavior; for instance, in evaluating dynamic similarities between MT populations.
In modeling the MT dynamic behavior, biological insights provide essential guidance. Similar behavior patterns are known to be shared between different experimental conditions, while MT populations within a cell may exhibit dissimilar patterns. Thus, modeling design should handle expected variations of behavior within each experimental condition, and similarities between different experimental conditions.
Formally, we denote each experimental condition by EC, consisting of groups of behavior patterns, w. All experimental conditions have a known label, while patterns making up a condition are unknown. The problem is to estimate a model λ for each pattern w, such that differences between EC i and EC j , i ≠ j, are emphasized, while each pattern may occur in different experimental conditions, w ∈ EC i and w ∈ EC j . Note that our formulation calls for a discriminative approach between EC, while descriptive models of w is the goal across different EC's.
A well known class of models used in representing activity is the Hidden Markov Models (HMMs). In the past, they have been used in numerous applications, most notably in speech recognition, , and in genomic sequence analysis, [25–28]. Particularly in activity context, HMMs were used in activity recognition , abnormal activity detection, gesture recognition, and American Sign Language recognition. In the next section we review the essentials of HMMs, while referring the reader to  for further details.
HMMs are probabilistic generative models estimating the statistics of a process from observation sequences generated by that process. The modeled process is assumed to be not directly observable, thus hidden states capture statistics of the process, subject to stochastic constraints. In practice, hidden states generally correspond to certain physical characteristics of the process. Detailed information on modeling with HMMs can be found in [24, 28]. Concisely, HMMs, denoted by λ, are described by parameters λ = (π, A, B), where π is the state priors, A is the transition, and B is the emission probabilities. Given an observation sequence O = (o1, o2,..., o T ), where t = 1..T denotes time, and a model λ = (π, A, B), the quantity P(O|λ) can be computed efficiently. Given a set of observation sequences, estimating the parameters of λ is generally performed using maximum likelihood methods, while discriminative techniques were suggested in classification tasks, [30, 31].
From the biological perspective, classification of tracks to respective EC is not the end goal for dynamics analysis since labels of EC are known a priori. However, estimated behavior models, λ, provide novel analytical capabilities. Furthermore, model parameters may reveal further insights into MT dynamic behavior. Our formulation of the problem aims to extract behavior patterns through estimating λ, while discriminating between different EC. In doing so, we employ the classification score as our measure of model reliability. The problem description motivates us to use a model based clustering approach to estimate a λ for each w. HMM based clustering methods are discussed in .
After parameter estimation, each EC is represented by a mixture of λ where dynamics variations within each EC are modeled by the components of the mixture. In this sense, each λ models the (pseudo-)center of a w, the component behavior patterns contributing to the resulting behavior in respective EC. The estimation of λ is primarily a modeling task, while discrimination between w is handled by clustering the observations, MT tracks, into behavior patterns represented by the respective w.
is maximized through
assign o to cluster C w such that
w = arg max w' log p(o|λ w' )
train λ w on C w , w = 1..W
In each iteration of the algorithm, observation o is assigned to maximally likely cluster C w , whose center λ w is re-estimated using the new members of C w . The iterations are terminated when no significant increase in the overall likelihood is observed.
Note that the decision is conditional on λ w , EC, representing contributions of each member λ w of EC.
We would like to thank to Motaz El-Saban and Emre Sargin for their contributions in the automated tracking of microtubules. This study was funded by Center for Bioimage Informatics under grants NSF-ITR 0331697(BSM, KR, LW, SCF), NIH-ROI NS35010(SCF), NIH-ROI NS13560(LW).
This article has been published as part of BMC Cell Biology Volume 8 Supplement 1, 2007: 2006 International Workshop on Multiscale Biological Imaging, Data Mining and Informatics. The full contents of the supplement are available online at http://www.biomedcentral.com/1471-2121/8?issue=S1
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.