Characterizing heterogeneous cellular responses to perturbations

Abstract

Cellular populations have been widely observed to respond heterogeneously to perturbation. However, interpreting the observed heterogeneity is an extremely challenging problem because of the complexity of possible cellular phenotypes, the large dimension of potential perturbations, and the lack of methods for separating meaningful biological information from noise. Here, we develop an image-based approach to characterize cellular phenotypes based on patterns of signaling marker colocalization. Heterogeneous cellular populations are characterized as mixtures of phenotypically distinct subpopulations, and responses to perturbations are summarized succinctly as probabilistic redistributions of these mixtures. We apply our method to characterize the heterogeneous responses of cancer cells to a panel of drugs. We find that cells treated with drugs of (dis-)similar mechanism exhibit (dis-)similar patterns of heterogeneity. Despite the observed phenotypic diversity of cells observed within our data, low-complexity models of heterogeneity were sufficient to distinguish most classes of drug mechanism. Our approach offers a computational framework for assessing the complexity of cellular heterogeneity, investigating the degree to which perturbations induce redistributions of a limited, but nontrivial, repertoire of underlying states and revealing functional significance contained within distinct patterns of heterogeneous responses.

Fig. 1.

Heterogeneous cellular populations can be characterized as mixtures of subpopulations with distinct phenotypes. For convenience of visualization, results are shown for k = 4 subpopulations. (A and B) Illustration of approach. (A) HeLa cells are treated with drugs, labeled with fluorescence markers, and imaged. Individual cellular regions are identified. (B) A high-dimensional vector of phenotypic features is extracted from each cell. Cells are assigned probabilities of belonging to subpopulations that have been identified by using a GMM. Scatter plot: cells are visualized as points via feature representation and PCA reduction to two dimensions; colored cell outlines colors in A and points in B correspond to most probable subpopulation assignment; ellipses represent unit standard deviation from the mean for each subpopulation. Shown are 3 representative cells taken from various locations in each ellipse. (C) Shown are subpopulation models for different marker sets. Models are built from “reference” sets of cells sampled from multiple conditions. Scatter plots are as in B. Stacked bar graphs represent (prior) probabilities of cell membership in each subpopulation. Sample cells from within 1 standard deviation of the mean phenotype are shown at right for each subpopulation (image colors correspond to markers in order; e.g., blue/green/red corresponds to DNA/SC35/anillin). (Scale bar: 20 μm.)

Fig. 2.

Subpopulation profiles allow quantitative comparison of heterogeneous cellular responses. (A) Subpopulation profiles reveal dose-dependent heterogeneous responses to drugs [DNA/SC35/anillin marker set, 20-h drug exposure (25)]. Thirteen concentrations are shown, increasing from left to right, at 3-fold titration. (B) Subpopulation profiles reveal time-dependent heterogeneous responses to drugs. Subpopulation profiles for selected drugs are combined in time series (DNA/p-p38/p-ERK marker set; concentrations are given in Table S1). Profiles in A and B are based on reference models computed with k = 4 subpopulations. Colors and subpopulations ids S1-S4 on left-most images correspond to Fig. 1C. (C) Drugs of similar mechanism often induce similar subpopulation profiles. Subpopulation profiles for 25 drugs at the 24-h time point and DNA/p-p38/p-ERK marker set are displayed in a heat map and grouped by category (D, DNA replication; H, histone deacetylase; M, microtubule; G, glucocorticoid receptor; T, topoisomerase). Drugs with <200 cells per well are marked with asterisks (mitomycin C and trifluoperazine).

Fig. 3.

A small number of phenotypic subpopulations can distinguish drug effects by mechanism. (A) As the number of subpopulations in the reference model increases (k = 2, 4, 10, and 100), MDS plots (40) reveal relative similarities among drugs (computed by using KL divergence) within 3 mechanism categories (H, blue diamonds; M, green triangles; G, orange stars). Axes are labeled in arbitrary units. (B) Drug classification performance is computed for each value of k (and drug category; see also Table S5). Performance is measured by −log₁₀ of the P value for a rank-sum statistic [Mann–Whitney U test (39)] obtained by comparing distributions of intra- and intercategory subpopulation profile similarity scores. Profile similarity is computed for 12- and 24-h time points of low drug concentration and combined by rank addition. Dotted lines indicate significance threshold corresponding to P = 0.05. Black arrowheads (Left) (k = 2, 4, 10, 100) correspond to MDS plots in A. Black arrowheads (Right) (k = 7) correspond to MDS plot at bottom right and uses reference model shown in Fig. S1D.

Fig. 4.

Exceptional phenotypes are detected by evaluating the fit of each drug to a reference model built by using only the low concentration of 1 drug in each functional category (control, aldosterone, aphidicolin, apicidin, colchicine, and doxorubicin) at all time points. (A) Reference models for DNA/p-p38/p-ERK and DNA/actin/α-tubulin marker sets are both computed with k = 4 subpopulations. Heat map reflects the ratio of model fit (log likelihood) of drug treatment data to training data. Four time points (3, 6, 12, and 24 h) are displayed for each of 2 marker sets, at both low and high concentrations of each drug. Gray boxes indicate poor fits due to low cell counts (<700 cells per well). (B) Shown are sample images of drug treated cells that fit the reference model well (a–c), as well as the 6 lowest scoring drugs (–6) as indicated in A.