Continuous clonal labeling reveals uniform progenitor potential in the adult exocrine pancreas.

The tissue dynamics that govern maintenance and regeneration of the pancreas remain largely unknown. In particular, the presence and nature of a cellular hierarchy remains a topic of debate. Previous lineage tracing strategies in the pancreas relied on specific marker genes for clonal labeling, which left other populations untested and failed to account for potential widespread phenotypical plasticity. Here we employed a tracing system that depends on replication-induced clonal marks. We found that, in homeostasis, steady acinar replacement events characterize tissue dynamics, to which all acinar cells have an equal ability to contribute. Similarly, regeneration following pancreatitis was best characterized by an acinar self-replication model because no evidence of a cellular hierarchy was detected. In particular, rapid regeneration in the pancreas was found to be driven by an accelerated rate of acinar fission-like events. These results provide a comprehensive and quantitative model of cell dynamics in the exocrine pancreas.

Introduction

The pancreas is a metabolic organ with key roles in the digestive and endocrine systems. These functions are performed by exocrine and endocrine cells, respectively. Endocrine islet cells include, for example, α cells that generate glucagon and β cells that produce insulin. Likewise, the exocrine pancreas consists of enzyme-producing acinar cells and ductal cells that mediate transport of acinar secretions to the intestine. Acinar cells and ducts are closely connected by centroacinar cells, which form the histological transition between these cell types. The vast majority of pancreatic tissue consists of exocrine cell types; in particular, acinar cells (∼90%) (Jiang et al., 2020; Shih et al., 2013).

Although we have extensive knowledge of the factors associated with lineage specification in the pancreas, our insights into the dynamic organization of tissue architecture and clonal relationships therein are much less developed. Only recently have studies provided detailed insight into the developing embryonic pancreas (Kopinke et al., 2011; Kopp et al., 2011; Larsen et al., 2017; Pan et al., 2013; Schaffer et al., 2010; Shih et al., 2013; Solar et al., 2009; Zhou et al., 2007). Most notably, it was established that, early in development, a pool of multipotent progenitor cells at the tips of growing ducts generate ductal and acinar cells (Sznurkowska et al., 2018). Stochastic branching of the embryonic pancreatic ducts leads to expansion of the ductal system by a process referred to as branching morphogenesis. Gradually, this phase is followed by the presence of more lineage-restricted progenitors that will generate acinar cells or duct/islet cells (Kopinke et al., 2011; Kopp et al., 2011; Larsen et al., 2017; Pan et al., 2013; Schaffer et al., 2010; Shih et al., 2013; Solar et al., 2009; Zhou et al., 2007). In part, these studies were made possible by the extensive proliferation and expansion of developing tissue, which allows fate mapping using lineage tracing on practical timescales. In contrast, the adult pancreas is much more static, as also reflected by its proliferation rates; around 40% of the cells express Ki67 directly postnatally compared with less than 2% in the adult organ (Houbracken and Bouwens, 2017).

As a result, our understanding of the clonal dynamics at play in maintenance of the adult pancreas is rather limited. Despite several lineage tracing studies performed in the adult pancreas, the renewal dynamics of its main cellular constituents, acinar cells, remain uncertain (Jiang et al., 2020; Ziv et al., 2013). In particular, it remains a topic of debate whether all acinar cells have the ability to generate new acinar cells in homeostasis and regeneration or whether specific cell populations can be recognized that harbor stem cell properties. For example, it has been demonstrated that Bmi1+ and Nestin+ cells partially actively proliferate and display prolonged self-renewal properties (Sangiorgi and Capecchi, 2009; Wollny et al., 2016). On the other hand, Stmn+ and Dclk1+ cells are largely quiescent and generate larger acinar clones following tissue damage (Westphalen et al., 2016; Wollny et al., 2016). This suggests that these stem cell marker-expressing cells demarcate separate and distinguishable cell populations with divergent self-renewal and clonal potential in homeostasis. In notable contrast, it has been demonstrated that cells expressing the terminal acinar differentiation marker Elastase were also capable of successful regeneration (Desai et al., 2007), which implies that all acinar cells have some degree of self-renewal and progenitor capacity. Interpretation of these and other studies is complicated by inherent technical caveats. For example, clonogenic outgrowth in vitro might not reflect in situ cell potential. Furthermore, in vivo lineage tracing from cell-type-specific promoters that drive tamoxifen-inducible Cre recombinase in combination with a Cre-activatable reporter gene (e.g., Bmi1-Cre) is fundamentally biased because it leaves the non-marker-expressing populations untested. In addition, tamoxifen has been reported to induce diffuse cell death in the intestine and stomach (Huh et al., 2012; Zhu et al., 2013) and proliferation in the acinar cells of the pancreas (Ahn et al., 2019). This implies that studies that have aimed to study homeostasis have, in fact, analyzed regeneration.

To circumvent these shortcomings and define the clonal dynamics of the exocrine pancreas in an unbiased fashion for the entire murine lifespan, here we employed a marker-agnostic and tamoxifen-independent lineage tracing method. We combined the obtained temporal in vivo clone size data with a stochastic modeling approach to establish that the acinar compartment is a self-replicating system. We concluded that, in the adult organism, all acinar cells have an equal ability to participate in self-renewal events and that there is no indication of a strict cellular hierarchy. Furthermore, we found that regeneration following severe pancreatitis is mediated by events resembling fission of acinar glands, which allows rapid restoration of pancreatic volume.

Results

Marker-free lineage tracing in the pancreas

We aimed to characterize the dynamics of pancreatic tissue in an unbiased fashion independent of specific promotor activity or tamoxifen induction. To this end, we employed a proliferation-dependent tracing system that relies on an out-of-frame enhanced yellow fluorescent protein (EYFP) reporter gene in the Rosa26 locus that can be stochastically activated during replication following strand slippage of a [CA]30 dinucleotide repeat (the R26 mouse model; Figure 1A). This rare event places the reporter gene in frame, leading to permanent EYFP expression that allows fate tracing of the labeled cell’s progeny. Previously, this model has been employed successfully to characterize stem cell dynamics in the intestine (Kozar et al., 2013) and mammary gland (Davis et al., 2016).

Figure 1

Marker-free lineage tracing in the pancreas

(A) Schematic of the R26 mouse model. All cells contain a [CA]30 dinucleotide repeat and a EYFP (enhanced yellow fluorescent protein) reporter gene in the Rosa26 locus. During replication, strand slippage of the [CA]30 repeat will result in expression of in-frame EYFP.

(B) Representative confocal images showing EYFP labeling (arrowheads) in different pancreatic cell types, co-staining with amylase for positive acinar cells (purple), Dolichos biflorus agglutinin (DBA; red) for ductal cells, and EYFP+ cells of staining-negative islet cells (dashed white line). Nuclear stain, Hoechst (blue). Scale bar, 100 μm. The pie chart shows quantification of EYFP+ clones per indicated cell type (n = 8,855 clones from N = 55 mice).

(C) UMAP (Uniform Manifold Approximation and Projection) embedding EYFP+ (green dots) and EYFP− (purple dots) acinar and ductal single-cell transcriptomes from dissociated pancreata.

(D) Gene expression feature plot of representative markers of ductal (Krt19, Sox9, and cftr) and acinar (Cela2a, Cela3b, and Try4) cells.

(E) Example microscopy image of a pancreas of a 521-day-old mouse. Colors as for (B). Right panels show magnifications of the indicated areas with corresponding numbers: 1, single EYFP+ cell; 2, single EYFP+ acinus; 3, cluster of multiple EYFP+ acini. Scale bars, 1 mm (main image) and 100 μm (magnification images).

(F) 3D image of a cleared pancreas of a 341-day-old mouse. EYFP+ cells are yellow, DBA+ duct cells are red. Magnifications as for (E). Scale bars, 200 μm (main image) and 100 μm (magnifications).

(G) Representative immunofluorescence images and diagrams showing different EYFP+ fractions (yellow) of acini, as defined by laminin staining for the basement membrane (purple) of pancreata of different ages as indicated. Scale bar, 50 μm.

(H) Histograms showing mean relative frequencies of EYFP+ fractions of acini at the indicated ages.

(I) Example confocal images of a pancreatic section and schematics of EYFP+ clusters, defined as clones of 2 or more EYFP+ acini (indicated by a black line in the diagram) of EYFP+ acini (purple line) at the indicated ages.

(J) Mean number of EYFP+ acini per cluster of multiple acini at the indicated ages.

(C and D) N = 7 mice; 3,617 EYFP− and 1145 EYFP+ cells were profiled in total.

(H and J) N = 8 mice (2 per time point); weighted average of 18 and 44 days, 129 and 227 days, 313 and 521 days, and 632 and 716 days; the number of analyzed acini (H) and clusters (J) are indicated per time point. Data represent mean ± SD.

In the adult pancreas of R26 animals, clones were readily detected in all three main cell lineages; i.e., acinar, ductal, and islet cells (Figure 1B). The relative number of labeled clone types directly reflected the relative abundance of these cell types in the pancreas (Figure S1A). A very small fraction of EYFP+ clones comprising more than one cell type was detected, indicating that lineage conversion in the pancreas is very rare. The fraction of such mixed clones was highest in young mice and decreased with age, which suggests that these clones originate from multipotent progenitor cells during development, whereas in adult mice, lineage-specified clones accumulate (Figure S1B). Alternatively, these mixed clones could arise from two independent tracing events occurring sporadically in close proximity. In this study, we predominately focused on the acinar compartment because this constitutes the majority of the exocrine pancreatic tissue. To confirm unbiased labeling of cells, we performed single-cell RNA sequencing and immunohistochemistry analyses to compare the EYFP+ and EYFP− cell populations (Figures 1C, 1D, and S1C–S1I). These analyses revealed that both populations are indistinguishable based on gene expression profiles, stem/progenitor marker expression, and proliferation rate. This is a key finding because it confirms that the EYFP+ clones are representative of the total acinar compartment of the pancreas and that no molecularly distinct population is labeled preferentially.

Importantly, the R26 mouse model allows adequate demarcation of acinar clones, which enabled us to accurately quantify clone size in 2D tissue sections (Figure 1E) as well as in 3D reconstructions of tissue whole mounts (Figure 1F). This is a vital observation to ensure adequate quantification of clone numbers and clone sizes. In adult mice, we detected a wide range of clone sizes ranging from individual labeled cells to large clones that span multiple acinar structures (Figures 1E and 1F).

More specific analysis of mice from a wide age range (36–650 days) revealed that in younger mice (36 days), the majority of clones were very small, comprising few cells only occupying a small fraction of acinar structures, whereas in older mice (184–650 days), the proportion of clones that occupy a complete acinus accumulated (Figures 1G and 1H). This confirmed that cells within an individual acinar structure are clonally related and that acinar compartmentalization is a constraint to further clonal expansion. Also, it suggests that, in homeostasis, acini slowly convert to monoclonality, resembling the dynamics in the intestine and other epithelial tissues (Leushacke et al., 2016; Van Keymeulen et al., 2011; Winton et al., 1988). Furthermore, clones that consist of multiple acinar structures were detected mainly in older mice (Figures 1I and 1J). The vast majority of these multi-acinar clones consisted of completely EYFP+ acini (Figures S2A and S2B), suggesting that fission-like events contribute to labeled acini in these clusters. This interpretation is further enforced by the finding that the proportion of single EYFP+ cells surrounding clones decreased with labeled clone size (Figures S2C and S2D). If spillover from one acinus to another is responsible for non-fully converted acini close to each other, then one would expect larger clones to have a higher chance to spill over labeled cells compared with small clones. Our data showed the opposite result, which makes spillage events from one acinus into a neighboring structure unlikely (Figures S2C and S2D).

These observations suggest that the adult pancreas is a dynamic tissue in which turnover of acinar cells results in clone purification in acinar structures and that expansion of clones beyond individual acini likely results from sporadic events resembling acinar fission.

A quantitative model of acinar cell dynamics

To quantitatively analyze the observed clone size distributions and obtain a more comprehensive understanding of acinar clone dynamics, we established a simple stochastic model (Figure 2A). The model captures the dynamics within acini as well as expansion of clones by acinar fission. Details regarding the model can be found in Methods S1. In brief, the intra-acinar dynamics are modeled by counting EYFP+ cells in an acinus as a continuous time-discrete random walk on 0,1,2,…,N, where N > 0 is the average number of clonogenic cells within an acinus, which we assume to be constant in time. The acinar cells display a turnover rate of ʎ (divisions/day, ʎ > 0), which mirrors loss of an acinar cell and replacement of that cell by division of another cell in the same acinus. Importantly, at every replacement event there is a probability μ (mutations/division, 0 < μ < 1) that the R26 locus mutates and that an EYFP+ clone is established from a previously EYFP− cell. Based on the experimental data showing that labeled and non-labeled cells have similar proliferation rates (Figures S1F and S1G), we assume that EYFP+ and EYFP− cells have an equal ability to replace vacant spaces within an acinus. Furthermore, acini are subject to fission at rate κ (fissions/day, κ > 0), by which cells present in the acinus are segregated randomly into two new acini. A final assumption is that the number of acini in the adult pancreas is stable over the lifespan of the mice; thus, a fission-like event not only created a new acinus but also removed a different random acinus. This was confirmed experimentally (Figures S2E and S2F).

Figure 2

A quantitative model of acinar cell dynamics

(A) Diagram explaining the mathematical model (progenitor model) for defining acinar cell homeostasis, with the following parameters: number of cells per acinus (N), turnover rate ʎ (divisions/day), mutation rate μ (probability/division), and fission rate κ (fissions/day). PPA, partly populated acinus; WPA, wholly populated acinus. See Methods S1 for details.

(B) Heatmap depicting goodness of the progenitor model fit for the shown combination of parameters for the turnover rate and fission rate (large heatmap) and the number of cells per acinus (right heatmap). Lower values indicate better fit. Error bars represent an area with a similar goodness of fit. The best fit for the stem cell model is also shown (green).

(C) Inference of the number of clonogenic cells per acinus of the model with the experimentally derived mutation rate, number of cells per acinus (black dashed line ± SEM), and turnover rate. The acinus fission rate was inferred for each x value independently. The goodness of fit was determined as the residual sum of squares of the number of WPAs and the number of acini per clone corresponding to the inferred fission rate. Specific cases are noted as the stem cell model, where 1 cell per acinus is capable of proliferation, and the progenitor model, where all cells are able to proliferate.

(D) Distribution of clone sizes for all time points of the training set (bars; N = 8 mice, two combined per time point) and the inference of the progenitor (purple line) and the stem cell (green line) model using the optimal inferred values as shown in (B). The difference between the model distributions and the training set distributions is quantified by the Kolmogorov-Smirnov distance (KSD) (inset), where lower values indicate a better fit.

(E) Comparison of experimentally derived number of cells per acinus inferred to 3D (black line; mean ± SEM, n = 640 acini from N = 8 mice) and model-inferred number of cells per acinus in the progenitor model (top panel, purple line) and stem cell model (bottom panel, green line), as determined in (B) (Methods S1). The shaded area indicates 95% confidence intervals of the model.

(F) Top panels: immunofluorescence for proliferating cells (Ki67+, yellow) in the pancreas during postnatal development (18 days) and homeostasis (744 days). Nuclear stain, Hoechst (purple). Scale bar, 100 μm. Bottom panel: quantification of the percentage of Ki67+ cells in the pancreas during the lifetime of mice (N = 34 mice). The mean percentage of proliferating cells in the adult pancreas (from ~100 days, dashed black line) and model predictions of the progenitor model (purple line) and the stem cell model (green line) are depicted (Methods S1). The difference in distance between the predicted percentage of Ki67 and the experimental value for the stem cell model and progenitor model was compared by the Mann-Whitney U test.

Next we performed 5 × 106 simulations capturing a wide range of the key parameters (N, ʎ, and κ) while fixing the mutation rate μ to one determined experimentally previously (Kozar et al., 2013). Using the data from the earlier-described set of mice (N = 8, 36–650 days, paired per two age-matched mice, training dataset) in combination with a goodness of fit method, we estimated the best combination of parameters that optimally describes the steady increase in wholly populated acini (WPAs) per 10,000 cells and the average number of acini per clone for the different age groups (Figures 2B and S3A–S3D; Methods S1). This analysis revealed that the optimal fit with the experimental data was N = 15 cells, ʎ = 0.036 divisions/day, and κ = 7.9 × 10−4 fission/day (Figures 2B and S3A–S3D; Methods S1). Indeed, using this combination of parameters, the accumulation of WPA clones with age and the increase in the number of acini per clone were reproduced accurately by the model (Figures S3E and S3F). This indicates that the model as well as the inferred parameters accurately describe the behavior of the clones present in pancreatic tissue.

For completeness of our analyses, we compared the stochastic model with competing models where we assume acinar cell heterogeneity. The two extremes of the model are cases where all cells in an acinus are clonogenic (the progenitor model) and where only a single cell has long-term clonogenic potential (the stem cell model). For both of these models, we inferred the parameters that best described the clone size distribution data we obtained in the training dataset (Figures 2B and S3G). When comparing the model fits with the experimental data, the progenitor model described the experimental data more accurately (Figures 2C and 2D). Most strikingly, the fit of the model improved with assuming more long-term clonogenic cells per acinus, arguing for a progenitor model (Figure 2C). The total number of cells per acinus inferred using the progenitor model was very similar to the number of cells per acinus we quantified, in contrast to the stem cell model, which underestimated the number of cells per acinus (Figure 2E; Methods S1). This implies that, in principle, all acinar cells serve as clonogenic progenitor cells with an equal probability to replace lost acinar cells. In addition, it indicates that the progenitor model is superior to the stem cell model in inferring parameter values that can be readily validated. Similarly, the replacement rate estimate of the progenitor model is also in agreement with the proportion of proliferating cells detected by Ki67 immunofluorescence and significantly more accurate than the prediction by the stem cell model (Figure 2F; Methods S1). This lends further support to the notion that key features of pancreatic homeostasis are captured accurately by the progenitor model. Finally, the estimated fission rate of acini was predicted to be 7.9 fission events per day per 10,000 acini without a significant difference between the progenitor or stem cell models. Therefore, fission is a rare event in the adult pancreas and could explain why acinar fission has not received much attention in the literature on pancreatic renewal mechanisms in homeostasis.

Next we evaluated how robust our model is to modifications of the inferred parameters. Fixing the number of cells per acinus to 13, which is the average number of acinar cells we determined visually (Figure 2E), resulted in highly similar estimates of the other parameters (Figure S4A). Similarly, inferring the mutation rate of the R26 locus using the model resulted in μ = 2.1 × 10−4 (Figures S4B–S4D), which is in the same range as the mutation rate we established experimentally previously in mouse embryonic fibroblasts and the murine intestine (Kozar et al., 2013). In combination, these analyses establish that, in contrast to the stem cell model, the progenitor model accurately infers parameters for which we have other (direct) means to determine with high fidelity and truthfully captures the clonal dynamics in the adult pancreas.

Validation of stochastic acinar cell dynamics in the adult pancreas

Our studies indicated that the vast majority of acinar cells in the adult pancreas are able to participate in stochastic replacement events to renew acinar cells. This suggests that there is no reason to assume that a specific and rare population of stem cells underlies renewal dynamics. In fact, as we demonstrated, a hierarchical model does not accurately describe the properties of the experimental data and infers parameters that are demonstrably incorrect. Especially in older mice, sizeable clones could be detected that most likely result from sporadic budding or fission of preexistent acini, which enables rapid clonal expansion. Also, in this process, there is no reason to assume a hierarchical organization of the cells that make up this tissue. To lend further support to our proposed model of acinar dynamics and avoid overfitting, we aimed to validate our findings so far (based on the training dataset) in a much larger panel of mice of a wide variety of ages (the validation dataset). We sacrificed R26 mice at various ages and quantified the number of acinar clones as well as their sizes. In this extended cohort of mice, we detected a steady increase in the average clone sizes with increasing age and a wide range of clone size distributions (Figures 3A, 3B, and S4E). These large clone size differences can be explained by the stochastic nature of our model, the different times when the CA30-EYFP label was turned on, in combination with the stochastic clone expansion dynamics. This clone size variation is equivalent to the wide range of clone size distributions detected previously in the intestine, skin, and esophagus and point to neutral and stochastic clone dynamics (Clayton et al., 2007; Doupé et al., 2012; Lopez-Garcia et al., 2010; Murai et al., 2018; Snippert et al., 2010).

Figure 3

Validation of stochastic acinar cell dynamics in the adult pancreas

(A) Schematic of the experimental procedure for clonal analysis in the pancreas.

(B) Heatmap showing the relative clone size frequency of EYFP+ acinar clones in the pancreas (indicated by color in the heatmap) per binned clone size (in rows) over time (columns). Ages of the mice are indicated in the top bar plot (range, 16–758 days). Number of acinar clones per time point are indicated below the heatmap. Clone sizes of the mice of the training set (Figure 2D) are also depicted (N = 55 mice, 1 mouse per time point, 2 mice for 18 days).

(C) Graph showing the fraction of EYFP+ acinar cells in pancreata of mice of different ages, determined experimentally (purple dots) and model predicted for the progenitor model (purple line) and the stem cell model (green line), using experimentally derived parameters for turnover rate and number of cells per acinus (Methods S1) (N = 47 independent mice). 95% confidence intervals (CIs; purple and blue shade) were determined computationally.

(D) Dissimilarity of the experimental clone size data (KSDs) and the progenitor model (purple) or stem cell model (green) on the validation set (N = 47 mice). Lower values indicate a better fit. The difference between distributions of distances was calculated by Mann-Whitney U test.

Most notably, the straightforward model we established before, based on few assumptions and only a limited number of mice, accurately predicted the behavior of the fraction of EYFP+ cells and the number of EYFP+ acini within clones (Figures 3C and S4F). This demonstrates that the derived model accurately predicts tissue dynamics beyond the mice that were used to derive the parameters and that features of the system that were not used in training the model are predicted faithfully. This is in contrast to the competing stem cell model, which described the experimental data significantly less accurately (Figures 3C and 3D). This confirms that we have established a robust and comprehensive framework to recapitulate acinar cell dynamics in homeostasis. Furthermore, it corroborates that all pancreatic acinar cells are bona fide progenitors that can maintain tissue integrity and that there is no reason to assume a strict tissue hierarchy, as suggested previously (Sangiorgi and Capecchi, 2009; Westphalen et al., 2016; Wollny et al., 2016).

Acinar fission underlies pancreatic regeneration

It has been described that, under conditions of repair after severe tissue damage, stem cells or progenitors drive regeneration of the pancreas (Pan et al., 2013; Westphalen et al., 2016; Wollny et al., 2016). This was, for example, concluded from the observation that cells expressing Stmn1 increased rapidly in a model for pancreatitis (Wollny et al., 2016) and that Dclk1+ cells following partial pancreatectomy give rise to large clones (Westphalen et al., 2016). To evaluate the effect of regeneration on pancreatic cell dynamics, we employed R26 animals, in which pancreatitis was induced using cerulein injections (Carrière et al., 2011; Lerch and Gorelick, 2013; Willemer et al., 1992; Figure 4A). Following this, we detected a clear inflammatory response and degeneration of acinar structures, which was largely restored after approximately 3 weeks (Figure 4B). This phase was accompanied by a transient increase in proliferation, as evidenced by Ki67+ immunofluorescence (Figure 4C). Tissue regeneration following induction of pancreatitis in R26 mice was accompanied by significant, and rapid expansion of EYFP+ clone sizes (Figures 4D and 3E). Importantly, the number of clones was not affected significantly (Figure S4G). These data suggest that extensive proliferation of a small number of stem cell-derived clones cannot drive tissue regeneration. If this was the case, we would have observed a small number of very large clones and a rapid increase in clone number that reflects extensive proliferation in a limited number of lineages. To formally evaluate the mode of tissue repair following pancreatitis, we inferred the most likely set of parameters that captures the rapid transition from the control situation to the post-pancreatitis condition (Figure 4F). These indicated that the key factor to drive tissue regeneration following cerulein-induced pancreatitis is a marked increase in acinar budding or fission-like processes. In the 2–3 weeks following pancreatitis, a minimum 46-fold increase in acinar fission was predicted. This increase in fission rate might be even higher because the model is likely to underestimate fission events (Methods S1). Indeed, numerous examples of adjacent labeled acini with structures that suggest them to be derived from a shared progenitor acinus could be detected in this period (Figures 4G–4I, S2C, and S2D). These data demonstrate that acinar fission drives regeneration instead of a specialized acinar subpopulation of cells.

Figure 4

Acinar fission-like events underly pancreatic regeneration

(A) Schematic depicting the experimental setup for clonal analysis after cerulein-induced pancreatitis. R26 mice received an intraperitoneal (i.p.) injection with cerulein (50 μg/kg) every hour 7 times a day for 2 consecutive days. At 10 different days post induction (DPI; range, 2–58 days), the pancreata were isolated for clonal analysis.

(B) Representative hematoxylin and eosin (H&E) images of pancreata during cerulein-induced pancreatitis (2–16 DPI, image is 4 DPI), post-pancreatitis (28–58 DPI, image is 28 DPI), and an age-matched control pancreas of a 129-day-old mouse. Scale bar, 50 μm. Right panel: quantification of the tissue area of H&E staining of pancreata during pancreatitis (2–16 DPI), post-pancreatitis (28–58 DPI), and age-matched control mice. 3, 6, and 3 mice are quantified, respectively.

(C) Immunofluorescence images of Ki67+ proliferating cells (yellow) of a pancreas 4 DPI, 28 DPI, and an age-matched control mouse (129 days old). Nuclear stain, Hoechst (blue). Scale bar, 100 μm. Right panel: quantification of the percentage of Ki67+ cells among all cells during tissue regeneration (2–16 DPI), after tissue regeneration (28–58 DPI), and in age-matched control mice. N = 4, 6, and 6 mice respectively.

(D) Relative clone size frequency of EYFP+ acinar clones in the pancreas of age-matched control mice during pancreatitis (2, 4, 14, and 16 DPI, columns) and post-pancreatitis (28, 30, 42, 44, 56, and 58 DPI, in columns) (indicated by color in the heatmap) per binned clone size (in rows) (pancreatitis, n = 734 clones from N = 4 mice; post-pancreatitis, n = 539 clones from N = 6 mice) and age-matched control mice (61–129 days; n = 201 clones from N = 4 mice; see also Figure 3A).

(E) Average clone size in pancreata post-pancreatitis (28–58 DPI, n = 539 clones from N = 6 mice) compared with age-matched control mice (61–129 days; n = 201 clones from N = 4 mice).

(F) Heatmap depicting goodness of fit to the progenitor model after recovery of pancreatitis (black square). Best fit during homeostasis as found in Figure 2B is indicated by white square. Lower values indicate better fit. Error bars represent area with similar goodness of fit.

(G) Violin plots depicting the model-predicted and experimental data (dots) mean of the number of cells per clone (control, n = 201 clones from N = 4 mice; post-pancreatitis, n = 539 clones from N = 6 mice) and the number of acini per clone (control, n = 43 clones from N = 1 mouse; post-pancreatitis, n = 539 clones from N = 6 mice) in age-matched control mice and regenerated mice.

(H) Histogram showing the relative frequency of the number of EYFP+ acini per cluster post-pancreatitis (28–58 DPI, n = 81 clusters from N = 6 mice) and healthy control mice (18–716 days, n = 66 clusters from N = 8 mice).

(I) Immunofluorescence images of EYFP+ progeny derived from one labeled acinar cell and possible acinar fission events; EYFP+ cells (yellow), laminin (purple), DBA+ cells (red). Nuclear stain, Hoechst (blue). Scale bars, 100 μm and 50 μm.

(B, C, E, and H) Data represent mean ± SD. Groups (control versus post-pancreatitis) were compared using unpaired two-tailed Student’s t test; p values are indicated in the graphs.

Discussion

We have described a marker-independent lineage tracing approach to comprehensively characterize the dynamics of the acinar cell compartment during the lifetime of the organism. Specifically, clonal labeling in conjunction with a simple stochastic model revealed that two distinct but related processes are at play to ensure tissue maintenance that operate at distinct spatiotemporal scales. First, our data indicate that acinar cells replace each other at a low but continuous rate during the lifespan of the mouse. Importantly, our results demonstrate that all acinar cells participate in this ongoing replacement and do not indicate the presence of an intrinsic cellular hierarchy. By comparing our progenitor model with a stem cell model we find that (1) the progenitor model more accurately describes the training data, (2) the progenitor model more accurately describes the data in the validation dataset, and (3) the parameters inferred by the progenitor model (number of cells per acinus, proliferation) are more in line with those detected experimentally. Second, at a much lower rate than during development and also much less frequent than intra-acinar replacements, our data suggest ongoing acinar fission. In the unperturbed pancreas, we estimate the rate of acinar gland fission to be 7.9 × 10−4 events per gland each day. This rate is much increased to enable rapid tissue regeneration following cerulein-induced pancreatitis.

Importantly, our unbiased clonal tracing data do not conform to a model in which rare stem-like cells drive tissue regeneration. Instead, they are in line with a self-replication model of differentiated acinar cells (Jiang et al., 2020) by intra-acinar replacement and gland duplication events.

Notwithstanding our conclusions regarding the functional homogeneity of acinar cells, previous studies have, in fact, reported significant phenotypic heterogeneity in acinar cell types. For example, single-cell sequencing studies have revealed molecularly distinct cell types within the acinar population that display different degrees of proliferation (Muraro et al., 2016; Tosti et al., 2021; Wollny et al., 2016). However, it is important to realize that cellular phenotypes and associated marker expression do not necessarily reflect stem cell or progenitor potential, as we uncovered previously in the intestine (Kozar et al., 2013; Vermeulen et al., 2013). Similarly, lineage tracing studies from distinct marker-expressing cells only reflect the properties of that cell population at the moment of induction. Hence, these analyses do not fully capture the (possibly widespread) cell plasticity that allows the changes in cellular phenotype which accommodate tissue renewal during homeostasis or repair. Especially in organs with a relatively slow cellular turnover, pulse-chase experiments do not reflect the actual long-term clonogenic potential of putative subpopulations. Indeed, our data, in combination with these earlier studies on phenotypical heterogeneity, suggest that extensive plasticity between putative acinar subsets should be present.

Limitations of the study

It will be interesting to further detail the conversions of cell states in follow-up research, which would require temporal integration of assessment of phenotypical and functional heterogeneity. A limitation of the current R26 mouse model to investigate cell lineage specification, and possibly conversions, in more detail is the low incidence rate of labeling. Therefore, our model is unlikely to capture multi-lineage clones that are defined early in development during a relatively short time frame. To capture such dynamics, a model that induces labeling at an accelerated rate (e.g., by using a CA30 model with a longer CA stretch) would be required. The low labeling rate of our model would also not be able to formally rule out that a very rare pool of quiescent stem cells exist. However, we can exclude that even if these exist, they do not contribute to pancreatic homeostasis in a numerically relevant fashion.

Furthermore, it will be important to investigate how our findings relate to formation of pancreatic ductal adenocarcinoma (PDAC) and, in particular, how the dynamics of the transforming cell populations change by activating Kras mutations during development of this malignancy. Similarly, at this moment it is unclear how our results relate to the dynamics of established PDACs. Also, in this case, conflicting results have emerged on the putative role of cancer stem cells in driving PDAC growth and progression (Hermann et al., 2007; Kong et al., 2011; Kopp et al., 2012; Li et al., 2007; Westphalen et al., 2016). Again, the inherent difficulties of testing functional properties of specific marker-expressing cell populations are equally relevant for the study of cancer tissue. Therefore, also in that context, marker-free, unbiased lineage tracing strategies are required to resolve the clonal dynamics of PDAC. This will form the focus of future investigations.

STAR★Methods

Key resources table

Resource availability

Lead contact

Further information and requests for resources and reagents should be directed to and will be fulfilled by the lead contact, Louis Vermeulen, l.vermeulen@amsterdamumc.nl.

Materials availability

This study did not generate new unique reagents

Experimental model and subject details

Genetically modified animals

The R26 (C57BL/6J background, RRID: IMSR_JAX:006148) mouse model has been previously described (Kozar et al., 2013). All animal experiments were performed in accordance with national guidelines and approved by the Animal Experimentation Committee at the Academic Medical Center in Amsterdam (LEX18-4667-2-01).

Method details

Lineage tracing and cerulein-induced acute pancreatitis

For clonal tracing experiments, mice (males and females) were kept in regular housing conditions until sacrificed by cervical dislocation, and pancreata were isolated at the indicated time points. Acute pancreatitis was induced in adult mice (60-70 days old) as previously described (Carrière et al., 2011). Briefly, mice received intraperitoneal injections with cerulein (Sigma-Aldrich, C9026) in PBS (50 ug/kg) every hour 7 times a day for 2 consecutive days. Metamizol (Sigma-Aldrich, 46232) in PBS (200 mg/kg) was orally administered as analgesia, 3 times daily on the days of cerulein treatment, without interfering with the cerulein-induced pancreatitis (Stumpf et al., 2016). Mice were sacrificed and pancreata were isolated at the indicated time points after the last cerulein injection.

Tissue isolation and processing for imaging

Immediately after isolation pancreata were fixed overnight using 4% paraformaldehyde followed by preservation in a 20% sucrose solution for at least 12 hours at 4°C before tissue was stored at −80°C. For image analysis 10-μm-thick frozen tissue sections were taken from different locations within the pancreas at intervals of > 100 μm, to prevent double sectioning of clones.

Immunostaining and immunohistochemistry

Fixed frozen pancreatic tissue sections were first incubated in Ultra-vision protein block (Thermo Scientific, TA-125-PBQ), to avoid aspecific staining. Next, the following primary antibodies were diluted in normal antibody diluent (KliniPath, VWR International, Radnor, PA, ABB999), applied on sections and incubated overnight in a humidified chamber at 4 °C; rabbit anti-amylase (Sigma, A8273, 1:500, AB_258380), chicken anti-GFP (Abcam, ab13970, 1:1000, RRID: AB_300798), rabbit anti-Ki67 (Sigma, SAB5500134, 1:200, RRID: AB_2892217), rabbit anti-laminin (Novus Biologicals, NB300-144SS, 1:1000, RRID: AB_921870), rabbit anti-Stathmin 1 (Proteintech, 11157-1-AP, 1:100, RRID: AB_2197114), rabbit anti-DCLK1 (Abcam, Ab31704, 1:500, RRID: AB_873537), mouse anti-nestin (Santa Cruz, 10C2 sc-23927, 1:100, RRID: AB_627994), rabbit anti-Ki67(SP6) (NeoMarkers, 9106S807E, 1:100, RRID:AB_149707). The following secondary antibodies were diluted in normal antibody diluent and incubated for at least 1 hour at room temperature; Alexa 647 goat anti-rabbit (Thermo Fisher, A-21244, 1:500, RRID: AB_2535812), Alexa 488 goat anti-chicken (Thermo Fisher, A11039, 1:500, RRID: AB_2534096), IR Dye 680 donkey anti- rabbit (Li-Cor Biosciences, 926-68073, 1:100, RRID: AB_10954442). Sections were counterstained with Hoechst 33342 (Sigma, 33342, 1:1000) to detect nuclei and Rhodamine labeled Dolichos Biflorus Agglutinin (DBA) (Vector Laboratories, RL-1032, 1:500, RRID: AB_2336396) was applied to visualize the ductal cells. Finally, sections were covered by ProLong Gold Antifade Mountant (Life Technologies, P36930) to ensure fluorescent signal preservation. Haematoxylin and eosin (Klinipath, 4085-9002) staining was performed on frozen tissue sections.

Frozen tissue section imaging

Frozen pancreatic tissue sections were imaged using an SP8-X confocal microscope (Leica). Whole pancreatic sections were scanned for Hoechst (405 nm laser), Alexa 488/ EYFP (488 nm laser), rhodamin red (553 nm laser), Alexa 647 (647 nm laser) and/or IR Dye 680 (670 nm laser).

Optical tissue clearing, immunostaining and imaging

Samples were cleared according to the iDisco protocol (Renier et al., 2014). Fresh pancreatic tissue was first fixed overnight in 4% paraformaldehyde. Next, the tissue was dehydrated with methanol/H2O series: 20%, 40%, 60%, 80%, 100%; 30min each, washed further in 100% methanol for 1hr and incubated overnight, with shaking, in 66% DCM/33% methanol at room temperature (DCM, Sigma 270997). Subsequently samples were washed in 100% methanol (2x), chilled at 4°C and bleached overnight in fresh 5% H2O2 in methanol at 4°C. Next, the tissue was rehydrated with methanol/H2O series: 80%, 60%, 40%, 20%, PBS; 30min each and washed twice in PBS/Triton-x 0.2% (PTx.2) for 1hr at room temperature. For immunolabeling tissue was incubated in permeabilization buffer (PBS/0.4% Triton X-100/0.3M Glycine/20% DMSO) for two days at 37°C, subsequently blocked in blocking buffer (PBS/0.2% Triton x-100/6% donkey serum/10% DMSO) for 4 days at 37°C and washed twice in PBS/0.2% Tween-20 with 10 μg/m Heparin (PTwH) (Sigma H3393) for 1hr at 37°C. Next, samples were incubated with primary antibody; chicken anti-GFP (Abcam, ab13970, 1:1000) and Rhodamine labeled Dolichos Biflorus Agglutinin (DBA) (Vector Laboratories, RL-1032, 1:500) in PTwH/5%DMSO/3% donkey serum for 6 days at 37°C, while slowly rotating. Then after washing several times in PTwH the tissue was incubated with the secondary antibody; Alexa 647 goat anti-chicken (Thermo Fisher, A21449, 1:100) in PTwH/3% donkey serum, slowly rotating, at 37°C for 4 days and again washed several times in PTwH. For clearing, the tissue was first dehydrated in methanol/H2O series: 20%, 40%, 60%, 80%, 100%, 100%; 1hr each at room temperature. Samples were left overnight in 100% methanol and the next day incubated for 3hr in 66% DCM/33% methanol at room temperature, 2x in 100% DCM for 15min and finally incubated in DiBenzyl Ether (DBE, Sigma, 108014) until samples were transparent. Samples were imaged using a light sheet microscope, Ultramicroscope II (LaVision BioTec). For image analysis Imaris 9.2 (Bitplane) was used.

Tissue isolation and preparation for Single-cell RNA-seq

To isolate single pancreatic cells for Single-cell RNA seq analysis freshly isolated pancreata of control (N = 5) and pancreatitis (N = 2), were first minced in small pieces and immediately transferred to 1ml cold digestion medium (0.5mg/ml Collagenase P, 10mM HEPES, 5mM glucose in HBSS). After processing all samples, the samples were incubated in a water bath at 37°C to acclimatize, followed by another 5min at 37°C for digestion. Next, 10ml cold 1% FCS DMEM was added to the samples to stop the digestion process. Cells were further dissociated by first pipetting them up and down with a P1000 pipet followed by filtering (70μm filter). After this, cells were centrifuged at 1200 RPM for 3min and the supernatant was removed to dissolve the pellet in 1ml cold 1%FCS DMEM. To collect single cell EYFP+ and EYFP- cells, the cells were sorted (FACSAria, BD) for EYFP expression (550-30nM) or no EYFP expression.

Single-cell RNA-seq library preparation and data analysis

Sorted EYFP+ and EYFP- cells were spun down at 720xg for 5min and supernatant was removed. Cells were diluted in 0.5% BSA/PBS to 500-2000 cells/μl and incubated on ice until further processing. Directly after this, libraries were prepared according to Chromium Next GEM Single Cell 3ʹ Reagent Kits v3.1 (PN-1000121). Briefly, Gel Beads-in-emulsion (GEMs) were generated by combining barcoded Single Cell 3′ v3.1 Gel Beads, a Master Mix containing cells, and Partitioning Oil onto Chromium Next GEM Chip. Immediately following GEM generation, the Gel Bead was dissolved, primers were released, and any co-partitioned cell was lysed. Primers were mixed with the cell lysate and a mastermix containing reverse transcription (RT) reagents. Next, Silane magnetic beads were used to purify the first-strand cDNA from the post GEM-RT reaction mixture. The full-length cDNA were amplified via PCR to generate sufficient mass for library construction. After end repair, A-tailing, adaptor ligation and PCR amplification, the final libraries were ready for sequencing (HiSeq).

Single cell sequencing data were aligned to the mouse reference genome (GRCm38) and processed using the CellRanger 3.1.0 software from 10x Genomics to generate unique molecular identifier (UMI) counts. The raw gene expression matrices were imported into R and further processed by the Seurat R package (Stuart et al., 2019) version 3.2.2 (filter < 20% of mitochondrial gene expression and > 200 unique gene counts (nFeature RNA) < 4000) and normalized by SCTransform (Hafemeister and Satija, 2019) with regression for nFeature RNA and the percent mitochondria. To correct the batch effect introduced by distinct batches during the sequencing process, we integrated and harmonized the datasets using the function FindIntegrationAnchors and IntegrateData from the Seurat package. Cell clusters were visualized using the UMAP algorithm with principal components as input. To perform function association between acinar EYFP- and EYFP+ cells, we conducted Gene Set Enrichment Analysis (GSEA) against the qvalue ranked list of genes in the experiment. Ranking was based on the qvalue generated by differentialGeneTest function from Monocle 2.18.0 (Qiu et al., 2017). Hallmark and Ontology (biological process) gene sets were downloaded from the Molecular Signatures Database (MSigDB) (Liberzon et al., 2011; Subramanian et al., 2005).

Quantification and statistical analysis

Image analysis

EYFP positive clone sizes were manually quantified on whole tissue sections imaged by confocal microscopy. The number of nuclei within an EYFP positive area was counted, the fraction of EYFP cells within an acinar structure (defined by a laminin staining) was determined and the cell type of EYFP positive cells was determined based on co-staining with either amylase (acinar cell), DBA (ductal cells) or negatively stained cells (islet cells). The total number of nuclei (Hoechst) and the total number of ductal cells (co-staining of Hoechst and DBA) of a whole tissue section was quantified automatically through intensity-based segmentation. The segmented areas were analyzed for overlap between different channels to determine the type of nuclei/cell. Nearest Euclidean distance between each EYFP cell and duct was calculated to classify as ductal or non-ductal EYFP. Noise reduction in images was carried out by eliminating isolated pixel groups smaller than area of 10μm2.

Quantification of Ki67+ cells and tissue+ area were performed using ImageJ software (Bethesda, MD) or Imaris (Bitplane). To determine the percentage of Ki67+ cells, the total number of nuclei (Hoechst) and Ki67+ cells of three representative areas of a whole tissue section (or a whole tissue section) per time point were determined with identical parameter settings for immunofluorescent intensity and particle size . The number of acini per mm2 (three representative images of 1 mm2 per mouse, N = 8 mice) and the nuclei within acini (80 acini per mouse, N = 8 mice) were manually quantified.

Quantification of the percentage of EYFP+ and EYFP- cells positive for STMN1, DCLK1, Nestin or Ki67 were performed manually. To quantify the total number of EYFP- cells, four representative areas of whole tissue sections per time point were determined with ImageJ software (Bethesda, MD) with identical settings for immunofluorescent intensity and particle size (N = 3 control mice and N = 3 pancreatitis mice). Imaris 9.2 (Bitplane) software was used to calculate the distance to a nearest duct for EYFP+ and EYFP- cells (N = 7 control mice and N = 4 pancreatitis mice).

Statistics and reproducibility

Statistical tests, sample sizes and definitions of error bars are indicated in the figure legends and calculated using Graphpad Prism 8 or MATLAB R2019a. All statistical tests were two-sided. A P value of < 0.05 was considered significant.

REAGENT or RESOURCE	SOURCE	IDENTIFIER
Antibodies
rabbit anti-amylase	Sigma	A8273; RRID: AB_258380
chicken anti-GFP	Abcam	ab13970; RRID: AB_300798
rabbit anti-Ki67	Sigma	SAB5500134; RRID: AB_2892217
rabbit anti-laminin	Novus Biologicals	NB300-144SS; RRID: AB_921870
Alexa 647 goat anti-rabbit	Thermo Fisher	A-21244; RRID: AB_2535812
Alexa 488 goat anti-chicken	Thermo Fisher	A11039; RRID: AB_2534096
IR Dye 680 donkey anti- rabbit	Li-Cor Biosciences	926-68073; RRID: AB_10954442
Rhodamine labeled Dolichos Biflorus Agglutinin	Vector Laboratories	RL-1032; RRID: AB_2336396
Alexa 647 goat anti-chicken	Thermo Fisher	A21449; RRID: AB_2535866
rabbit anti-Stathmin 1	Proteintech	11157-1-AP; RRID: AB_2197114
rabbit anti-DCLK1	Abcam	Ab31704; RRID: AB_873537)
mouse anti-nestin	Santa Cruz	10C2 sc-23927; RRID: AB_627994
rabbit anti-Ki67(SP6)	NeoMarkers	9106S807E; RRID:AB_149707
Chemicals, peptides, and recombinant proteins
ProLong Gold Antifade Mountant	Life Technologies	P36930
Haematoxylin and eosin	Klinipath	4085-9002
Hoechst	Sigma	33342
Ultra-V Block	Thermo Scientific	TA-125-PBQ
Normal antibody diluent	KliniPath, VWR International, Radnor, PA	ABB999
cerulein	Sigma-Aldrich	C9026
Metamizol	Sigma-Aldrich	46232
DCM	Sigma	270997
DiBenzyl Ether	Sigma	108014
Collagenase P	Merck/Roche	11213857001
Critical commercial assays
Chromium Next GEM Single Cell 3ʹ Reagent Kits v3.1	10x genomics	PN-1000121
Deposited data
Mathematical model code	This paper	https://zenodo.org/record/5017711
Single cell RNA sequencing	This paper	GEO: GSE171731
Experimental models: Organisms/strains
The R26[CA]30EYFPmouse strain	Doug Winton; Kozar et al., 2013	https://doi.org/10.1016/j.stem.2013.08.001. RRID: IMSR_JAX:006148
Software and algorithms
Imaris 9.2	Bitplane	N/A
ImageJ software	Bethesda, MD	https://imagej.nih.gov/ij/
Graphpad Prism 8	Graphpad Software	https://www.graphpad.com/prism
MATLAB R2019a	MathWorks	https://www.mathworks.com/products/matlab
Other
Ultramicroscope II	LaVision BioTec	https://www.lavisionbiotec.com/
SP8-X confocal microscope	Leica	https://www.leica-microsystems.com/