Michael W Held

Michael W Schmitz

Bernd H Fischer

Beate Thomas Walter

H Neumann

Michael H Olma

Matthias H Peter

Jan H Ellenberg

Daniel W Gerlich
email: daniel.gerlich@bc.biol.ethz.ch

Michael H A Schmitz

Beate Neumann

CellCognition: time-resolved phenotype annotation in high-throughput live cell imaging

Keywords: Live cell imaging, RNAi screening, GFP, Machine Learning, Image analysis, hidden Markov model, Cdc20, mitotic exit RNAi: RNA interference; siRNA: small interfering RNA

come   L'archive ouverte pluridisciplinaire

Introduction

The availability of RNAi technology for high-throughput gene inactivation experiments, fluorescent protein labeling, and automated microscopy has opened a new era of screening possibilities in higher eukaryotes 1 . Indeed, imaging-based RNAi and chemical compound screening has become one of the most important discovery tools for the identification of new gene function, for example in the regulation of DNA damage and repair 2 , endocytosis [START_REF] Collinet | Systems survey of endocytosis by multiparametric image analysis[END_REF] , mitosis [START_REF] Sonnichsen | Full-genome RNAi profiling of early embryogenesis in Caenorhabditis elegans[END_REF][START_REF] Goshima | Genes required for mitotic spindle assembly in Drosophila S2 cells[END_REF][START_REF] Neumann | Phenotypic profiling of the human genome by time-lapse microscopy reveals cell division genes[END_REF] .

Imaging-based screens typically assay altered incidence of cells with specific features within a population of fixed, fluorescently labeled cells. The development of computational methods for the automated annotation of high-throughput imaging data was key to establish microscopy-based screening as a routine technology in a wide research community. Especially machine learning for supervised classification of cellular morphologies is one of the most powerful annotation strategies [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF][START_REF] Loo | Image-based multivariate profiling of drug responses from single cells[END_REF][START_REF] Conrad | Automatic identification of subcellular phenotypes on human cell arrays[END_REF][START_REF] Glory | Automated subcellular location determination and high-throughput microscopy[END_REF][START_REF] Harder | Automatic analysis of dividing cells in live cell movies to detect mitotic delays and correlate phenotypes in time[END_REF][START_REF] Zhou | A novel cell segmentation method and cell phase identification using Markov model[END_REF] .

Many biological processes depend on stochastic events and occur in an unsynchronized and transient manner, which limits the applicability of single time point assays. Particularly, complex dynamic processes such as cell division or intracellular trafficking demand for time-resolved live cell imaging [START_REF] Gerlich | 4D imaging to assay complex dynamics in live specimens[END_REF] . Automated microscopes now enable live imaging with high throughput and spatio-temporal resolution 1,[START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF][START_REF] Schmitz | Automated live microscopy to study mitotic gene function in fluorescent reporter cell lines[END_REF] . Computational analysis of such data is challenging and existing machine learning and classification approaches do not provide sufficient accuracy to correctly annotate cellular trajectories with multiple time points. Published live imaging-based RNAi screens scored phenotypes either exclusively at the cell population level [START_REF] Neumann | Phenotypic profiling of the human genome by time-lapse microscopy reveals cell division genes[END_REF][START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF] , or relied on visual evaluation of single cell dynamics [START_REF] Sonnichsen | Full-genome RNAi profiling of early embryogenesis in Caenorhabditis elegans[END_REF] . However, cell population analysis cannot detect stochastic and transient phenotypes, and visual interpretation of morphological dynamics is very time consuming and often unreliable.

To improve the classification accuracy of machine learning methods, the temporal context can be taken into account. For example, if the biological process underlying an assay is well known, a biological model can be explicitly defined in an error correction scheme that suppresses illegitimate stage transitions. This has been applied to the pattern of mitotic chromatin morphology changes [START_REF] Harder | Automatic analysis of dividing cells in live cell movies to detect mitotic delays and correlate phenotypes in time[END_REF][START_REF] Zhou | A novel cell segmentation method and cell phase identification using Markov model[END_REF] . However, temporal error correction based on biological a priori models limits the detection of unexpected phenotypic variations, and the adaptation to different biological questions requires reimplementation of the underlying models by the user for each new assay.

Here, we present CellCognition, an integrated computational strategy that combines machine learning methods for supervised classification and hidden Markov modeling to measure morphological dynamics in live cell microscopic movies. Our error correction method does not require a priori definition of the temporal progression, which enables its application to a wide range of assays and phenotypic variations. We demonstrate efficiency and sensitivity of the methodology in various assays and perturbation conditions.

Results

High-throughput imaging of cellular dynamics

To visualize morphological dynamics of various cellular structures, we generated a collection of human HeLa reporter cell lines stably expressing different combinations of fluorescent markers. All cell lines expressed a red chromatin marker (core histone 2B (H2B) fused to mCherry). In this background, we co-expressed markers for microtubules (mEGFP-α-tubulin), the Golgi apparatus (Galactosyl transferase (GalT) fused to EGFP), or DNA replication factories (proliferating cell nuclear antigen (PCNA) fused to mEGFP). This diverse set of secondary markers (Fig. 1a) provides a well-suited test case for the implementation of a generic annotation method. With these cells, we performed multi-location time-lapse imaging on an automated widefield epifluorescence microscope [START_REF] Schmitz | Automated live microscopy to study mitotic gene function in fluorescent reporter cell lines[END_REF] . We typically recorded 96 movies in parallel, with a temporal resolution less than 5 min over a total duration of 24 h, generating datasets of about 100,000 images, or 200 Gigabyte, per day and microscope. The analysis of such a single experiment requires annotation of up to 25 million cellular morphologies derived from about 260,000 objects per movie with a 10x microscope objective.

Machine learning and classification of morphologies

Timing measurements in live cell imaging data are often based on the progression through distinct morphologies that relate to specific biological states. An excellent example for this is mitosis, for which the chromatin morphology can be used to annotate the canonical mitotic stages (Fig. 1b; Supplementary Movie 1). We decided to use this classic assay as a test case to measure timing events at the single cell level.

We first implemented a canonical strategy for automated annotation of morphological classes [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF][START_REF] Loo | Image-based multivariate profiling of drug responses from single cells[END_REF][START_REF] Conrad | Automatic identification of subcellular phenotypes on human cell arrays[END_REF][START_REF] Boland | the patterns of all major subcellular structures in fluorescence microscope images of HeLa cells[END_REF] , based on object detection, multivariate feature extraction, and supervised machine learning (Fig. 1c). We used local adaptive thresholding [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF] , followed by a watershed split-and-merge segmentation error correction [START_REF] Wahlby | Combining intensity, edge and shape information for 2D and 3D segmentation of cell nuclei in tissue[END_REF] to detect individual cells at an accuracy of 95.7% (n = 1876 objects; 2.6% over-segmented (falsely cut objects); 1.7% under-segmented (falsely merged objects)). A set of 186 quantitative features [START_REF] Walker | Statistical geometric features-extensions for cytological textureanalysis[END_REF][START_REF] Haralick | Dinstein & Shanmugam Textural features for image classification[END_REF] (Supplementary Table 1, Supplementary Fig. 1) describing texture and shape was then calculated for each object. Next, a support vector machine classifier [START_REF] Boser | A training algorithm for optimal margin classifiers[END_REF] was trained for the discrimination of 8 different object morphologies (Fig. 1b; interphase, six different mitotic stages, and apoptosis). These classes were defined by manual annotation of 28 to 195 example objects. The match between human and computer annotation was 94.6% (mean of all classes; five-fold cross-validation), ranging between 75.0% for the early anaphase class, and 99.0% for interphase class (Fig. 1d). This performance is similar to a number of previously reported supervised machine learning applications [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF][START_REF] Conrad | Automatic identification of subcellular phenotypes on human cell arrays[END_REF][START_REF] Harder | Automatic analysis of dividing cells in live cell movies to detect mitotic delays and correlate phenotypes in time[END_REF][START_REF] Wang | [END_REF] . Next, individual cells were tracked over time by a nearest-neighbor algorithm that supports trajectory splitting (e.g., cell division) and merging (e.g., cell-to-cell fusion). The automated tracking matched 99.8% of the human annotated object-over-frame connections (n = 1942), again comparable to the performance of previous studies on cell tracking [START_REF] Harder | Automatic analysis of dividing cells in live cell movies to detect mitotic delays and correlate phenotypes in time[END_REF][START_REF] Chen | Automated segmentation, classification, and tracking of cancer cell nuclei in time-lapse microscopy[END_REF] .

The overall accuracy of the individual computational steps appears high. However, considering >500 frames per cell trajectory for our time-resolved datasets, almost no error-free trajectories were obtained by this approach (Fig. 1e, Supplementary Movie 2).

Detecting scarce events in long-term movies

Mitotic events are scarce in comparison to the much longer duration of interphase (Fig. 1e). To improve the sensitivity for mitotic stage annotation, we automatically selected mitotic events based on a morphology class sequence motif of prophaseprometaphase. This yielded a sub-graph highly enriched for mitotic events (Fig. 2a; Supplementary Movie 3; 81.5% of all mitotic events were automatically extracted; n = 294 mitotic events in three movies). This set of trajectories contained 2.1% misclassifications per object (a posteriori compared with human annotation).

Untrained biological users may annotate the classifier training set less reliably. To test the sensitivity of the support vector machine towards annotation errors, we randomized the labels on an increasing fraction of training objects, and measured the overall classification accuracy (Supplementary Fig. 2). Surprisingly, randomization of the labels on 50% of the training objects reduced the overall annotation accuracy only slightly below 90%. This demonstrates that classification by support vector machine is relatively insensitive to annotation errors.

Hidden Markov model for time-lapse imaging

Single object-based machine learning and classification does not take the temporal context into account. However, objects with ambiguous morphologies occur within a typical context of preceding and following morphologies, which could help to derive correct annotation. This could be particularly relevant for gradual morphology changes at stage transitions, where single object-based classification is relatively inaccurate (e.g., interphase -prophase -interphase -interphase -prophaseprometaphase, see Fig. 2b or Supplementary Movie 3)

We reasoned that taking the history of a cell into account might provide a means to correct for such noise at stage transitions, as well as confusion between closely related morphology classes. We assumed that the true state of a cell at a given time point (the mitotic stage in this assay) is not known, but that it correlates with an observed state (the morphology class prediction probabilities). We further assumed that the progression to the next state entirely depends on a given present state. This fulfils the criteria for a hidden Markov model, which can be used for error correction in timeresolved data [START_REF] Durbin | Biological sequence analysis: probabilistic models of proteins and nucleic acids[END_REF] .

We built a model with five components: 1) hidden states, representing the true morphology classes (for example, mitotic stages), 2) observed states (the class prediction probability vectors of the support vector machine), 3) probabilities of hidden state transitions, 4) observation probabilities, and 5) initial probabilities of hidden states. All elements of this model were computationally derived from the data without further user interaction. The hidden states were defined by the initial class annotation, as described above (Fig. 1b). The observed states were derived from the support vector machine as a vector of class prediction probabilities for each time point. The hidden state probabilities were initialized at the first time point by the support vector machine predictions. Transition probabilities between hidden states were calculated based on the support vector machine prediction probabilities of all cellular trajectories per experimental condition (Fig. 2c,d), and the observation probabilities between hidden and observed states were estimated based on the confusion matrix of the support vector machine. We derived the overall maximum likelihood path for the progression through mitosis by the Viterbi algorithm [START_REF] Viterbi | Error bounds for convolutional codes and an asymptotically optimum decoding algorithm[END_REF] (thick black line in Fig. 2e). This increased the overall per-object accuracy to 99.0%.

Iterative learning of transition probabilities by the expectation-maximization algorithm [START_REF] Baum | A maximization technique occurring in the statistical analysis of probabilistic functions of Markov chains[END_REF][START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF] did not improve prediction accuracy (98.1% after five iterations). We suspected that the confusion matrix overestimates observation probabilities, as classes that are difficult to discriminate (prophase and early anaphase) were over-represented in the annotation data. We therefore tested the performance of temporal error correction with lower error rates in the observation probabilities (0.1% for all transitions). Indeed, this eliminated noise at state transitions and corrected single frames of misclassified objects even more efficiently, yielding overall accuracy of 99.4% per object, and 91% completely error-free trajectories (n = 100 trajectories; 4,000 objects; Fig. 2f, Supplementary Fig. 3; Supplementary Movie 4).

We next tested if incorporation of a priori biological knowledge on state transitions further increases the annotation accuracy. Specifically, we constrained the state transition graph to the forward direction of three consecutive classes, and defined apoptosis as a terminal state (Supplementary Fig. 4a,b). The probability matrix for constrained state transitions improved the error correction performance of the hidden Markov model to 99.7% per object, yielding 94% completely error-free trajectory annotations (n = 100 trajectories; 4,000 objects; Supplementary Fig. 4c).

Temporal error correction by the hidden Markov model is expected to depend on good estimates of the predicted morphology classes. We therefore investigated the robustness of temporal error correction towards simulated classification noise. We randomized the class prediction probability vectors of an increasing fraction of objects, then learnt the hidden Markov model on the noisy trajectories, and applied it to correct classification errors (Supplementary Fig. 5). Comparison with manually annotated data demonstrated that the hidden Markov-based error correction improved the overall accuracy at all noise levels.

We also tested if the temporal error correction was sensitive to changes in the timelapse interval by generating trajectories sampled to every 2 nd up to every 6 th time point (Supplementary Fig. 6). Comparison with the manually annotated labels showed that the hidden Markov model increased the overall annotation accuracy at all sampling intervals.

In conclusion, hidden Markov modeling provides a robust and efficient means to eliminate misclassifications and noise at morphology state transitions. The combination of mitotic event selection and hidden Markov error correction reduced the per-object error rate about 10-fold below single time point-based classification.

Generic strategy for annotation of cellular dynamics

We next used our tools for other assays and fluorescent markers. We were particularly interested in simultaneous analysis of multiple markers in the same cell, for example, to address temporal coordination of mitotic processes. We defined cytoplasmic areas based on their relative position to the chromatin marker, using non-overlapping region growing of the contours derived from the chromatin channel (Supplementary Fig. 7a,b). While this may be less precise than segmenting in the secondary channel, it proved to be robust over many different assays and was insensitive to temporal dynamics (see Fig. 3,4, and below). Tracking results of the primary channel were applied to the secondary channel, and all subsequent analysis of temporal dynamics was performed independently for primary and secondary channels, as outlined above (see Supplementary Fig. 7c for workflow).

We first applied our methods to movies from cells expressing mEGFP-α-tubulin to annotate mitotic spindle assembly and disassembly (Fig. 3a and Supplementary Movie 5), and to movies from cells expressing GalT-EGFP to study mitotic breakdown and reassembly of the Golgi apparatus (Fig. 3b, and Supplementary Movie 6). We trained classifiers for six (α-tubulin), or five (GalT) distinct morphology classes. The mean accuracy of object class predictions was 96.5% for mEGFP-αtubulin, and 97.3% for GalT-EGFP (5-fold cross-validation, computational versus visual scoring). This yielded 55% (α-tubulin), or 38% (GalT) completely error-free trajectories. By hidden Markov model error correction, the accuracy increased to 89% completely error-free trajectories for α-tubulin (Fig. 3d and Supplementary Movie 7), and 90% for GalT (Fig. 3e and Supplementary Movie 8; n = 100 for both assays; corresponding H2B-mCherry annotations are shown in Fig. 3g,h).

To apply our methods to non-mitotic cellular dynamics, we next annotated the timing of S-phase progression. We imaged a HeLa cell line stably expressing H2B-mCherry and EGFP-PCNA, a marker for DNA replication foci, which visualizes a characteristic pattern of morphology changes during S-phase progression (Fig. 3c and Supplementary Movie 9). We trained classifiers for six distinct PCNA morphology classes, and established a hidden Markov model for error correction. This yielded 98.2% correctly annotated objects and 90% completely error-free trajectories (n = 100 trajectories containing 15,000 objects; Fig. 3f and Supplementary Movie 10, see Fig. 3i for H2B annotations of same cells). The high performance in this diverse set of assays demonstrates a generic applicability of our computational methods.

Quantitative phenotyping and kinetic measurements

Our methods were designed for the detection of timing phenotypes. We therefore established perturbation conditions that are known to delay or shorten particular stages of mitosis. First, we used the microtubule-depolymerizing drug Nocodazol, which arrests cells in prometaphase by permanent activation of the spindle checkpoint (Fig. 4a; Supplementary Movie 11). This was reliably detected by our computational tools (96.2% completely error-free annotated trajectories, n = 154; Fig. 4b).

Next, we depleted the essential spindle checkpoint component Mad2 by RNAi, which is known to accelerate the timing from mitotic entry until anaphase onset in HeLa cells by about two-fold [START_REF] Meraldi | Timing and checkpoints in the regulation of mitotic progression[END_REF] (Fig. 4a; Supplementary Movie 12). We evaluated the accuracy of automated timing measurements, scoring the time from prometaphase until anaphase onset based on the chromatin marker (cells that did not segregate chromosomes were omitted). Automated measurements of 47.2 ± 20.0 min (mean ± s.d.; n = 195) in control cells did not significantly differ from manual annotation of the same dataset (48.5 ± 18.0 min; two-sided Mann-Whitney-Wilcoxon test: p = 0.12). Automated timing measurements in Mad2 RNAi cells demonstrated mitotic acceleration (13.0 ± 3.6 min), again matching well measurements by manual annotation (12.4 ± 3.4 min; two-sided Mann-Whitney-Wilcoxon test: p = 0.23). As expected from the known biological function of Mad2, the mitotic acceleration in Mad2 RNAi cells was mainly due to a shortened metaphase stage (1.6 ± 1.1 min in Mad2 RNAi cells; 36.5 ± 16.6 min in control; Fig. 4b).

Simultaneous measurements of morphological dynamics and the state of regulatory factors provide a powerful approach for mechanistic dissection of perturbation phenotypes. Here, we combined the annotation of mitotic stages with kinetic measurements of Securin degradation, which is required for anaphase initiation [START_REF] Hagting | Human securin proteolysis is controlled by the spindle checkpoint and reveals when the APC/C switches from activation by Cdc20 to Cdh1[END_REF] (Fig. 4a; Supplementary Movies 11-13). In the normalized degradation kinetic profiles (Fig. 4c), we found that the Securin-mEGFP degradation in control cells initiated briefly before anaphase (compare Fig. 4b andc), consistent with spindle checkpoint inactivation at this stage. In nocodazol-arrested cells, almost Securin-mEGFP remained stable within the measurement period of 138 min, consistent with an efficient and permanent activation of the spindle checkpoint. Securin-mEGFP degradation in Mad2 RNAi cells initiated directly after mitotic entry, at a stage where chromosomes were still in prometaphase configuration, indicating that the anaphasepromoting complex was activated before complete chromosome congression, as expected for a compromised spindle checkpoint function. In conclusion, these experiments demonstrate accurate timing phenotype annotation in RNAi-and drugperturbed cells.

RNAi screen for mitotic exit regulators

To test the sensitivity and performance of our computational methods in a highthroughput application, we performed a screen for regulators of mitotic exit. Specifically, we aimed to identify regulators of post-anaphase stages of mitosis, for which RNAi phenotypes have not been reported so far. Mitotic exit control is well understood in budding yeast, yet it is unclear if homologues of the yeast factors also control mitotic exit in higher eukaryotes [START_REF] Bollen | Mitotic phosphatases: from entry guards to exit guides[END_REF] . We therefore designed a library of 283 siRNA targeting 93 candidate regulators, including all known human genes with homology to budding yeast mitotic exit regulators and some additional genes known to be involved in mitotic regulation (see Supplementary Table 2). As an assay for mitotic exit timing, we scored the timing from anaphase onset, based on the chromatin marker H2B-mRFP, until postmitotic nuclear envelope reassembly, based on the nuclear import substrate IBB-EGFP (Fig. 5a; Supplementary Movie 14).

For solid-state transfection of siRNAs into HeLa cells, we used a high-density transfection array with 300 spots of different siRNA transfection solutions printed to the glass surface of a chambered coverslip [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF] . We seeded the cells onto this array and 20 h later started parallel imaging of 108 movies per experiment, for a total duration of 46 h and with 3.7 min time resolution. We automatically annotated the mean mitotic exit timing per experimental condition within the 1.6 TeraByte data containing 646'754 images and 16'314 mitotic events. Only one siRNA delayed mitotic exit above a z-score threshold of 3.0 (Fig. 5b, Supplementary Fig. 8a; 6.8 ± 2.0 min mean ± sd; n = 50 mitotic events). This oligo depleted the anaphase promoting complex co-activator Cdc20, as validated by Western Blotting (Supplementary Fig. 8b). The specificity of the phenotype was confirmed in two additional replicas with standard liquid phase transfection, and with an additional siRNA (Fig. 5c).

To test if Cdc20 was required for other cellular reorganization processes during mitotic exit, we assayed chromosome decondensation and mitotic spindle disassembly. High resolution confocal time-lapse imaging of cells co-expressing H2B-mCherry and mEGFP-α-tubulin (Fig. 5d,e, and Supplementary Movies 15 and 16) showed that 100% (n = 30) of control cells started chromosome decondensation within 14 minutes after chromosome segregation, whereas only 54% (n = 36) did so after Cdc20 depletion. 31% (n = 36) of Cdc20-depleted cells started kinetochore fiber spindle disassembly 7 minutes post anaphase onset, in contrast to 87% (n = 30) in control cells. These data suggest a requirement of Cdc20 for various cellular processes leading to postmitotic reassembly of interphase cells. This is unexpected given that Cdc20 has so far been thought to act mainly at pre-anaphase stages of mitosis, and it has not been noticed in previous phenotypic analysis of Cdc20 RNAi cells [START_REF] Wolthuis | Cdc20 and Cks direct the spindle checkpoint-independent destruction of cyclin A[END_REF] .

Discussion

In this study, we present CellCognition, a computational framework for time-resolved single-cell assays in high-throughput imaging applications. Building on existing machine learning methodologies, the design of a generic workflow for annotation of morphological dynamics faced two main challenges. First, the classification noise at continuous morphology stage transitions impairs coherent trajectory annotation.

Second, some biologically distinct classes appear morphologically similar, which leads to high classification confusion. By hidden Markov modeling, our methods efficiently correct both types of errors based on the temporal context. The hidden Markov models are learned individually for each experimental condition, without any human supervision. This allows the software to automatically adapt the error correction scheme to phenotypic deviations. Biological a priori knowledge to suppress state transitions that are assumed to be impossible can also be used to improve annotation accuracy [START_REF] Harder | Automatic analysis of dividing cells in live cell movies to detect mitotic delays and correlate phenotypes in time[END_REF][START_REF] Zhou | A novel cell segmentation method and cell phase identification using Markov model[END_REF] . Such explicit error correction schemes cannot be applied to new markers or assay systems without adaptation, and they may not apply to phenotypes with potentially altered stage progression. We find that the gain in accuracy by biological a priori constraints on the temporal progression is only minor. Our hidden Markov implementation models time series analysis in a high dimensional feature space with an intrinsic class-discriminant dimensionality reduction. This preserves context-specific structures, in contrast to principle component analysis as used in [START_REF] Zhou | A novel cell segmentation method and cell phase identification using Markov model[END_REF][START_REF] Wang | Context based mixture model for cell phase identification in automated fluorescence microscopy[END_REF] , which may explain the large gain in accuracy compared to the previous implementations (see Supplementary Tables 3 and4). Compared to the models by [START_REF] Harder | Automatic analysis of dividing cells in live cell movies to detect mitotic delays and correlate phenotypes in time[END_REF][START_REF] Zhou | A novel cell segmentation method and cell phase identification using Markov model[END_REF] , our model is the only one able to handle arbitrary relationships between phenotypic cell classes, providing a powerful and generic solution for time-resolved cellular phenotyping.

Using a variety of different structural markers, we demonstrate that our analysis methods can be used for a broad range of biological assays. We are not aware of any constraints that would preclude the use of our methods in other biological context, e.g., apoptosis or cellular differentiation. However, the texture and shape features implemented into our software do not enable assays relying on absolute object counts, for example in centrosome duplication assays. Also, assays scoring rapid intracellular dynamics would require integration of motion feature extraction methods into our published software source code. Supervised machine learning as in this study requires user-defined morphology classes. It is therefore not possible to detect aberrant phenotypic morphologies that do not occur in the control conditions used for annotating the classifier training set. This limitation may be overcome in future studies by implementing unsupervised machine learning methods for the analysis of image time series.

In conclusion, we present a powerful computational strategy for high-throughput phenotyping of single cell dynamics. Our methods are integrated into the platformindependent software package CellCognition, with graphical user interface and supporting high-throughput batch processing on computer clusters. CellCognition is published as open source software (current version 1.0.7 in Supplementary Software), along with high quality reference image data on http://www.cellcognition.org/. With the increased availability of live cell screening microscopes, we anticipate that timeresolved imaging assays will soon dominate a significant fraction of high content screening and systems biology applications.

Methods

Cell culture, RNAi and cell transfection arrays, and Western Blotting

HeLa 'Kyoto' cells were cultured in DMEM (Gibco) supplemented with 10% fetal calf serum (PAA Laboratories) and 1% Penicillin/Streptomycin (Invitrogen), and grown on LabTek chambered coverslips (Nunc) for live microscopy. All experiments were performed with monoclonal cell lines stably expressing combinations of the fluorescent markers as indicated throughout the manuscript. Live imaging was in DMEM containing 10% fetal calf serum and 1% Penicillin/Streptomycin, but without phenolred and riboflavin to reduce autofluorescence of the medium. Cell transfection arrays for live cell RNAi screening were produced and used as described in [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF][START_REF] Erfle | Reverse transfection on cell arrays for high content screening microscopy[END_REF] . All other RNAi interference experiments were performed using single RNAi duplexes (Qiagen) that were liquid phase transfected with either Oligofectamine (Invitrogen) or HiPerfect (Qiagen) as transfection reagent according to the manufacturers protocols.

Fluorescent reporter plasmid constructs

For efficient generation of cell lines stably expressing fluorescently tagged marker proteins, the genes were subcloned into pIRES-puro2 and pIRES-neo3 vectors (Clontech) that allow expression of resistance genes and tagged proteins from a single transcript. For details on the plasmids, see Supplementary Table 5.

Stably expressing cell lines

For generation of stably expressing cell lines, HeLa Kyoto cells were first transiently transfected using FuGENE6 (Roche) following the manufacturer's instructions. Cells were then seeded to clonal density and grown in culture medium supplemented with 500 µg/ml Geneticin (Invitrogen) and/or 0.5µg/ml Puromycin (Merck/Calbiochem) for three weeks. Individual colonies of resistant cells were picked, expanded, and validated for homogeneous expression levels and correct sub-cellular localization of fluorescent proteins. All cell lines used in this study had a normal morphology and cell cycle progression as compared to the maternal line. For details on the stable cell lines, see Supplementary Table 6.

Live microscopy

Automated microscopy with reflection-based laser auto focus was performed on a Molecular Devices ImageXpressMicro screening microscope equipped with 10x 0.5

Image analysis

Cell nuclei were detected by local adaptive thresholding [START_REF] Neumann | High-throughput RNAi screening by time-lapse imaging of live human cells[END_REF] , which is robust towards variable expression levels of the fluorescent chromatin marker in individual cells, and inhomogeneous illumination typical for wide-field microscopy. To improve segmentation accuracy, we implemented a split-and-merge approach. First, we split objects containing directly adjacent nuclei, using watershed transformation based on object contours. In some cases, this incorrectly split single objects. Thus we implemented object merging based on a priori definition of size and circularity criteria [START_REF] Wahlby | Combining intensity, edge and shape information for 2D and 3D segmentation of cell nuclei in tissue[END_REF] . Regions of interest for the secondary marker were derived by region growing of the chromatin segmentation to a fixed size, but constrained by regions of neighboring cells. Depending on the marker, we defined nuclear, cytoplasmic, or total cellular areas. This segmentation strategy turned out to be more precise than direct segmentation in the secondary channel, as many secondary markers dramatically changed in intensity levels or pattern throughout the time course of the experiment.

Texture and shape features [START_REF] Walker | Statistical geometric features-extensions for cytological textureanalysis[END_REF][START_REF] Haralick | Dinstein & Shanmugam Textural features for image classification[END_REF] (see Supplementary Table 1) were extracted from the two channels and all regions individually. For secondary region classification, only texture features were used since the shape information only depended on the chromatin segmentation.

Samples for morphology classes were manually annotated on the original images overlaid with the segmentation contours, to establish a training set for supervised classification. Support vector classification with radial-based kernel and probability estimates [START_REF] Wu | Probability estimates for multi-class classification by pairwise coupling[END_REF] was then computed with libSVM. Classification performance was calculated with five-fold cross-validation. Samples and feature plots for all classifiers used in this study can be accessed online through a web browser interface (see resource section).

Tracking cells over time was achieved by a constrained nearest-neighbor approach based on the Euclidian distance between objects [START_REF] Chen | Automated segmentation, classification, and tracking of cancer cell nuclei in time-lapse microscopy[END_REF] . Since tracks might be lost due to segmentation errors or migration of cells into the field of view the tracking must be able to create new tracks for all objects without incoming edges. To detect cell division events, or potential cell-to-cell fusion events, the tracking algorithm needed to support both splitting and merging. This yielded a hierarchical directed graph of isolated tracks for each cell over time. Tracking errors resulted mostly from segmentation errors and lead to wrong edges between the cell tracks. Secondary objects are tracked indirectly by the primary objects associated with them. Mitotic motifs were detected in this graph structure by the transition from prophase to prometaphase. Sub-graphs (mitotic trajectories) were extracted by considering a predefined number of frames preceding and following this mitotic motif, resulting in synchronized mitotic trajectories of equal length, as displayed in the figures.

Hidden Markov model and statistical analysis

A hidden Markov model λ is defined as λ = (X, A, Y, B, π), where X is the set of hidden states, A is a matrix of transition probabilities from one state to another, Y is the set of observable variables per state, B is a matrix of observation probabilities storing the probability of observation k being produced from state j (also termed emission or observation probability), and π is a vector of probabilities of the initial state (first time point) in the trajectory.

The hidden states X are the true cellular stages expressed by the class labels (8 classes for fluorescent H2B, see Fig. 1b). The hidden Markov model is learned by maximum likelihood estimates from the aligned trajectories of estimated prediction probabilities of the support vector machine, which is a three-dimensional array over trajectories, time points, and classes. Transition probabilities A are learned from the prediction probabilities along the trajectories on the underlying graph structure. In a free model all transitions between morphology classes were allowed (Fig. 2c). In a constrained model some transitions were suppressed based on biological a priori knowledge (transition probabilities were set to 0 for edges missing in the graph; Supplementary Fig. 6a). For the initial probabilities π the prediction probabilities of all trajectories at the first time point are considered. The observables Y are the class labels. The observation probabilities were either set to an error rate of 0.1%, or derived from the confusion matrix of support vector machine training.

Using the Viterbi algorithm, each trajectory was corrected based on its sequence of support vector machine probability estimates and the trained hidden Markov model for a given experimental condition (decode problem). This correction scheme was calculated individually for each marker and experimental perturbation condition.

To detect the onset of nuclear envelope breakdown and nuclear envelope reformation the time series of IBB-EGFP intensity ratios of individual cells were analyzed. We computed the ratio by a shrunken area of the chromatin object and a ring around. The onset was defined as the time point where the ratio was 1.5 fold increased above the ratio at the time point of chromosome segregation.

For data normalization of Fig. 5b we computed the z-scores of mitotic exit timing for all siRNA conditions (mean over all values of one condition). The z-score was computed by the mean of negative controls and the standard deviation of the entire data set.

Implementation and performance

The basic image processing was implemented in C++ using VIGRA (http://hci.iwr.uni-heidelberg.de/vigra) and in house-developed extensions. The C++ code was then wrapped for Python, which is a programming language particularly well suited for handling complex data structures and integration of external modules.

Statistical analysis and plots were performed with the R-project (http://www.rproject.org). The entire software package is platform-independent, and was compiled for Mac OS X and Windows environments.

Computation of each movie required 4-20 s per image and processor node, consuming 500-1500 MB RAM, depending of the number of frames and objects per frame. As an example, a single movie of Fig. 2 with 206 frames and ~37,000 objects required a total processing time of 34 min on a single processor node. For high-throughput analysis, we implemented distributed computing on a farm of desktop computers (four MacPro 2.2GHz, 28 cores total).

Software and data resources

CecogAnalyzer is a platform-independent graphical user interface, which covers the entire workflow presented in this paper. The software is publicly available in source and binary versions and was tested on MacOS X Leopard/SnowLeopard and Windows XP/7. We use a subversion repository for concurrent software development by remote contributors, and tracking of software changes. Our website is based on the project management tool TRAC (http://trac.edgewall.org/), which allows coordination of this open-source project by milestones, tickets, wiki pages and browsing of code changes.

The software, a subset of raw images presented here, the classifiers and parameters used for generating the figures are available online at http://www.cellcognition.org.

The classifiers data sets consisting of annotated samples and extracted features are interactively visualized by Adobe Flex and can be browsed online at http://flex.cellcognition.org. Normalization was per trajectory to the first prometaphase frame. Scale bar: 10 µm. Each gene was targeted by three different siRNA oligos (For full list of oligos, see Supplementary Table 1).

Final

  siRNA concentrations were 50 nM for Oligofectamine or 10 nM for HiPerfect. Cdc20 siRNA validation oligos were obtained from Qiagen with the following target sequences: AACCTTGTGGATTGGAGTTCT (Cdc20_1), CACCACCATGATGTTCGGGTA (Cdc20_2). Total HeLa cell lysates for SDS/Page analysis were prepared according to standard procedures. Rabbit-anti-human Cdc20 antibody (diluted 1:5000) was from Bethyl laboratories.

Figure 1 .

 1 Figure 1. Supervised machine learning and classification of morphologies. (a) Confocal images of live HeLa cells stably expressing a chromatin marker (H2B-mCherry), together with GalT-EGFP to visualize the Golgi apparatus, with mEGFPα-tubulin, or with the replication factory marker EGFP-PCNA. The images show maximum intensity projections of five z-sections. (b) Live imaging of HeLa cells expressing H2B-mCherry at different cell cycle stages, or apoptosis (2D time series imaged with wide field epifluorescence 20x dry objective, see Supplementary Movie 1). The color scheme relates to H2B-mCherry morphology classifications of subsequent figures. (c) Object detection (contours) and classification (colors) of cellular morphologies corresponding to predefined mitotic stages as shown in (b). Cells were tracked over time (arrows). See Supplementary Movie 2. (d) Classification performance of support vector machines with radial basis functions. The confusion matrix displays the matching of human versus machine annotation, identical annotations are on the diagonal. (e) Automated annotation of cell trajectories over time by the workflow shown in (c). 80 randomly selected trajectories (rows) over 40 time frames (columns) are displayed (time-lapse: 4.6 min). Colors refer to morphology classes as labeled in (b). Tick marks indicate sampled time points. Mitotic events are rare, and the trajectories contain many single frames of mitotic annotations, likely due to classification errors. Scale bars: 10 µm.

Figure 2 .

 2 Figure 2. Hidden Markov modeling of progression through morphology stages. (a) Automated extraction of mitotic events. Cells were synchronized in silico to the prophase -prometaphase transition. The plot displays a random selection of 100 mitotic events (from a total set of 172 mitotic events out of 8 movies; time-lapse: 4.6 min; see Supplementary Movie 1). Predicted morphology classes were color-labeled as in Fig. 1b. Asterisks: classification errors. Black frame indicates region of interest displayed by contour overlays on image data. For complete data, see Supplementary Movie 3. (b) Single cell and corresponding trajectory of class labels. Asterisks: classification errors. (c) Graph for all possible transitions between classes. Node 0 is start node, all other nodes are color-labeled as in Fig. 1b. (d) Learned class transition

Figure 4 .

 4 Figure 4. Timing phenotypes and kinetic measurements. (a) Mitotic progression assayed by H2B-mCherry morphology, and degradation of Securin-mEGFP. Examples are shown for untreated control cell (larger region of original data shown in Supplementary Movie 13), a cell with Mad2 RNAi-inactivated spindle checkpoint (siMad2; larger region of original data shown in Supplementary Movie 12), and a cell arrested in prometaphase by a Nocodazol (Noc; larger region of original data shown in Supplementary Movie 11). Time-lapse: 2.7 min. (b) Automated classification of mitotic stage progression as in Fig. 2f for the three experimental conditions shown in (a). (c) Securin-mEGFP degradation kinetics for the same cells shown in (b).

Figure 5 .

 5 Figure 5. RNAi screen for mitotic exit regulators. (a) Assay for mitotic exit timing. Live imaging of a cell line expressing H2B-mCherry and IBB-EGFP. The timing from anaphase onset (red bar) until onset of nuclear accumulation of IBB-EGFP (green bar) was used to define mitotic exit timing (arrow). Time is in min:s. Larger region of original data shown in Supplementary Movie 14. (b) Mitotic exit timing in an RNAi screen for 300 different RNAi conditions. 108 movies of different siRNA transfections were recorded in parallel over 46 h, to collect the entire dataset in four experiments. Time-lapse: 3.7 min; see Supplementary Movie 14. Each point in the graph indicates the z-score for one siRNA (for calculation of z-scores, see methods). Dashed lines indicate z-score threshold, solid line indicates mean of the entire dataset.

  (c) Cumulative percentage of cells exiting mitosis after onset of chromosome segregation (t = 0 min). The curves represent all mitotic events from two experimental replica. Cells were transfected in liquid phase with two different siRNA targeting Cdc20, or a non-targeting oligo for control, as indicated in the legend. (d) Confocal time-lapse imaging of a cell stably expressing H2B-mCherry and mEGFP-α-tubulin. Time is in min:s, maximum intensity projection of five zslices. See Supplementary Movie 15. (e) Confocal imaging as in (a) for a Cdc20 RNAi cell. See Supplementary Movie 16. Scale bars: 10 µm.

  Figure 1

  Figure 2

  Figure 3

  Figure 4

  94

  Effect of support vector machine class prediction noise on hidden Markov model error correction(a) Support vector machine prediction probabilities (8-vector; color corresponds to class label with highest probability) of Fig.2awere randomly sampled and replaced by random numbers (uniform distribution; sum 1.0). Fractions of the 100x40 objects were randomized from 0% to 100% in 10% steps and 8 repetitions. The overall prediction accuracy (correct predictions / total predictions) was measured based on the manual annotation. Trajectory panels are shown for 0%, 10%, 20%, and 50% randomized data w/o correction (upper panel) and with HMM correction (lower panel). (b) A hidden Markov model was trained on the partially randomized data shown in (a) and errors corrected on the randomized trajectories. (c) Plot visualizing the overall prediction accuracy w/o correction (red line) and with hidden Markov model correction (green line) relative to the percentage of randomized data. Average and standard deviation of 8 repetitions shown.

  Every n th frame (from 1 to 6) of the trajectories from Fig. 2a was selected to simulate different time-lapse intervals (4.6 min to 27.6 min). The overall accuracy was measured by comparison with manual annotation. (b) Same data and procedure as in (a) but with hidden Markov model trained on and applied to the subsampled trajectories.(c) Plot visualizing the overall accuracy w/o correction (red line) and with hidden Markov model correction (green line) for differently sampled time-lapse. Note: A higher time-lapse might decrease tracking accuracy, which is not reflected by this plot.

The MetaMorph journals developed for fast and robust acquisition of the time-lapse experiments presented here are available on our group website: http://www.bc.biol.ethz.ch/people/groups/gerlichd.

Acknowledgments

The authors thank Christian Conrad and Wolfram H. Gerlich for critical comments on the manuscript, F. Oliver Gathmann for helpful discussions about software engineering, Nicolas Graf for outstanding IT support, and Gabor Csucs, the ETHZ Light Microscopy Center (LMC), and the ETHZ RNAi Screening Center (RISC) for excellent technical support. We thank Jack Rohrer (University Zurich, Switzerland) for providing GalT-EGFP plasmid, and Jonathon Pines (Gurdon Institute, Cambridge, UK) for providing Securin-EYFP and cyclin B1-EGFP plasmids, and Katja Beck and Ulrike Kutay for providing images of cells expressing fluorescent α-tubulin and H2B, and Qing Zhong for generating the plot for Supplementary Fig. 1. Work in the Gerlich laboratory is supported by Swiss National Science Foundation (SNF) research grant 3100A0-114120, SNF ProDoc grant PDFMP3_124904, a European Young Investigator (EURYI) award of the European Science Foundation, an EMBO YIP fellowship, and a MBL Summer Research Fellowship to D.W.G., an ETH TH grant, a grant by the UBS foundation, a Roche Ph.D. fellowship to M.H.A.S, and a Mueller fellowship of the Molecular Life Sciences Ph.D. program Zurich to M.H. M.H. and M.H.A.S are fellows of the Zurich Ph.D. Program in Molecular Life Sciences. B.F. was supported by European Commission's seventh framework program project Cancer Pathways. Work in the Ellenberg laboratory is supported by a European Commission grant within the Mitocheck consortium (LSHG-CT-2004-503464). Work in the Peter laboratory is supported by the ETHZ, Oncosuisse, SystemsX.ch (LiverX) and the SNF.

invariant) h1_ASM h1_CON h1_COR h1_COV h1_DAV h1_ENT h1_IDM h1_PRO h1_SAV h1_SET h1_SHA h1_SVA h1_VAR h1_average h1_variance h2_ASM h2_CON h2_COR h2_COV h2_DAV h2_ENT h2_IDM h2_PRO h2_SAV h2_SET h2_SHA h2_SVA h2_VAR h2_average h2_variance h4_ASM h4_CON h4_COR h4_COV h4_DAV h4_ENT h4_IDM h4_PRO h4_SAV h4_SET h4_SHA h4_SVA h4_VAR h4_average h4_variance h8_ASM h8_CON h8_COR h8_COV h8_DAV h8_ENT h8_IDM h8_PRO h8_SAV h8_SET h8_SHA h8_SVA h8_VAR h8_average h8_variance Haralick 1 (32 gray levels, normalized to min/max gray values per object, distances 1,2,4 and 8 pixels rotation invariant) h1_2ASM h1_2CON h1_2COR h1_2COV h1_2DAV h1_2ENT h1_2IDM h1_2PRO h1_2SAV h1_2SET h1_2SHA h1_2SVA h1_2VAR h1_2average h1_2variance h2_2ASM h2_2CON h2_2COR h2_2COV h2_2DAV h2_2ENT h2_2IDM h2_2PRO h2_2SAV h2_2SET h2_2SHA h2_2SVA h2_2VAR h2_2average h2_2variance h4_2ASM h4_2CON h4_2COR h4_2COV h4_2DAV h4_2ENT h4_2IDM h4_2PRO h4_2SAV h4_2SET h4_2SHA h4_2SVA h4_2VAR h4_2average h4_2variance h8_2ASM h8_2CON h8_2COR h8_2COV h8_2DAV h8_2ENT h8_2IDM h8_2PRO h8_2SAV h8_2SET h8_2SHA h8_2SVA h8_2VAR h8_2average h8_2variance Statistical Geometric Features 2 ls0_CAREA_avg_value ls0_CAREA_max_value ls0_CAREA_sample_mean ls0_CAREA_sample_sd ls0_DISP_avg_value ls0_DISP_max_value ls0_DISP_sample_mean ls0_DISP_sample_sd ls0_INTERIA_avg_value ls0_INTERIA_max_value ls0_INTERIA_sample_mea n ls0_INTERIA_sample_sd ls0_IRGL_avg_value ls0_IRGL_max_value ls0_IRGL_sample_mean ls0_IRGL_sample_sd ls0_NCA_avg_value ls0_NCA_max_value ls0_NCA_sample_mean ls0_NCA_sample_sd ls0_TAREA_avg_value ls0_TAREA_max_value ls0_TAREA_sample_mean ls0_TAREA_sample_sd ls1_CAREA_avg_value ls1_CAREA_max_value ls1_CAREA_sample_mean ls1_CAREA_sample_sd ls1_DISP_avg_value ls1_DISP_max_value ls1_DISP_sample_mean ls1_DISP_sample_sd ls1_INTERIA_avg_value ls1_INTERIA_max_value ls1_INTERIA_sample_mea n ls1_INTERIA_sample_sd ls1_IRGL_avg_value ls1_IRGL_max_value ls1_IRGL_sample_mean ls1_IRGL_sample_sd ls1_NCA_avg_value ls1_NCA_max_value ls1_NCA_sample_mean ls1_NCA_sample_sd ls1_TAREA_avg_value ls1_TAREA_max_value ls1_TAREA_sample_mean ls1_TAREA_sample_sd Gray level features n2_avg n2_stddev n2_wavg n2_wdist n2_wiavg Gray level features (min/max normalized) n_avg n_stddev n_wavg n_wdist n_wiavg

Author contributions

Per class prediction performance compared with manual annotation of data without error correction (Fig. 2a) and with HMM error correction (Fig. 2f).

Supplementary Movie legends