Spatial protein analysis in developing tissues: a sampling-based image processing approach.

Advances in fluorescence microscopy approaches have made it relatively easy to generate multi-dimensional image volumes and have highlighted the need for flexible image analysis tools for the extraction of quantitative information from such data. Here we demonstrate that by focusing on simplified feature-based nuclear segmentation and probabilistic cytoplasmic detection we can create a tool that is able to extract geometry-based information from diverse mammalian tissue images. Our open-source image analysis platform, called 'SilentMark', can cope with three-dimensional noisy images and with crowded fields of cells to quantify signal intensity in different cellular compartments. Additionally, it provides tissue geometry related information, which allows one to quantify protein distribution with respect to marked regions of interest. The lightweight SilentMark algorithms have the advantage of not requiring multiple processors, graphics cards or training datasets and can be run even with just several hundred megabytes of memory. This makes it possible to use the method as a Web application, effectively eliminating setup hurdles and compatibility issues with operating systems. We test this platform on mouse pre-implantation embryos, embryonic stem cell-derived embryoid bodies and mouse embryonic heart, and relate protein localization to tissue geometry. This article is part of a discussion meeting issue 'Contemporary morphogenesis'.

Introduction

One of the main areas of focus in developmental biology is how tissue complexity is generated and how diverse cell types are consistently organized into a complete organism. A common challenge in the field is the three-dimensional spatial analysis of protein distribution within tissues and cells [1,2]. The function of a protein can vary depending on the cellular compartment it is localized to and differential expression of proteins across tissues during development is important in pattering those tissues. Cutting edge optical microscopy and image analysis tools are now routinely used for spatial analysis of developing tissues from various organisms [3-9]. Given the variability of different tissues, existing software often must be optimized to work for specific cases. For example, full cell segmentation tools like MorphoGraphX [3] or Packing Analyzer [10] work extremely well with samples such as plants and Drosophila, where high-quality fluorescence signals and low noise make it possible to segment complete cell shapes.

As the signal-to-noise ratio of the image drops or features become indistinct, it becomes more difficult to accurately segment and analyse individual cells. One of the most reliable solutions that is still widely used is manual segmentation [11,12], because humans are adept at feature recognition. However, manual outlining is time consuming, which is particularly a problem for volume or time-lapse data. To automate image segmentation, there has been a shift towards machine learning segmentation approaches. A notable example is RACE [9], which uses machine learning for tissue segmentation into individual cells. Arguably, this is one of the most successful approaches yet for automated segmentation, especially from membrane signals. However, it is still often impossible to correctly segment cells based on their three-dimensional outlines alone, because of limitations in image quality.

Dealing with multiple spatial and fluorescence variables raises additional challenges in terms of signal normalization, which is crucial in correcting for variations due to image acquisition or optical heterogeneity within tissues. A commonly used approach involves taking ratios of nearby cellular compartments, such as nuclei and cytoplasm. This has been previously used for example in investigating the Hippo pathway, where YAP localization could be quantified by taking the ratio of nuclear and cytoplasmic fluorescence while at the same time normalizing the raw signals [13].

Even in cases where image quality prevents reliable segmentation, image data often contain a great deal of useful quantitative information. To extract such information, we have developed a sampling-based approach to measure nuclear, cytoplasmic and plasma membrane fluorescence levels. Importantly, the localization of these samples within the image volume is retained, which can then be used for detailed statistical analysis of the spatial distribution of the detected proteins. This approach works across the scale of individual cells to entire tissues. At the level of the cell, it allows us to quantify relative protein levels in different cellular compartments, while at the level of the tissue it allows us to identify cell populations based on relative protein expression, statistically test the difference between these populations and analyse how their location within the tissue correlates with protein levels.

Our approach employs well-defined features such as nuclei to perform partial but robust segmentation using well-established techniques [6,14]. Alternatively, well-resolved membrane signals or manually defined regions can also be used. If nuclei are used as points of reference, nearby voxels are used to sample nearby nuclear and cytoplasmic regions, identified based on nuclear stain intensity. Simultaneous sampling of different cell compartments allows one to normalize the data with respect to experimental variation and tissue optical heterogeneity. Regions of membrane sampling are located by using a combination of detected nuclei positions and intensity of a membrane stain. Membrane signal quantification can be further improved by using commercially available software such as Imaris to manually outline membrane regions, which are then automatically partitioned by our software into contacting (basolateral), exposed (apical) and junctional domains. Since the spatial distribution of proteins is of particular interest to developmental biologists, we have also incorporated the ability to define points or regions of interest in the image volume, so that the relative positions of these sampled regions with respect to these points of interest can be recovered.

We have implemented this approach as a software package called ‘SilentMark’, designed for general use. The package is a standalone GUI-based application written in Matlab, available at https://process.innovation.ox.ac.uk/software/p/13299a/silentmark-academic/1. The core algorithm has also been implemented in Python as a Web application and is hosted in the analysis section at https://dropletgenomics.com. Some limitations of our sampling-based approach are that it does not provide information on cell shape, or the total sum of fluorescence, or facilitate cell tracking. However, the unique aspect is that it allows one to process routine lower-quality images, investigate multiple cell compartments and quantify the spatial distribution of fluorescence within cells and tissues. Here we demonstrate the use of our algorithm on mouse pre- and post-implantation embryos, as well as mouse embryonic stem cell-derived embryoid bodies.

Results and discussion

To quantify protein levels in different cellular compartments and determine the relative proportion of the plasma membrane in different domains (apical, junctional and basolateral), we developed a sampling-based algorithm named ‘SilentMark’, implemented in Matlab and Python. The overall strategy relies on detecting well-resolved tissue elements, such as nuclei or membranes, and then analysing their closest environment to take samples of subcellular region fluorescence. Depending on which cell marker is available, the first step is to detect well-resolved objects such as nuclei or membranes, after which their closest environment is analysed. Care has been taken to implement efficient algorithms for three-dimensional convolution and searching, which enables work with images containing millions of voxels.

Nucleus-based fluorescence sampling

Available nuclear detection techniques differ in their accuracy, throughput and noise tolerance. One of the fastest and most robust methods for spherical object detection is based on the Hough transform [14], which is a tensor voting technique for the most probable circle positions within an image. The original algorithm is developed for two-dimensional image analysis. However, it can be extended for three-dimensional measurements by slicing the image into sections, analysing them as separate images and then deconvolving into three-dimensional objects, as done previously [9]. In this case, confocal image stacks were analysed as separate images to detect circular nuclear cross-sections, the coordinates of which were then merged to generate a three-dimensional point cloud of nuclear centres of mass (figure 1a). Using these nuclear coordinates as reference points, a region from the confocal image is cropped around each detected point and classified according to the intensity of the nuclear stain as either a cytoplasmic or nuclear region. The fluorescence from both regions is then averaged, with cytoplasm pixels being weighted according to their distance from the detection point (figures 1b and 4a,b). This was done to represent the probability that pixels closer to the nucleus are more likely to belong to the same cell as the nucleus. The fluorescence levels along with the nuclear coordinates were recorded and this formed the basis of segmentation and quantification.

Figure 1.

Schematic illustrating strategy for detecting the coordinates of the nuclear centre of mass (a) and for extracting nuclear and cytoplasmic fluorescence based on nuclear coordinates (b).

Figure 4.

Quantifying lineage determinants during mouse pre-implantation development. (a) Sampling algorithm for detecting three subcellular regions; (b) three-dimensional rendering of different sampling regions. Red and green are exposed and contacting membranes, respectively; grey is cytoplasm; and cyan is nucleus. (c) Three-dimensional rendering of a mouse embryo (grey, actin; red, YAP; cyan, DAPI); (d) 8-cell embryo in which cell membranes have been manually segmented and then membrane regions automatically sub-segmented into ‘exposed’ (blue), ‘junctional’ (white) and basolateral (red). (e) Error distribution of automated sampling-based quantitation compared with manual segmentation-based quantitation. (f) Levels of key transcription factors during development. (g) YAP distribution in the early pre-implantation embryo. (h) Spatial YAP patterning during pre-implantation development. (i) YAP and CDX2 relationship in the morula and early blastocyst. Scale bars in all panels represent 50 µm. N/C, ratio of nuclear to cytoplasmic fluorescence.

If a membrane stain was available in addition to the nuclear stain, it was possible to add another fluorescence reading for exposed and contacting membrane regions, which is particularly useful for looking at membrane-localized proteins. To detect the outside surface, the membrane fluorescence marker was thresholded and the outer pixels were assigned to the closest nucleus as an exposed membrane sample. To detect the contacting membrane samples, a line was drawn between two neighbouring nuclei, and the highest intensity region for membrane fluorescence along the line was designated as a contacting membrane sample and was assigned to both nuclei.

Membrane-based fluorescence sampling

In cases where nuclear markers are not available as a reference, membrane fluorescence, if available, can be used to create the points of reference for subsequent analysis. To detect membrane regions in space, a three-dimensional stack was processed using contrast-limited adaptive histogram normalization and Gaussian blur to reduce noise. Following this, a three-dimensional Canny edge detector [15] was applied to find membrane edges, to which spheres were fitted using an outlier-tolerant RANSAC algorithm. The detected membranes were represented by a location of the fitted sphere and a vector representing membrane curvature and direction (figure 2a). This information would then be the starting point for sampling membrane and cytoplasmic fluorescence intensities. The advantage of this method is that it looks for fragments of membranes, which is more tolerant to noise than other approaches that rely on total membrane segmentation. These detected membrane regions can now be used as sampling masks for membrane and cytoplasmic fluorescence. The mask considers a narrow region in space around the detected point to sample membrane fluorescence and draws a normal to the membrane surface to sample cytoplasm on either side. These masks are used on three-dimensional image stacks to analyse membrane and cytoplasmic fluorescence distribution in space (figure 2b). Membrane signal quantification can be improved by using commercially available software Imaris to manually outline membrane regions, which SilentMark then automatically partitions into contacting, exposed and junctional domains (figure 4d), where proteins are likely to play different roles.

Figure 2.

Schematic illustrating the membrane fragment detection algorithm (a) and analysis of membrane and cytoplasmic fluorescence based on membrane fragment detection (b).

Development of SilentMark software

We have developed these sampling-based algorithms into a software package designed for general use. The package is a standalone application, based on a GUI designed on Matlab. The algorithm has also been implemented in Python as a Web application. It accepts three-dimensional images as stacks of TIFF files or the Zeiss confocal microscopy ‘.lsm’ format. The user is required to enter two intuitive and robust parameters—an estimate of the diameter of the nuclei and an estimate of nuclear channel brightness for thresholding. These two minimal parameters allow one to work with different types of tissues and image qualities. All the remaining parameters needed for the analysis are derived and optimized by the software.

The software will output a list (CSV format) of detected objects (nuclei or membrane fragments), each of which will have a measure of nuclear and cytoplasmic fluorescence as well as membrane fluorescence, if it is labelled. The data will also contain object coordinates to describe spatial protein distribution, automatically calculated distance to exposed surface and distance to regions of interest, if they are designated manually in three dimensions. This information can then be processed with statistical software tools for further analysis.

Method validation with CAD images

To validate this probabilistic sampling-based approach, we first used computer-generated three-dimensional images to compare the accuracy of our automated sampling-based approaches against standard manual outlining (using Imaris; figure 3). Mock-up tissues composed of several spherical objects representing cells were created on Google SketchUp and converted to a stack of images, which were then filled with pixels of known brightness. Unsurprisingly, manual outlining was the most accurate method of image segmentation and, in addition to quantifying fluorescence intensity, could also be used to estimate the proportion of ‘exposed’ cell surface, that is, cell surface not in contact with surrounding cells. In general, intensity measurement errors were below 10%, except for feature detection relying on membrane fragments. The most likely source of error for membrane-based feature detection came from the detection algorithm while fitting spheres, which is responsible for identifying both membrane segment position and orientation.

Figure 3.

Accuracy of sampling-based image segmentation. (a) A section of a computer-generated image for testing image analysis methods. (b) A three-dimensional rendering of computer-generated objects for testing image analysis methods. Method accuracy was assessed by comparing measured fluorescence levels with known intensity values in computer-generated images. Four models, each consisting of four to six cells, were used in this test.

Method validation with pre-implantation embryos

The mouse pre-implantation embryo is a model system where cell fate is strongly influenced by tissue geometry [16,17]. The underlying basis for cell-type specification in the pre-implantation embryo has been studied in great depth and this large amount of existing information makes it ideal for validating our analytical approach.

During the process of the first cell fate decision in mammals, the protein YAP shuttles between the cytoplasm (where it is inactive) and the nucleus (active), in a tightly regulated manner. The differential localization of YAP protein in outside and inside cells ultimately determines cell fate in this context, and the ratio of nuclear to cytoplasmic YAP reflects the position of cells. In outside cells, YAP is in the nucleus, where it can interact with nuclear TEAD4 to drive the expression of tissue-specific genes such as Cdx2 to give rise to the trophectoderm. In inside cells, YAP is excluded from the nucleus and these cells go on to become the pluripotent inner cell mass [18-20].

To estimate the errors resulting from the sampling-based approach (figure 4a), we used manual cell outlining (figure 4d) to create a ground-truth dataset of pre-implantation embryos stained for YAP (figure 4b,c), (n = 67 embryos, comprising 20 2-cell, 15 8-cell, 17 16-cell, 12 32-cell and three 64-cell embryos; examples in figure 5 and electronic supplementary material, movie S1). For nuclear and cytoplasmic fluorescence measurements, mean error was lower than 5% (figure 4e), which indicates a low bias during the automated sampling. Most of the errors result from poor contrast of nuclear stain (which particularly occurs deeper into the tissue), as this makes it harder to correctly classify voxels as nuclear or cytoplasmic. Another source of errors is sampling outside the tissue, as the software does not calculate bounding tissue surfaces.

Figure 5.

(a) An example of a manually segmented 16-cell mouse embryo, with cell outlines in green and nuclei in yellow. (b–d) Example of a 16-cell mouse embryo stained for DAPI, YAP and CDX2. In all images, cell outlines are in yellow, stained with phalloidin (F-actin). The scale bar represents 50 µm.

We used our automated sampling-based approach to explore the relationship between YAP and CDX2 levels and the extent to which this correlates with a blastomere's position, as measured by exposure to the outside or distance of a cell from the surface. We used embryos stained for YAP, CDX2 and DAPI. Published single-cell RNA-seq data on mouse pre-implantation embryos [21] (GEO accession: GSE45719) indicate that Cdx2 becomes upregulated at the 8-cell stage (figure 4f, N = four 2-cell, four 8-cell, four 16-cell and three 32-cell embryos). Our analysis shows that at the 8-cell stage, blastomeres begin to show signs of segregating into two populations, with high or low nuclear to cytoplasmic YAP ratios (figure 4g and example in figure 5; N = 20 2-cell, 15 8-cell, 17 16-cell, 12 32-cell and three 64-cell embryos). More robust segregation of these two populations occurs at the 16-cell stage, approximately 12 h after Cdx2 transcript appearance at the 8-cell stage. Analytical detail makes it possible to demonstrate that populations are separated gradually, as YAP translocates to the cytoplasm in inside cells (figure 4h and electronic supplementary material, figure S1). It is also evident, that the CDX2–YAP correlation is stronger at the 32-cell stage than the 16-cell stage (figure 4i and example in figure 5; N = seven 8-cell, nine 16-cell and six 32-cell embryos). This is consistent with the presence of additional stabilizing mechanisms for robust lineage specification [22,23].

Method validation with embryonic stem cell-derived embryoid bodies

Similar to mouse pre-implantation embryos, mouse embryonic stem cell aggregates also display inside–outside patterning during embryoid body (EB) formation in vitro [24]. The outer layer of cells forms an endoderm-like layer around an inner core that is epiblast-like. These aggregates tend to be composed of crowded cells (figure 6a,b) and represent a more challenging task for both manual and automated segmentation or sampling approaches. To test its capabilities, we applied our sampling-based image approach to investigate endoderm formation during EB differentiation. We stained EB for the endoderm marker GATA6 (figure 6b), which specifically stains outside cells, and quantified the relative levels of nuclear GATA6 expression with respect to the surface of EB (n = 6). Quantitative statistical analysis (figure 6c) shows preferential GATA6 expression at the outer layer of the EB (distance to outer surface is smaller than two cell radii).

Figure 6.

Automated quantitative analysis of differentiated mouse stem cell-derived embryoid bodies and mouse embryonic heart. (a,b) Optical sections of a differentiating mouse stem cell EB stained with DAPI (a, cell nuclei) and GATA6 (b, marking outer endoderm layer). (c) GATA6 localization in the EB as a function of distance from exposed surface (n = 6 EBs). The approximately single-cell-wide outer layer is distinct from the remaining EB. (d–f) Optical section of an embryonic heart stained for F-actin (magenta in all panels, showing cell outlines), DAPI (d, grey, all cell nuclei), NKX2-5 (e, grey nuclei in cardiac crescent) and SRF (grey nuclei in cardiac crescent). (g) Maximum intensity projection of image volume of embryonic heart, with F-actin in magenta and NKX2-5 in grey. The green dots and line are regions of interest marking the head folds and cardiac crescent, respectively. (h) Localization and correlation of SRF and NKX2-5 transcription factors, with boxplots comparing proteins in different embryo regions. The blue and red dots were clustered by the Clara algorithm based on SRF and NKX2-5 levels and their distance from the heart and head folds. p-values are less than 0.001 (t-test). Scale bar in panel (a) represents 50 µm and in panel (f) 100 µm. N/C, ratio of nuclear to cytoplasmic fluorescence.

Method validation with embryonic hearts

While mouse pre-implantation embryos present an ideal case for testing automatic image segmentation, because it is possible to also manually outline the cells for comparison, they do not necessarily represent a typical mammalian tissue. A more representative mammalian tissue is the embryonic heart (figure 6d–g), with a larger number of cells and closely spaced nuclei.

To test the ability of our image analysis approach for this problem, we examined the localization of NKX2-5 and SRF, two transcription factors important for cardiomyocyte differentiation [25] in 8.0 days post coitum embryos (figure 6e,f). In order to represent protein levels of these two transcription factors, we used a ratio of nuclear to cytoplasmic fluorescence (N/C), as the latter is expected to represent background signal for internal normalization. To reduce the dynamic range of the measurements and normalize the distribution of the values, we used a natural logarithm of the ratio. In addition to these normalized fluorescence measurements, we also measured cell distance to two regions, the head folds and cardiac crescent, manually designated using our software (green dots and line in figure 6g).

The four variables formed a dataset, which was clustered using partitioning around medoids (PAM, R ‘Clara’ package) [26]. The method was chosen over other clustering approaches [27] owing to computational efficiency in dealing with large datasets. The analysis verified as expected that both proteins were preferentially expressed in the heart region and there was a region-specific correlation between nuclear levels of NKX2-5 and SRF (figure 6h).

In summary, we have presented a sampling-based three-dimensional image analysis approach, designed for relative quantification of protein levels in microscopy images of complex tissues. Our approach is the first to simultaneously quantify nuclear, membrane and cytoplasmic fluorescence, while also marking regions of interest to extract spatial information relating to the quantified signals. Automation of analyses enabled by our program will enable high-throughput quantification and statistical analysis of spatial protein organization. As three-dimensional microscopy is versatile and widespread, we expect that our publicly available open-source software will be useful not only for developmental biology but also more broadly, in the context of cell biology.

Methods

Mouse strains, husbandry and embryo collection

All animal experiments were carried out according to UK Home Office project licence PPL 30/3155 and 30/2887 compliant with the UK Animals (Scientific Procedures) Act 1986 and approved by the local Biological Services Ethical Review Process. All mice were maintained in a 12 h L : 12 h D cycle. Noon of the day of finding a vaginal plug was designated 0.5 dpc (days post coitum). To obtain embryos, C57BL/6 males were crossed with CD1 females (Charles River). Embryos of the appropriate stage were dissected in M2 medium (Sigma-Aldrich) at room temperature. Pre-implantation embryos were collected by oviduct and uterus flushing.

Embryonic stem cell culture

Embryonic stem cells were cultured in T-25 flasks coated with 0.1% w/v gelatin and used at passage no. 50–60. Cell culture medium consisted of 90% DMEM medium (Gibco), 10% fetal calf serum (FCIII, Hyclone), 1 mM sodium pyruvate, 1 mM non-essential amino acid mix, 0.1 mM beta-mercaptoethanol, 2 mM l-glutamine and 1000 U ml−1 leukaemia inducible factor (LIF). Cells were passaged every 2 days and seeded at approximately 400 000 cells per T-25 flask. Embryoid bodies were formed by aggregating cells in hanging drops and culturing for 2 days before fixation.

Antibodies

Vectashield with DAPI (H-1200) was purchased from Vector Laboratories, phalloidin atto647 was purchased from Sigma (65906-10MMOL). Mouse anti-YAP monoclonal (sc-101199) and goat anti-Nkx2-5(N-19) (sc-8697) polyclonal antibodies were purchased from Santa Cruz Biotechnology Rabbit anti-CDX2 (no. 3977S) monoclonal and rabbit anti-pERM(T567) (no. 3149P) monoclonal antibodies were purchased from Cell Signaling. Rat anti-uvomorulin/E-cadherin (U3254-100UL) monoclonal antibody was purchased from Sigma. Rabbit anti-SRF (PA5-27307) polyclonal antibody was purchased from Thermo Fisher. We used the following secondary antibodies from Life Technologies: goat anti-rabbit AlexaFluor 488 (A11008), donkey anti-goat AlexaFluor 488 (A11055), donkey anti-mouse AlexaFluor 555 (A31570), donkey anti-rabbit AlexaFluor 647 (A31573).

Immunofluorescence on pre-implantation embryos and embryoid bodies

Embryos and EBs were fixed for 20 min at room temperature in 4% paraformaldehyde (PFA). Samples were then washed three times in PBT-0.1% (phosphate-buffered saline with 0.1% Tween) for 10 min, permeabilized in PBT-0.25% for 40 min and washed again in PBT-0.1%. The tissues were transferred to a blocking solution for 1 h at 4°C. Primary antibodies were then added to the solution and incubated overnight at 4°C. The embryos were washed in PBT-0.1% and incubated for 1 h at 4°C in PBT-0.1% with the secondary antibodies, then subsequently washed two times in PBT-0.1% for 15 min and mounted in Vectashield with DAPI at least 6 h prior to imaging.

Immunofluorescence on post-implantation embryos

Dissected embryos were fixed for 1 h at room temperature in 4% PFA. The embryos were then washed three times in PBT-0.1% for 15 min, permeabilized in PBT-0.25% for 40 min and washed again three times in PBT-0.1%. The embryos were transferred to blocking solution overnight at 4°C. Primary antibodies were then added to the solution and incubated at 4°C. The embryos were washed three times in PBT-0.1% and incubated overnight at 4°C in PBT-0.1% with the secondary antibodies, then subsequently washed three times in PBT-0.1% for 15 min and mounted in Vectashield with DAPI at least 24 h prior to imaging.

Software and statistics

The software and graphical user interface were written in MATLAB (v. 2015b, Mathworks) and compiled using the in-built compiler for Mac OSX and Windows operating systems (https://process.innovation.ox.ac.uk/software/p/13299a/silentmark-academic/1). The algorithm has also been implemented in Python as a Web application (hosted by Droplet Genomics at http://dropletgenomics.com). Source code for SilentMark is available upon request. Three-dimensional image visualization was done using Volocity software (Perkin Elmer), and Imaris software (v 6.1) was used for manual cell outlining and analysis. R statistics package alongside ‘cluster’ library was used for automatic normalized cluster detection (PAM algorithm) and cell population analysis.

The first step of the algorithm is nuclei or membrane detection. This initial two-dimensional segmentation step can be done using a variety of approaches (including convolutional neural nets), but we found Hough circle detection to be the most efficient and robust. Detected nuclear or membrane feature coordinates on two-dimensional stacks are then convolved into three-dimensional space to build an interconnected graph. The second step then relies on efficient spatial searching, which segments stack voxels based on this graph. Delaunay triangulation was used for this search algorithm. Once the voxels are assigned to nearest nuclear or membrane centres, normalized intensity of nuclear or membrane stain is used to classify them as cytoplasmic. Then the algorithm measures all fluorescence markers in the voxels designated by type (nuclear, cytoplasmic or membrane).