Reconstructing temporal and spatial dynamics from single-cell pseudotime using prior knowledge of real scale cell densities

Modern cytometry methods allow collecting complex, multi-dimensional data sets from heterogeneous cell populations at single-cell resolution. While methods exist to describe the progression and order of cellular processes from snapshots of such populations, these descriptions are limited to arbitrary pseudotime scales. Here we describe MAPiT, an universal transformation method that recovers real-time dynamics of cellular processes from pseudotime scales by utilising knowledge of the distributions on the real scales. As use cases, we applied MAPiT to two prominent problems in the flow-cytometric analysis of heterogeneous cell populations: (1) recovering the kinetics of cell cycle progression in unsynchronised and thus unperturbed cell populations, and (2) recovering the spatial arrangement of cells within multi-cellular spheroids prior to spheroid dissociation for cytometric analysis. Since MAPiT provides a theoretic basis for the relation of pseudotime values to real temporal and spatial scales, it can be used broadly in the analysis of cellular processes with snapshot data from heterogeneous cell populations.

www.nature.com/scientificreports www.nature.com/scientificreports/  www.nature.com/scientificreports www.nature.com/scientificreports/ Next, we examined if MAPiT can also capture the heterogeneity within the different zones in the spheroid, typically a proliferative layer followed by quiescent and finally necrotic cells.
The DNA content of cells in the outermost spheroid layers exhibited a bimodality typically for proliferating cells with subpopulations in G1 (2N), S and G2/M (4N) phases (Fig. 5a). In contrast, the majority of cells in the inner layers remain in a quiescent G1/G0 (2N) state. The distinct distributions of additional markers, Ki-67, p27 and RNA content corresponded to this pattern (Fig. 5a).
We then studied the robustness of MAPiT, by comparing its performance to recover the spatial position using different pseudotime algorithms, markers and spheroid sizes. We first used the Wanderlust and Diffusion Maps algorithms 4,5 to obtain a pseudotemporal ordering. The order of cells was largely conserved in both algorithms (Fig. 5b, monotonously increasing data points), however the pseudotime values were clearly different between both algorithms (Fig. 5b, deviation from the diagonal). Transforming the DPT and Wanderlust axis to a distance scale with MAPiT resulted in almost identical spatial profiles of Ki-67, RNA and p27 (Fig. 5c), thereby proving that MAPiT results are not affected by the pseudotime algorithm. MAPiT should likewise provide identical marker trajectories unaffected by the choice of markers used to generate the pseudotime order. To validate this, we conducted "leave one out cross validation", where we took subsets of the markers as input to the Wanderlust and DPT algorithms, and compared the results obtained by MAPiT. MAPiT robustly provided correct distance profiles of all markers in all combinations, demonstrating that MAPiT is not influenced by the choice of inputs to the pseudotime algorithms (Fig. S1). Thus, MAPiT is not affected by the choice of markers or pseudotime algorithms, provided that the order of cells on the pseudotime scale reflects the sequence or directionality of the biological processes.
Spatial reconstruction by MAPiT also provides scope to significantly accelerate high throughput studies that make use of advanced 3D culture-based cell screens, such as spheroid based viability assays and drug effect screens. Indeed, spheroids are considered superior models in comparison to conventional cell cultures 24 , www.nature.com/scientificreports www.nature.com/scientificreports/ but spatiotemporal analyses still require cumbersome and manual labour-intensive work for the analysis of cross-sectional slices. MAPiT performed reliably in 3D reconstruction when processing flow cytometry data from dissociated spheroids of different sizes and ages, based on Ki-67, RNA and p27 amounts ( Fig. 5d-f). The profiles for Ki-67, RNA and p27 were conserved throughout spheroid sizes, indicating a dependence of these markers solely on the distance from the surface (Fig. 5f) and therefore their suitability for spatial reconstruction. Overall, robust spatial markers together with MAPiT therefore allow the rapid and versatile reconstruction of spheroids, providing a basis for studying spatiotemporal changes in other measurable variables.

Discussion
In summary, MAPiT provides a solution to a fundamental problem, namely the transformation of pseudotime to real-time or the true scale of a biological process. For example, snapshot single-cell data of cell populations can be converted to extract real-time kinetics of cellular processes and responses, which otherwise could only be obtained by live-cell microscopy, which is more complex, time consuming and limited by the availability of suitable live-cell reporters. By reconstructing spatial and temporal spheroid compositions from single-cell data, MAPiT provides insights to the evolution of cellular heterogeneity within tumour-like microenvironments and allows to understand how responsiveness to therapeutics manifests within spheroidal environments. Since MAPiT provides the means to not only employ flow cytometric data but data from any other single-cell based high-throughput multiplex measurement, such as CyToF or single-cell RNA-seq, it provides a foundation for high-throughput and high-content studies of 3D-spheroid models by recovering the spatial information lost during spheroid dissociation. By extension, applying MAPiT to other single-cell snapshot data, such as single-cell transcriptomics and proteomics data, might significantly improve the inference of complex regulatory processes and networks by recovering real temporal and spatial dynamics. = spheroids in three independent replicates. (f) Median profiles for Ki-67, RNA and p27, as obtained by MAPiT, were conserved throughout spheroid sizes (median of = n 3 spheroids, in three independent replicates). Notable deviations towards the ends of the profiles (dashed) and thus within the centre of the spheroids arise from inaccuracies due to the low number of cells available for analysis at these locations.
www.nature.com/scientificreports www.nature.com/scientificreports/ In the present paper we apply MAPiT to two examples where we obtained the distribution on the desired scale by theoretical view on the source of heterogeneity. Constructing a target distribution for other processes i.e. the most common application of pseudotime single-cell analysis on stem cell differentiation, might be more challenging.
In the case of differentiation process, it might be necessary to not only measure cellular components for single cell analysis, but also cell death and cell division rates at every stage of the process (further discussed in the Supplementary Information). A time scale separation argument can be used to regard cell death and cell division as processes perpendicular to the cell differentiation. A suitable choice of markers may then allow inference of cell death and division rates from single cell data. We envision that combination of ergodic theory and mathematical models (as discussed in the Supplementary Information) to map single-cell differentiation data with MAPiT on a real time scale. Label data (e.g. sampling timepoints) may likewise be used to support or validate MAPiT results. However, knowing sources and sinks along the process like cell death or cell division is critical for obtaining the right real-time distribution and subsequent correct transformation with MAPiT.
Recently, pseudotime algorithms were further developed to robustly recognise also branching processes in differentiation pathways 8,25,26 , providing scope to apply MAPiT to study differentiation dynamics in individual branches. MAPiT could be applied to all paths from the root to the end points on the respective branches treating other the flow of cells into other branches as sinks.
Overall, MAPiT is a robust and universal tool to recover temporal or spatial cellular trajectories from high-throughput, high-dimensional single-cell experiments. MAPiT can be combined with pseudotime algorithms, and a MATLAB implementation is available through GitHub (https://github.com/karstenkuritz/MAPiT) 16 .

MApit theory.
MAPiT is based on the measure-preserving transformation which states that one must conserve the area under the curve when transforming a probability distribution to another scale.
Consider a measure space λ L X ( , , ), where x is a set, L is a σ − ring of measurable subsets of X, and λ is the measure. Given a map τ from a measure space x s The mapping τ → s x : from pseudotime to real-time was obtained by solving equation (2) for τ, which then depends on the cumulative distributions of cells on both scales x s 1 Thus, by definition, the transformation τ requires knowledge of the distribution (or cumulative distribution) of cells on the pseudotime scale and the desired scale. If these distributions are positive (larger zero) over their support, then the cumulative distributions are monotonically increasing and the inverse exists. Once the mapping τ is known, one can apply the transformation to the joint densities of pseudotime and the observed quantities p s y ( , ) s to obtain the desired joint distribution of the true scale x and measured markers y Density in pseudotime for cell cycle data was estimated by kernel density estimation with linked boundary conditions to account for doubling of cell density during cell division 27 . Joint densities were calculated as sum of the product of the individual kernels Bandwidths h s and h y were derived from Silverman's rule 28 .
Distributions on the real scales. The distributions on the real scales can either be derived from theoretical considerations or from empirical measurements (discussed in Supplementary Information).
Cell density with respect to cell age in proliferating populations. For analysis of cell cycle-dependent processes with single-cell measurements, MAPiT requires the distribution of cells related to their cell cycle stage or equivalently their age a. Cell age, refers to the time since cell birth via cytokinesis. Our analysis is based on the following assumptions: (1) population is unperturbed and in its exponential growth phase, (2) no cell death in the population, (3) cell cycle progression is homogeneous. We thus restricted our analysis to unperturbed cell populations in their exponential growth phase with growth rate γ and cell cycle length T related by In such a case, the theoretical steady state age distribution of a cell population is given by 19 : The cumulative distribution and its inverse can be obtained in closed form: a a P y y ( ) 1 ln 1 2 (11) Thus, in case of an unperturbed cell population it is sufficient to know the growth rate of the population to recover cellular age with MAPiT and thus obtain the temporal changes related to cell cycle progression of measured markers from one single-cell experiment. We verified all assumptions by live cell imaging. In addition, light scattering characteristics in the flow cytometric datasets were used to probe the population for cell death.
Cell density in tumour cell spheroids with respect to distance from surface. Cell density depending on the distance from the surface was obtained from sphere geometry (Fig. 2b, Supplementary Information Fig. S2). Our analysis is based on the following assumptions: (1) Spheroids are radial symmetric, (2) all cells in the spheroid have the same size. We verified both assumptions by visual inspection of whole spheroids and spheroid slices. The volume of a sphere with radius r is given by The volume of a spheroid with necrotic core with radius = − r r d max( , 0) with d N being the distance from the surface where the necrotic core begins. The volume of a spherical shell at distance x from the surface of the spheroid with necrotic core is then given by M N M N Normalising (14) with the total spheroid volume V M results in the normalised volume with respect to the distance to the surface of the spheroid which represents the cumulative distribution of cells related to the distance from the surface MAPiT is furthermore based on the probability density function and the inverse of the cumulative distribution which can be calculated analytically www.nature.com/scientificreports www.nature.com/scientificreports/ Spheroid radius r and the radius of first appearance of a necrotic core r N were inferred from experiments, described in more detail in the Supplementary Information. For the present study we used r m 270 and a spheroid radius based on the linear regression = . + . r t t ( ) 19 2 22 4 (18) for spheroid growth as shown in Fig. 5e. cell culture. The human colon carcinoma cell line HCT116 was obtained from the Banca Biologica e Cell Factory of the IRCCS Azienda Ospedaliera Universitaria San Martino in Genoa (ICLC HTL95025). The geminin expressing non-small-cell lung cancer NCI-H460 cells have been described previously 13 . NCI-H460/geminin cells were maintained in RPMI 1640 medium (Gibco, 21875034) supplemented with 5% heat-inactivated fetal calf serum (FCS, Pan -Biotech GmbH; P303309) and HCT116 cells were cultured in RPMI 1640 medium with 10% heat-inactivated FCS at 37  C in a humidified incubator with 5% CO 2 . For generation of tumour cell spheroids, cells were transferred into Terasaki multiwell plates (Greiner bio-one; 653180) in a volume of 25 µl RPMI 1640 medium with a concentration of 4000 cells/ml. Thereafter, the plates were inverted to allow spheroid formation at the bottom of the emerging hanging drops and placed in humid chambers located in the incubator. Two to three days after seeding, formed spheroids were transferred to agarose-coated 96-well cell culture plates (Greiner bio-one; 655180 coated with 1.5% agarose (Carl Roth GmbH & Co. KG; 3810.3) in RPMI 1640 medium). Spheroid growth was monitored with an EVOS FL Cell Imaging System (Thermo Fisher Scientific Inc.) and spheroid diameters were determined from generated pictures using the Fiji software (distribution of ImageJ 29 ).
time-lapse imaging. NCI-H460/geminin cells were imaged for their total cell cycle length and length of G1 or S/G2/M phases by time-lapse fluorescence microscopy using the Cell Observer system (Carl Zeiss, Oberkochen, Germany) equipped with a humidified imaging chamber at 37  C and 5% CO 2 . Randomly chosen cells were manually tracked for at least one full cell cycle and geminin signal intensity was recorded. Cell trajectories were obtained using the Tracking Tool (tTt) and qTfy for single-cell tracking and quantification of cellular and molecular properties in time-lapse imaging data 30 . For comparison to MAPiT derived results, background signal θ 1 of cell trajectories and scaling factor 2 θ for single-cell trajectories were chosen to maximise the log-likelihood between MAPiT density and individual traces, where p is the density obtained from MAPiT and y i are live cell imaging values for geminin of cell i at time-points t k after cell division.

Generation of spheroid sections and immunofluorescence staining for Ki-67.
To generate spheroid sections, spheroids were fixed with 4% paraformaldehyd for 10 min at room temperature. Thereafter, spheroids were washed three times with PBS and finally kept in a sucrose solution (30% sucrose (Carl Roth GmbH & Co. KG; 4661) in PBS) at 4  C. After 48 h, the sucrose solution was removed and replaced with Tissue Freezing Medium (A. Hartenstein GmbH; TTEK). Such embedded spheroids were stored at 20 −  C and finally cut into 10 m µ slices with a CM30505 cryostat. Generated sections were mounted on Polysine Microscope Adhesion Slides (Thermo Fisher Scientific Inc.; 10219280), fixed with 4% PFA for 10 min at room temperature (RT) and washed twice with PBS for 5 min. Thereafter, permeabilization was carried out with 0.1% Triton X-100 in PBS for 10 min before the slices were incubated with blocking solution (5% FCS and 0.1% Triton X-100 in PBS) for 30 min at RT. Incubation with the primary antibody against Ki-67 (1:400, Cell Signalling Technology; #9449) was carried out in blocking solution for 1 h at RT. Thereafter, sections were washed two times with blocking solution before they were incubated with the secondary antibody Alexa Fluor 647-conjugated goat anti-mouse IgG (1:500, Thermo Fisher Scientific Inc.; A-21236) solved in blocking solution for 1 h at RT in the dark. Next, sections were washed with blocking solution and incubated with DAPI (1 µg/ml in PBS, Thermo Fisher Scientific Inc.; D1306) to stain DNA for 10 min at RT before they were covered with coverslips using Fluoromount-G ® (SouthernBiotech; 0100-01). Samples were dried and fluorescence was analysed with a LSM 710 laser scanning microscope (Carl Zeiss, Oberkochen, Germany) and the blue edition of the ZEN software. flow cytometric analysis. Tumour cell spheroids were dissociated with trypsin/EDTA (Gibco; 59418C) and single cells were fixed with 4% PFA in PBS for 10 min. Thereafter, permeabilization was carried out with 90% ice-cold methanol in PBS for 30 min on ice. Subsequently, cells were washed two times with a washing solution (5% FCS, 0.05% BSA (Sigma-Aldrich; A2058), 0.02% NaN 3 (Carl Roth GmbH & Co. KG; Ca4221) in PBS) before they were incubated with the primary antibodies diluted in washing solution for 1 h at RT (anti-Ki-67, Cell Signalling Technology; #9449 and anti-p27, Cell Signalling Technology; #3686). Next, cells were washed and www.nature.com/scientificreports www.nature.com/scientificreports/ incubated with the respective secondary antibodies, Alexa Fluor 488-conjugated goat anti-rabbit IgG and Alexa Fluor 647-conjugated goat anti-mouse IgG diluted in washing solution for 1 h (1:100, Thermo Fisher Scientific Inc.; A11008 and 1:500, Thermo Fisher Scientific Inc.; A21236). After washing, cells were incubated with 100 l µ Hoechst 33342 solution (10 g ml / µ Hoechst 33342 (Thermo Fisher Scientific Inc.; H3570) in PBS with 0.05% BSA and 0.02% NaN 3 ) for 1 h at 37  C, followed by the addition of 5 l µ Pyronin Y (stock solution: 100 µg ml / solved in ddH 2 O, Sigma-Aldrich; 83200) and incubation for 15 min at 37  C. Finally, cells were pelleted, dissolved in PBS and fluorescence was measured by flow cytometry with the MACSQuant analyser 10. NCI-H460/geminin cells were analysed as described in 13 . Data size and statistical reporting. This section reports statistics of the described experiments. We abbreviate the number of biological replicates, representing replicates under the same protocol but in different experiments, with n b . Technical replicates abbreviated by n t are obtained in identical experimental conditions, e.g. flow cytometric measurements of = n 3 t spheroids from the same culture. Furthermore, we denote the samples size, e.g. number of single cells in a single cell experiment, by n s .
Flow cytometry measurements of HCT116 spheroids were performed in three independent experiments (n 3 b = ). Each experiment was carried out with three spheroids (n 3 t = ). A samples size of = n s 10,000 cells was subsequently used for analysis with MAPiT. Ki-67 profiles in sliced spheroids are reported for one representative spheroid, with few defects caused by the slicing procedure ) wherein we recorded the geminin intensity of n 12 s = randomly chosen cells over one full cell cycle.
Signal intensities in the spheroid experiments were normalised to the median intensity of all cells up to a distance of m 150 µ from the surface. MATLAB code for comparison of MAPiT analysis with imaging data is available in the MAPiT repository (https://github.com/karstenkuritz/MAPiT) 16 . www.nature.com/scientificreports www.nature.com/scientificreports/