Inferring statistical properties of 3D cell geometry from 2D slices.

Although cell shape can reflect the mechanical and biochemical properties of the cell and its environment, quantification of 3D cell shapes within 3D tissues remains difficult, typically requiring digital reconstruction from a stack of 2D images. We investigate a simple alternative technique to extract information about the 3D shapes of cells in a tissue; this technique connects the ensemble of 3D shapes in the tissue with the distribution of 2D shapes observed in independent 2D slices. Using cell vertex model geometries, we find that the distribution of 2D shapes allows clear determination of the mean value of a 3D shape index. We analyze the errors that may arise in practice in the estimation of the mean 3D shape index from 2D imagery and find that typically only a few dozen cells in 2D imagery are required to reduce uncertainty below 2%. Even though we developed the method for isotropic animal tissues, we demonstrate it on an anisotropic plant tissue. This framework could also be naturally extended to estimate additional 3D geometric features and quantify their uncertainty in other materials.

Introduction

Over the past decade, improved live-imaging techniques including multi-photon confocal [1] and light sheet microscopy [2] have dramatically altered our ability to quantify tissue architecture in in vivo and in vitro biological systems. In tandem, there has been an increased focus on developing mathematical models that can help organize and drive hypotheses about these complex systems.

Quite a bit of analysis and modeling has focused on confluent monolayers, where there are no gaps or overlaps between cells. These two-dimensional sheets of tissue are often studied in cell culture systems [3-5] and can also be found during embryonic development [6, 7]. Much of that work focuses on understanding how cellular properties (interfacial tensions, adhesion, adherens junctions) give rise to local cellular shapes and also how they help to generate the large-scale, emergent mechanical properties of tissue.

For example, researchers have developed a suite of mechanical inference techniques to estimate interfacial tensions and pressures from detailed images of cell shapes [6, 8, 9]. Others have quantified precisely the deformation mechanisms in the developing fruit fly using dynamical shape changes [10]. These methods rely heavily on automated watershed algorithms to segment membrane-labeled cell images in order to identify cell-cell interfaces in a network of many cells [11-16]. Existing segmentation algorithms have largely been optimized to work on two-dimensional cell sheets.

Another set of experiments and models has focused on the statistics of cell shapes as a metric to quantify global mechanical tissue properties. Specifically, studies of 2D cell vertex models (VMs) have found that cell shape may determine mechanical properties of confluent tissues (tissues with no gaps between cells) [17-19]. The models predict that when cells have a compact shape, so that their cross-sectional perimeter is small relative to their cross-sectional area, the tissue as a whole is solid-like in the sense that cells cannot migrate. In contrast, when cells have an elongated shape, so that their perimeter is large relative to their area, then the tissue is fluid-like in the sense that cells can easily exchange neighbors and migrate. The transition from solid-like to fluid-like behavior is predicted to occur at a specific value of the dimensionless 2D shape index, p2, which is defined as the ratio of the perimeter to the square root of the area. This prediction was shown to be precisely realized in human epithelial lung cell culture [5].

Given that many biological tissues are fully three dimensional, it is natural to wonder whether any of this work can be extended to 3D. From a modeling perspective, it is straightforward, although technically challenging, to develop 3D simulations. We have recently developed a 3D vertex-like model, called the 3D Voronoi model, and demonstrated that it, too, has a fluid-to-solid transition governed by cell shape [20]. In this case, the governing shape parameter is p3, which is the dimensionless ratio of the surface area of each cell S to its volume V: p3 = S/V2/3. It also appears fairly straightforward to generalize mechanical inference methods to 3D [21].

Although advances in imaging techniques have allowed much clearer and deeper imaging of 3D structures, it remains a technically difficult challenge to reconstruct the full network of cellular contacts in 3D [22-25]. For example, watershed algorithms for segmentation will fail if there is even one 2D slice where the membrane structure is poorly resolved, and so in general they have a very large error rate in 3D. In addition, many 3D structures of interest lie deep inside tissues where optical scattering makes live-imaging techniques difficult. In some cases, such as histological sections for staging of cancer tumors, only individual 2D images are available. Finally, a 3D reconstruction requires that all of the cells must remain sufficiently stationary while an image stack is acquired. Therefore, to our knowledge, very little of the exciting work in 2D can be robustly extended into live-imaged 3D experimental data.

This suggests that there may be an unexplored opportunity to use statistics of 2D images, which are standard in the field, to infer something about the statistics of 3D structures, an idea which has been exploited previously in materials science. Methods to estimate the grain size distribution within poly-crystalline materials have been proposed that use processed 2D imagery and assume 3D grain shapes [26-28]. Statistical reconstruction of 3D structure from 2D imagery has also been investigated for porous two-phase random media [29], particulate media [30], and media with shaped inclusions [31]. Typically, these methods start with a random 3D structure and have a process for evolving that structure to reduce differences between its 2D projections and 2D experimental data.

In our case, we would like to understand whether we can infer useful 3D shape information from 2D slices. Such an approach will not be directly helpful for mechanical inference methods, which rely on precise reconstructions of angles between junctions in 3D. However, it could prove very useful for testing predictions of vertex-like models where tissue mechanics is predicted to depend on cell shape, or perhaps for testing models for studying constrained cell migration through complex networks. Such migration can lead to DNA damage that depends sensitively on the shapes and sizes of pores in the constraining environment [32].

Therefore, the goal of this manuscript is to test whether information about 3D cell shapes can be reconstructed from randomly selected 2D image slices. Many experiments on mechanics and migration of cells in 3D focus on prepared tissues in collagen matrix or in centrifuged cell aggregates, and on other tissues, including organoids, certain tumors, and certain embryonic tissues, which appear isotropic and have relatively simple structure. We therefore perform this analysis in the context of a 3D Voronoi model [20], which is perhaps the simplest model for confluent cell bodies in 3D. In contrast to previous 3D inference methods, which typically sample only a fraction of the model phase space using Monte Carlo or dynamic minimization techniques [29-31], the simplicity of this model allows us to quantify the relationship between the generated 2D and 3D data across the entire relevant phase space.

We focus on determining whether the mean 3D shape index, which models suggest is strongly correlated with tissue mechanics, can be inferred from the 2D shape index, although we also explore other 2D and 3D shape descriptors as well. We find that there is a robust correlation between the 2D and 3D shape, and quantify the sensitivity of this correlation to sample size, experimentally relevant systematic errors, and tissue heterogeneity. We find that relatively few cells are required to converge on the correct mean 3D cell shape index, and that the estimates are quite robust with respect to moderate errors in 2D cell perimeter measurements, dropping cells with a small cross-sectional area, and cell size heterogeneity. This general framework may be extended to other contexts, ranging from more complicated or anisotropic tissues to extracellular matrix to disordered inorganic materials; all that is required is to substitute the 3D Voronoi model with reasonable 3D models of such systems or full 3D cell shape reconstructions and to replace p3 with the relevant 3D shape descriptors.

Results

2D cell shapes in slices of 3D cell packings

As described in the methods section, we use the 3D Voronoi model to create cell packings with a specified, homogeneous 3D shape index, given by p3 = S/V2/3. We cut the packings to obtain parallel slices of randomly oriented systems, yielding 2D networks of edges and vertices, as illustrated in Fig 1. For simplicity, we assume the slices to be extremely thin, although we can also vary the thickness of the optical section as described later. Our systems are small compared to experimentally-obtainable slices; for a large slice of an isotropic tissue, the results would be the same if one analyzed many different cells from the same slice. We then calculate the 2D shape index p2 for each cell in the 2D slices.

Fig 1

Visualization of the simulated tissue, confocal cross-sections, and individual 3D cell shapes.

(a) A 3D tissue geometry in which all 3D cells have the same shape index, p3 = 5.40, from the 3D Voronoi model. (b) The geometry seen in simulated ideal confocal imagery of that tissue. (c) Samples of the cell shapes from the model with specified p3 value. A 3D tissue geometry created using the 3D Voronoi model and the geometry seen in simulated ideal confocal imagery. In this example, all 3D cells have the same shape index (p3 = 5.45).

The distributions of p2 values for 20, 000 cells extracted from 500 3D packings are shown in Fig 2, for several different values of p3. Each curve demonstrates that even though each of the 3D cells has an identical 3D shape index, p2 exhibits a broad distribution of values. The lowest possible value of p2 is that of a circle, , although compact shapes in tessellated patterns are typically higher. Fig 1 shows a compact cell (labeled “a”) with a shape index of p2 = 3.75. In contrast, very elongated shapes occasionally arise when polyhedral edges are nearly parallel to the slicing plane, as shown by the cell labeled “b” in Fig 1. These very elongated shapes are a purely geometric effect that contain little information about the underlying 3D shape index, as shown in Fig 2(b). They contribute to a long tail in the distribution of p2 that strongly affects both the mean and variance. Therefore, we choose to focus on quantities such as the location of the peak of the distribution, , and the half-width at half max (HWHM) that are less sensitive to these tails.

Fig 2

2D cell shape index distributions.

(a) The distribution of 2D shape index, p2, in slices from 3D cell packings provide a signature of the 3D shape index, p3. p3 values vary from 5.4 (solid black) to 6.4 (dashed red) in increments of 0.2. (Data available online [33].) (Inset) The location of the distribution peak varies smoothly with p3. (b) The p2-distribution tail decays with a similar form for all values of p3. Example cross sections illustrate that arbitrarily large p2 values can be produced when the cross section is nearly parallel to a nearby edge—this produces nearly-rectangular slim shapes, for which p2 diverges when the cutting plane is near the edge. In contrast, general cross sections near a vertex are triangles, and the shape index, which is insensitive to the size of the triangle, retains a constant value independent of proximity of the cutting plane and vertex.

A first observation is that the peak in the 2D shape index () shifts dramatically with p3, suggesting that it should be possible to infer the 3D shape index from 2D data. This is quantified by the inset to Fig 2(a), which indicates that scales linearly with p3. More comprehensive measures of sensitivity such as the K-S test, discussed below, confirm that the p2-distribution sensitively reflects the p3 value.

While vertex models suggest p3 is an important parameter governing tissue rheology, this approach could also be useful for inferring quantities that are important for other models. Similarly, one might wonder whether there are better 2D shape descriptors than p2. In a later section, we study alternative geometry descriptors and find similar sensitivity, meaning that the inference can be made as easily. Therefore in the following sections we focus on the shape indices in 2D and 3D (p2 and p3), since they are simple to calculate and they control tissue rigidity in 2D and 3D cell vertex models.

Precision of p3 estimation depends on sample size

We ask how our estimate of p3 depends on the sample size (the number of p2 values extracted from cells in the 2D imagery, N). In experiments on 2D lung epithelia [5], the mean shape index of cells only varies about 4% in total. Therefore, the uncertainty in p3 must be smaller than a few percent.

It can be arduous to obtain the p2 values upon which the p3 estimate is based. In some biology experiments, cell boundaries are digitized by an analyst who hand-traces the outlines of labeled membrane proteins on a computer. Alternatively, automated computer algorithms may segment a 2D image into cell (and non-cell) regions, with trade-offs between speed and accuracy, to record the 2D shapes [12, 14, 15]. It is therefore important to establish how many cells must be analyzed to achieve the necessary level of accuracy in the estimate of p3.

To study how the best estimate of p3 and its uncertainty depends on the sample size, for each fixed value of p3 we segment a very large number of 2D cell images (30, 000 cells) and generate a fixed reference distribution for p2. We repeat this process for 20 values of p3 spaced equally between 5.4 and 6.4, to generate a library of reference distributions. These reference distributions are publicly available [33].

Next, for a fixed value of p3 we segment N randomly selected cells in 2D slices and generate a histogram of p2. This histogram can be compared to the reference library. The standard Kolmogorov-Smirnov (K-S) test (described in the section below), identifies the most likely distribution that produced the histogram, and correspondingly the estimated p3 value, . For this purpose we also created a publicly available online tool to compare data to the reference distributions [33]. We repeat this process 1000 times to measure the spread of estimates that occur. The fractional random error of the p3 estimate is , where is the standard deviation of . The fractional systematic error is . Fig 3 shows the random and systematic error as a function of number of cells traced. Note the difference of y-axis scale; Δp3 is usually insignificant compared to . After tracing 50 cells at random, the random error in the p3 estimate is less than 2%. The estimate converges to the actual p3 with 1% uncertainty after 100 random samples for p3 = 5.4 and after 200 random samples for p3 = 5.8. These results show that one may infer a sufficiently accurate estimate of p3 after imaging only a moderate number of cells.

Fig 3

Required number of measured cell shapes.

The random error, (upper), and the systematic error, Δp3 (lower), in the p3 estimate both fall rapidly with the number of cells traced, N. The systematic error is much smaller than the random error. The total error can be reduced to 1% by calculating the shape index of about 100 cells in this idealized case where the 2D measurements are exact. Curves are for p3 = 5.5 (black solid), 5.8 (blue dashed), 6.0 (red dotted), 6.1 (black dash-dotted), 6.2 (red long dashed).

Estimation of p3 is sensitive to systematic tracing errors

For the analysis in the previous section, we assumed that the 2D measurements of cell shape were exact, but data from experiments will have additional sources of noise and error. Therefore, we identify several likely sources and assess their impact on 3D shape estimation.

Cell shapes in 2D imagery can be measured by manual tracing on a computer or by automatic image segmentation and analysis. Some programs that measure cell shape report the perimeter length as the number of pixels that the cell perimeter passes through, and this artificially raises the length of a line segment by an amount that depends on the angle relative to the pixel axes. Attempting to infer the shape of isolated cells rather than a compact cluster of cells that tessellates space can also artificially inflate cell perimeters. Imprecise tracing of shape may result from other factors including inconsistent fluorescent dye saturation, uncertainty in identifying the cell borders, or limited image resolution, for example. We can model these sources of error by distorting the images before processing them and then estimating p3.

For certain types of noise, the measured distributions of p2 can be directly computed from those without noise. For example, a systematic overestimate of perimeter by a fraction ΔL has the effect of scaling the distributions of Fig 2(a) from P(p2) to P(p2/(1 + ΔL))/(1 + ΔL). A random mis-estimate of perimeter L by a fraction σ is modeled by multiplying the perimeter by a Gaussian random factor (with standard deviation Lσ). This convolves the distributions of Fig 2(a) with a Gaussian kernel.

Using this information, we can, as in the previous section, compute the errors in the estimated p3, plotted in Fig 4. Fig 4(a) shows that , the random error in the p3 estimate, is nearly independent of the perimeter overestimation ΔL and decreases with the number of traced cells, to less than 1% with N = 256 cells. Δp3, the systematic error, is approximately linear in ΔL and the error in the estimate can be reduced below 1% if the perimeter error is less than 3%.

Fig 4

Error propagation in estimate of p3.

The fractional random error, , and the systematic error, Δp3, in the p3 estimate increase with mis-measurement of cell perimeter, but are not significantly influenced by neglecting small cell areas. Curves are for N = 16 (solid), 64 (dashed), and 256 (dotted) cells, and p3 = 5.6 (thick lines) and p3 = 6.0 (thin lines). There is little sensitivity to p3, so results for other p3 values are not shown. (a,b) A systematic overestimate of cell perimeter by a fraction ΔL produces a nearly-proportional systematic error in p3 estimate. Tracing additional cells reduces uncertainty but does not reduce systematic error. (c,d) A fractional random error of standard deviation σ in the cell perimeter measurement also produces an error in the p3 estimate. (e,f) Accidentally neglecting cells with small area causes only a small error.

Fig 4(c) and 4(d) show the error accrued from random tracing errors that result in mis-estimate of cell perimeter by a fraction σ. For large values σ, the estimated p3 can differ systematically from the true p3, but if σ < 0.06 (less than 6% uncertainty of each L measurement) the systematic error is less than 1%, and tracing a large number of cells reduces the random error in .

Another possible source of systematic error is that the smallest shapes in microscopy may be accidentally overlooked and not traced. We generate reference distributions for this error by filtering the small-area shapes from the distribution functions. Neglecting shapes of area below a certain threshold, a, introduces relatively little error into the p3 estimate as shown in Fig 4. When cells with area 10% of the mean are accidentally neglected, the result is only an error of 1% in p3. To avoid errors when comparing against the distributions generated here, cell shapes of all sizes should be included. If only large cell shapes can be traced reliably in experiment, then the resulting histograms may be compared to reference distributions that are generated using only the large shapes.

The above results show that generally an error in a tracing measurement generates a similar order of magnitude fractional error in the p3 estimate. If multiple errors occur (e.g. both ΔL and σ are significant), and they are independent and uncorrelated, systematic errors Δp3 are approximately summed while random errors are approximately added in quadrature.

Effect of heterogeneity

So far, we considered homogeneous tissues where all cells have the same cell volume V and 3D shape descriptor p3. Of course, this is an idealization, as in real tissues there will be variations in these quantities. Now, we study the influence of such variations by generating 3D cell packings with Gaussian distributions of shapes p3 or volumes V.

We first focus on variations in cell shape. We again generate the cell packings for an additional 200,000 cells as energy-minimized states of the 3D Voronoi model (see Methods), where states are included only if the cells achieve the target mean and standard deviation of the 3D shape index; cases near extreme values of (near 5.4 or 6.4) targeting high- were not able to achieve the targets as seen in other simulations [19, 20] and so are excluded.

As discussed earlier, we focus on the location of the peak in the 2D shape index () and the half-width-at-half-max (HWHM) to minimize contributions of the universal tail. Fig 5 shows and HWHM for various values of and . While the peak and width of the distribution are strongly correlated with , these quantities are much less sensitive to . A sensitivity analysis based on the K-S test confirms this result for both these and other shape descriptors.

Fig 5

Properties of p2 distributions from heterogeneous tissues.

The peak and half-width-at-half-max (HWHM) of the p2-distribution accurately reflect the mean constituent 3D cell shape index, but are relatively insensitive to variation in shape. Curves vary for (black), 0.1 (blue, dashed), and 0.15 (red, dotted), and vary for (magenta, long dash), 5.65 (cyan, long dash-dot), and 5.8 (yellow, dash-dot).

It is natural to ask why and how the 2D shape distribution is mostly insensitive to the variance in 3D shape index. To gain an intuition for this observation, it is useful to look at a specific example. Fig 6(a) illustrates the 2D shape distribution for a binary mixture of approximately 50% p3 = 5.6 and 50% p3 = 6.2 cells. The distribution P from the mixed system (dashed green line) is distinct from P of the homogeneous systems, the most similar of which is p3 = 5.9 (solid black line). This emphasizes that the 2D shape distribution from a mixed system is distinct from any of the homogeneous-system distributions.

Fig 6

p2-distributions for heterogeneous 3D cell shape or volume.

(a) The p2-distribution for cell packings of homogeneous 3D shape index p3 = 5.9 (solid black) is compared to a packing of a binary mixture (dashed green). The binary mixture is 50% of p3 = 5.6 and 50% of p3 = 6.2 cells, so that the mean 3D shape index is 〈p3〉 = 5.9. The superposition (dotted violet) of the p2-distributions of packings with homogeneous p3 = 5.6 and p3 = 6.2, respectively, are a close approximation of the heterogeneous system. (b) The p2-distribution for cell packings of homogeneous 3D shape index p3 = 5.9 (solid black) and of heterogeneity in 3D shape index (dashed blue) and in cell volume σ = 0.15V0 (dotted red).

Next, we plotted the superposition of the homogeneous distributions, P ≈ 0.5P5.45 + 0.5P6.0 (Fig 6(a), violet dotted line) and found good agreement with the distribution of the mixed system. This suggests that the mixed system is composed of approximately the same 3D shapes that are generated in homogeneous simulations of the constituent p3 values. Under this approximation, the p2-distribution of the heterogeneous system is where f is the fraction of the cells observed in a cross section of the system that are in subpopulation β. The value of f is affected not only by the relative frequency of the cell subpopulation β in the 3D tissue, but also from the increased number of intersections of the imaging plane with large or elongated 3D cells than with small or compact 3D cells. It can be connected to the mean 2D cross-sectional area of cells of that subpopulation, 〈a〉, as well as to the fraction of the total cell volume in the tissue from subpopulation, . In terms of these quantities, it is where the sum occurs over each constituent cell subpopulation α.

Since the 2D shape distribution from heterogeneous systems is approximately given by the superposition of pure-p3 distributions, the sensitivity to p3-heterogeneity is linked to the variation with p3 itself. The lack of change of the distribution with heterogeneity simply arises from two facts. First, many similar constituent distributions underlie the distribution of the heterogeneous system, so for instance, Gaussian-distributed p3 heterogeneity with , 68% of the distribution weight comes from the narrow range . Second, the smooth variation of the distribution with p3 (e.g. Fig 2) means that Gaussian heterogeneity (which is symmetric above and below the mean p3) produces somewhat canceling differences to the homogeneous distribution.

We also show explicitly in an example that variations in volume have even less impact on the estimations of cell shape. This can be seen by comparing (Fig 6(b)) distributions from a homogeneous case, σ = 0 (solid black), a p3-heterogeneous case (dashed blue), and a V-heterogeneous case (dotted red). The mean V and p3 are the same in all cases, and either V is the same for all cells while p3 is drawn from a Gaussian distribution (standard deviation ) or p3 is the same for all cells and V is drawn from a Gaussian distribution (σ = 0.15μ).

The curve corresponding to heterogeneous volume is much more similar to the homogeneous case than the one corresponding to the similarly heterogeneous distribution of cell shapes. The K-S distance quantifies this; the distance between the dotted red curve and the black curve is , which is about 10% of the distance between the black curve and the dashed blue curve. Practically speaking, this means that the 2D shape data presented in this manuscript Figs 2, 5 and 7 provide accurate estimates of 3D cell shape even if there are fluctuations in preferred cell volume.

Fig 7

Distributions from alternative shape descriptors are similarly distinct.

The distribution of p2 values depends on the volume fraction of the largest inscribed sphere f, and distributions of the 2D anisotropy index m reflect the 3D shape descriptors p3 and f. Curves correspond to f value, varying in 0.1 increments from 0.15 (dashed, red) to 0.65 (solid, black) or p3 value, varying from 5.4 (solid black) to 6.4 (dashed red) in increments of 0.2.

Taken together, these results suggest that 2D shape analyses are very good at estimating the mean 3D shape index, even with heterogeneity, but are not very useful for estimating the variance of the 3D shape index. Although the p2-distributions vary with heterogeneity (and so the method is sensitive to heterogeneity in principle), the variation is relatively small. Future modeling work should therefore focus on understanding what features of heterogeneous systems are important for understanding global tissue mechanics.

Alternative geometry descriptors

The above analysis focuses on cell shape indices p2 and p3 as 2D and 3D shape descriptors, respectively. However, the method is quite general and other shape descriptors can also be considered. Here we demonstrate this point with two other quantities. First we show that the p2 distribution can be used to extract the mean value of f, the ratio of the max inscribed sphere volume to Voronoi cell volume, which provides another characterization of a 3D shape. To study the distribution of p2 values, we use the same set of simulations of 20,000 cells as before and bin cells by their f value. Fig 7 shows that the distribution of p2 varies with the value of f, which ranges between 0.1 and 0.7. The distributions in Fig 7 resemble those in Fig 2, because both p3 and f characterize compactness of the 3D shapes.

We additionally ask if a different 2D shape descriptor can be used to extract the mean value of p3. We consider distributions of m, the anisotropy index [34]. The shape anisotropy, m, is determined using the 2 × 2 moment of inertia tensor of the shape, Here, the integral runs over the area of the shape in the x-y plane. is the 2D vector from the shape center of mass. The anisotropy index, m, is the difference of the eigenvalues of the matrix G divided by their sum. Again, a confocal slice through a collection of cells of a specified p3 produces a distribution of these 2D shape descriptors.

The additional distributions of m, for varied values of p3 and f, are also shown in Fig 7 As m remains bounded within the unit interval, the distribution lacks a long tail and may be preferable in some cases to p2 as a 2D shape descriptor. The presence of the peak at low-m is reflective of the correlations between m and p2.

Multi-dimensional (joint) distributions , where a is cell 2D cross-sectional area, may also be used for more complete comparisons between the 3D model and the 2D geometry in slices. Additional 2D descriptors provide more sensitivity and points of comparison between experiment and model, whereas additional 3D descriptors provide more information about the 3D geometry.

Shape estimation from experimental data

Here, we demonstrate our method on experimental data. While we have developed our method for 3D isotropic disordered animal tissues, segmenting these tissues in order to precisely measure p3 remains a difficult task, and not surprisingly, we have indeed found no publicly available segmentation data of such tissues. In contrast, many plant cells are more readily imaged in 3D, but, conversely, they are typically ordered with a more regular structure. Here, we use segmented images of cells in the Arabidopsis plant, in data provided in Ref. [35], (labeled “Plant 1, Hour 0”), to illustrate the use of our method. Unlike the our 3D Voronoi model data set that our method is based on, these tissues are anisotropic. Despite these differences in geometry, we use the CVM-derived p2 distributions provided above to estimate 3D cell shape.

We sample the 3D-segmented dataset on a cubic grid of spacing 0.66 × 0.66 × 0.74 μm. The marching cubes algorithm [36, 37] identifies the 3D surface of each cell. The mean p3 is 5.5 and standard deviation is 0.2. Cells are then sorted according to p3 into bins of width 0.1. For each p3-bin, we slice all cells within the bin and find the p2 values of shapes in the 2D slices. The K-S test indicates the most similar distribution.

Fig 8 compares distributions of p2 from the plant cells with those from the CVM. Although the roots are not the isotropic aggregates represented by the model, the distributions for a given p3 are quite similar to the ones derived from the CVM. The estimate is near the actual value of p3. Over most of the range of p3 values, . At the low end, p3 = 5.2 cells are estimated to have , since that is the lowest value of p3 considered in the model geometries, producing a mis-estimate of 0.2. At the high end, p3 = 6.2, there are only four cells in the dataset analyzed, producing a noisy p2-distribution, and a best-fit estimate of . Additional experimental 3D cells with large p3 would be required to smooth the distribution. When slices are taken from the full 3D dataset without separating cells, the result is (K-S distance 0.07), a difference of 0.05 from the mean p3 = 5.50.

Fig 8

Analysis from experimental geometries.

The p2-distribution from slices of cell geometries extracted from experimental imagery (red) are compared with the reference distributions generated from the 3D Voronoi model (gray). The closest reference distribution (blue), identified by the K-S test, gives the p3 estimate.

Notably, the plant cells are not too dissimilar from the CVM geometry, approximately having convex polyhedral shapes. Tissues with the same geometry as the CVM may compare directly to its p2 distributions. Tissues with dramatically different cell shapes may expect large K-S distances from all CVM distributions, and would require distributions from an appropriate model. Because the plant cells are reasonably well-described by the CVM, the 3D shape parameters from 2D slice data that we estimate using this method are in excellent agreement with the actual 3D shape parameters obtained from the 3D cell reconstructions.

Connection to volume or surface area determination

This analysis has focused on determination of the average dimensionless cell shape index in confluent tissues. In certain applications, it may be desired to determine properties of cell volume or surface area independently.

Techniques to estimate average cell volume can be comparatively simple, and there are several possibilities. First, from a 3D image, average cell volume for confluent tissues is simply the imaged volume divided by the number of cells. This requires counting cells (e.g. by counting the nuclei for mononucleated cells), which is much simpler than a 3D segmentation of cell membrane boundaries, as needed for cell shape.

Second, it is tempting to infer mean cell volume from dimension-full 2D cell geometry such as areas a or perimeters, in line with the method presented in this paper. Dimensional analysis suggests, 〈V〉 = c1〈a〉3/2, where 〈a〉 may can come from a simple count of confluent cells, and c1 still depends on cell shape. In the homogeneous 3D Voronoi geometry, we find that c1 varies between 1.52 and 1.85 and depends on the 3D cell shape parameter p3 approximately as c1 ≈ 0.32 ⋅ p3 − 0.2. We find no significant change in this relation when we introduce heterogeneity using a Gaussian-distributed p3 or V up to standard deviation σ < 0.3V0 and .

If cell volume can be determined, cell surface area can be approximately recovered from p3, the dimensionless surface area. With heterogeneous shapes, because the estimate follows the mean p3, the mean surface area is —where this result is exact for homogeneous-volume cells.

Discussion

We have investigated a method for determining statistical information about the 3D shapes of cells using only 2D slices through a tissue. We focused on the 3D shape index p3 and found that it can be inferred directly and reliably from the distribution of measured p2 values. We quantified the number of p2 measurements (and therefore cell traces) required to reduce the random and systematic errors in p3, showing that typically only of order 100 cell traces are required to reduce uncertainty in p3 below a given threshold. The method is reasonably robust against the modeled sources of error and we suggest that the peak and width of the p2 distribution are easily-accessible quantities that biologists can compare directly to our reference distributions [33]. We also study the effect of tissue heterogeneity, and find that 2D shapes can be used to estimate the average 3D shape, but provide little information about the variance in 3D shape.

It is possible to envision many extensions of this general framework. For example, confocal images occasionally have rather limited resolution along the z-axis, so that the effective thickness of the 2D slice is significant compared to the diameter of the cell. If the slice thickness is well-characterized, it is possible to use ray-tracing as discussed in the methods section to generate reference distributions relating 2D shapes to 3D shapes corresponding to a specific value of the slice thickness. Another extension would be to study different models for cellular structure, such as vertex models [6, 7, 38] where the degrees of freedom are cell junctions instead of the cell centers, cellular Potts models [39, 40] where each cell is composed of a grid of points, or even models for non-biological materials such as foams [41].

We note that we are assuming that the tissue is isotropic on average. For such tissues, the observed distribution of 2D shapes does not depend on the orientation of the 2D slice. Conversely, for anisotropic tissues, a dependence on orientation will generally be observed, which could be used to infer the degree of alignment. We therefore presume that it is not an experimental challenge to determine whether a tissue is anisotropic. If it is not possible to study the dependence on orientation, then tissue anisotropy will show up in the 2D shape distributions, which will typically differ from those presented here.

We further note that while it has been found that Voronoi tessellations approximate 2D tissue geometry reasonably well for many purposes [42], the degree to which the 3D Voronoi model mimics a specific 3D biological tissue can be expected to vary with the type of tissue. On the other hand, because the high sensitivity of the p2 distribution to p3 originates from the fact that more compact cells produce more compact cross-sections, we expect the p3 estimate to be robust to minor variations in the geometry. Indeed, this is seen in the demonstration on experimental cell geometry for Arabidopsis, even though the tissue is anisotropic and is presumably not well-described by the cell vertex model.

It is also possible to envision extensions of this work beyond confluent cellular structures. For example, it may be possible to perform a similar analysis on particulate models and compare to nuclei-labeled images. Alternatively, cell migration in fiber networks is conjectured to be limited by the rate at which cells can squeeze their nuclei through the pores in the mesh [43]. It would be interesting to see if one could estimate typical pore sizes from statistics of 2D slices of fiber network models.

Although great progress is being made in developing techniques to fully reconstruct 3D cell shapes [22-25], their use is restricted to specimens that are optically transparent, or are compatible with light-sheet microscopy, and they can not be used in situations where cells exchange neighbors faster than a 3D scan can be completed. In contrast, 2D images of cell structure are ubiquitous in medicine and biology, from histological sections of cancer tumors for use by pathologists to standard brightfield microscopy techniques, and they are also fast to obtain.

Our results suggest that the techniques described here could be easily used by biologists and clinicians to obtain information about 3D shapes in many cases where full 3D reconstruction is not possible or extremely time-consuming. We have shown that our methods are quite robust, and that the distributions developed here work quite well even on cell shapes in an Arabidopsis tissue. In addition, this paper highlights a different way of thinking about how to use 2D cell images. Currently, cell shape and shape polarization are often described in terms of simple 2D quantities such as the axes of the best-fit ellipse or the length and width, and this naturally leads to a focus on how those quantities correlate with motion in a 2D plane. For cells in 3D tissue, our work suggests that 2D quantities provide information about 3D shape, which could be used to drive and test hypotheses about cell migration and tissue mechanics in 3D.

Methods

To investigate the connection between 3D cell shape and 2D cross-sectional cell shapes in a 3D geometry that is representative of simple isotropic 3D tissues, we use the recently developed 3D Voronoi model [20]. Based on a Voronoi tessellation of cell center positions [44], we divide a simulation volume with periodic boundaries into polyhedra, representing cells. Cell centers move to minimize a Hamiltonian, in which each cell i has a target volume V0 and surface area S0, and N is the number of cells. We choose the target cell parameters either to be equal for all cells (S0 = S0 and V0 = V0) to simulate a homogeneous system, or such that target volumes V0 and target shape indices are independently drawn from Gaussian distributions. The parameters k and k are stiffness constants which we set to 1 here without loss of generality, since we will focus on cases in the fluid regime of the model [20] where all cells attain their target parameters and E ≈ 0. The simulation volume is fixed at V = ∑ V0. Starting from uniformly randomly distributed cell centers, we use the FIRE minimizer [45] to minimize the energy with respect to all cell positions, generating an isotropic ensemble of 3D cell shapes.

For homogeneous target shapes, we observe that when , all cells satisfy their optimal values, S = S0 and V = V0. This allows generation of a large ensemble of disordered cells all with the same shape index, p3 = S/V2/3. When , the geometry becomes pinned [20] such that the mean shape index is 〈p3〉 ≈ 5.4, because it appears that disordered packings with smaller shape index do not exist [20]. For values of p0 > 6.4, multifold vertices become common and minimization algorithms face challenges [19] and so we restrict simulations to p0 ≤ 6.4. We created packings with shape indices p3 from 5.4 to 6.4 with 0.05 increments.

Once 3D packings with a defined distribution of p3 are generated, we simulate the acquisition of 2D cross sections by generating images of intersections of cells with a specified plane (Fig 1) using the software POV-Ray [46]. Based on a segmentation of these cellular cross sections, we quantify the cell 2D shape index p2 as the quotient of perimeter divided by the square root of the cross-sectional area, which provides us with a histogram of 2D shape indices p2 for the given ensemble of cell packings with its predefined p3 distribution.

To compare different p2 distributions with each other, we use two kinds of measures. As a practical measure for experimentalists to extract the 3D shape index from a p2 distribution, we propose to use peak and HWHM of the p2 distribution (Figs 1 and 5). As a complementary measure to compare p2 distributions in more detail, we also use the K-S test, which measures the maximum distance between the two respective cumulative distributions [47] (used for Figs 3–5 and 8).

K-S test for goodness of fit and sensitivity

The Kolmogorov-Smirnov (K-S) test characterizes the likelihood of a set of measurements, given a probability distribution for the measured quantity. [47] For our purposes, K-S assigns a distance, D, to the dissimilarity of a measured histogram and one of the reference distribution, P(p2). Specifically, D is computed from the maximum separation between the empirical cumulative distribution function (EDF) of the data and the cumulative distribution function (CDF) of the reference distribution. The CDF of a distribution P(X) is . The corresponding quantity for the set of p2 values, {X}, that go into the histogram is the EDF, . N is the number of samples in the histogram, and Θ(x − X) is the Heaviside step function, equal to 0.0 for x less than X and otherwise equal to 1.0. The K-S distance D is the maximum separation between the two, . This distance allows identification of the model distribution that best fits the sampled data in the main text.

The K-S test quantifies the goodness of fit between data and model in the presence of random noise [47]. This is done by comparing the K-S statistic with a published table and indicates in our context whether the 2D geometry from the experiment is consistent with the 3D geometry of the model considered.

Furthermore, the K-S distance D can be used to quantify the sensitivity of a distribution to p3 by characterizing how quickly a p2-distribution changes upon changing p3. For a small change of p3 by δ, we find that the K-S distance D increases from zero linearly. The sensitivity is given by dD/dδ. The p2-distributions of Fig 2 are most sensitive to p3 at p3 ≈ 5.4 where dD/dδ ≈ 1.5, while near 6.4 the distributions are less than half as sensitive with dD/dδ ≈ 0.5. Thus, while the location of the distribution peak varies linearly with p3 with slope 0.4, the changes in distribution shape produce additional sensitivity near p3 = 5.4.

Similarly, the sensitivity of the m-distribution to p3 is found to be at p3 = 5.4 and falls to by p3 = 6.4. Thus distributions of m and p2 are about equally sensitive to p3, and the uncertainties in Figs 3 and 4 in the main text would be similar for m. Here is the K-S distance between m-distributions from p3 values that differ by δ. Analogously, being the K-S distance between p2-distributions from f values that differ by , the sensitivity to f is found to be and near f = 0.2, rising to and at f = 0.5. This quantifies the level of sensitivity that is evident in Fig 7.