article: New approaches to lineage tracking have allowed the study of differentiation in multicellular organisms over many generations of cells. Understanding the phenotypic variability observed in these lineage trees requires new statistical methods. Whereas an invariant cell lineage, such as that for the nematode Caenorhabditis elegans, can be described by a lineage map, defined as the pattern of phenotypes overlaid onto the binary tree, a traditional lineage map is static and does not describe the variability inherent in the cell lineages of higher organisms. Here, we introduce lineage variability maps which describe the pattern of second-order variation in lineage trees. These maps can be undirected graphs of the partial correlations between every lineal position, or directed graphs showing the dynamics of bifurcated patterns in each subtree. We show how to infer these graphical models for lineages of any depth from sample sizes of only a few pedigrees. This required developing the generalized spectral analysis for a binary tree, the natural framework for describing tree-structured variation. When tested on pedigrees from C. elegans expressing a marker for pharyngeal differentiation potential, the variability maps recover essential features of the known lineage map. When applied to highly-variable pedigrees monitoring cell size in T lymphocytes, the maps show that most of the phenotype is set by the founder naive T cell. Lineage variability maps thus elevate the concept of the lineage map to the population level, addressing questions about the potency and dynamics of cell lineages and providing a way to quantify the progressive restriction of cell fate with increasing depth in the tree.
  The lineage map is a concept born from the study of invariant lineages, such as that for C. elegans, where the fixed pattern of phenotypes can be found, in principle, by tracking the family tree, or pedigree, from a single founder cell. When pedigrees are highly variable, however, seemingly identical founder cells can give rise to different patterns of descendants. Which of these represents the lineage map? Simply averaging the phenotype at each lineal position by pooling across multiple pedigrees may not give a representative pattern since the averaging will suppress correlations between lineal positions. Furthermore, the variability between pedigrees, which reflects the potency of founder cells, is an important quantity itself and cannot be represented in a lineage map. While lineage variability is minimal in simple organisms such as C. elegans [9,10] and leech [11], it is greater in higher organisms such as insects and vertebrates [1,12] and is significant in mammalian cells of clinical importance such as stem cells [13] and lymphocytes [14,15]. Such increased variability likely plays a role in being able to respond effectively to environmental changes. Given the additional variation inherent in molecular-level measurements [16] it is becoming increasingly important to extend the concept of the lineage map to account for variability. A fundamental property of any lineage tree measurement, which is crucial when studying variability, is whether it is ordered or unordered. We define a lineage tree to be ordered if the labels used to distinguish the lineal positions of two daughter cells (sisters) are meaningful. This gives each daughter a unique identity. The tree is unordered if these labels are arbitrary, making the daughter positions unidentifiable. This distinction has significant consequences in a statistical analysis. For example, in an unordered tree consisting of a mother A and its arbitrarily-labeled daughters B and C, we should not distinguish between the properties of positions B and C when comparing different pedigrees; nor should we distinguish between the mother-daughter relationships A-B and A-C. The C. elegans lineage is an example of an ordered lineage. Here each daughter is labeled by its orientation, at the time of division, with respect to the developing organism [4,5, 17]. Daughters thus have meaningful labels. For example, a mother labeled ‘Epl’ divides into ‘Epla’, the anterior daughter, and ‘Eplp’ the posterior daughter. In higher organisms, it can be difficult or impossible to find a meaningful way to label each daughter in a pair. For example, in the in vitro murine T-lymphocyte pedigrees discussed in this paper, the orientation of the mouse, even if it were known, can hardly be expected to be a useful way to distinguish between daughter cells. Now, if the phenotype pattern is invariant, a single complete pedigree measurement represents the lineage map and being unordered (or not) does not matter. However, when comparing variable pedigrees, the unidentifiability of daughter positions can significantly affect the ability to detect phenotypic patterns. For example, simply averaging different pedigrees to enhance subtle bifurcation patterns risks instead canceling these patterns if the pedigrees are not ordered the same way [18]. Since the majority of pedigree measurements from higher organisms are both variable and unordered [1] (see Table 1), an important question is whether the concept of a lineage map is even useful anymore. How can we associate fate specification with fixed lineal positions when the pattern of descendants varies from one apparently identical founder to the next? Clearly a statistical approach is required. A number of statistical methods have been developed to analyze variable, unordered lineage trees. Though these approaches do not directly address the question of how to construct a lineage map, many of them address aspects of the problem. The bifurcating autoregressive model [22,23] was developed to estimate mother-daughter and daughter-daughter correlations using a sample of unordered pedigrees from either E. coli or tumor cultures. The model was later used to analyze data from ordered pedigrees to test for lineage asymmetry [24,25]. This stationary, parametric model allowed for daughters to be conditionally dependent (with respect to their common mother) but forced cousins and more distant relatives to be conditionally independent (with respect to their most recent common ancestor). The subsequent discovery that cousins could be conditionally dependent motivated a theory of cellular inheritance involving deterministic, nonlinear dynamics in lymphoblasts [26]. However, such distant intragenerational dependence might instead be interpreted as a delay between fate specification and expression, where a phenotype that has been specified in a mother and its daughters is not expressed until its four granddaughters. These analyses illustrate the importance of having lineages that are large enough, and a model that is general enough, to examine correlations of distant relatives [27]. They also remind us that simple branching process models, which we here define to be those assuming conditional independence of daughters, do not properly represent the correlations in a lineage, a fact that was established in early lineage analyses [28,29]. Although population numbers can be satisfactorily modeled using branching processes [30], allowing for sibling correlations can affect the inferred population dynamics [14,31]. As we indicated earlier, identifying a subtree, or subclone, of shared phenotypes is the first step to inferring where fate is specified. This idea forms the basis of methods to study cell state transitions in bacterial cells or mouse embryonic stem cells [32,33], where phenotypic similarity among relatives in the same generation was used to infer how much earlier in the pedigree a transition occurred. A similar idea was used in hematopoietic stem cells to infer the multi-generational delay between when an invisible molecular decision occurred and when its effect was expressed as a surface marker [34]. These models assume that cell states transition over timescales that are slow compared to the cell cycle duration; alternatively they could be synchronized to cell divisions [35]. Note that, in a lineage map, the generation of a cell is a meaningful quantity, representing the number of divisions since the founder cell, whether that be a zygote, a naive lymphocyte, or some progenitor initiated with a particular stimulus. Thus any model of a developing lineage must be non-stationary. Several other approaches to statistical lineage analysis have been reported recently. A factor graph method was used to model conditional dependence between daughters [36], with the goal of testing whether pre-programmed instruction or differential cell death was responsible for differentiation of hematopoietic progenitor cells; direct inference of Nanog expression, a pluripotency factor, was used to understand its dynamics in embryonic stem cell lineages [37]; and, a parametric characterization of lineage patterns has been applied to achieve early identification of hematopoietic stem cells [38]. However, these methods are of less relevance to our question of how to build a statistical lineage map. Major efforts are underway to improve the throughput and quality of lineage measurements (see reviews [13,39–42] and commentary [43–45]). Recent breakthroughs have resulted in a wealth of data from automated microscopy-based [19,46–48] and sequencing-based [20,49–57] techniques. While the technological barriers for these measurements are severe, there are significant barriers to the analysis of the data as well. As we have discussed, there is currently no way to construct a useful lineage map from variable, unordered pedigrees. Yet, according to Shapiro [40], “Central unresolved problems in human biology and medicine are in fact questions about the human cell lineage tree: its structure, dynamics, and variability during development, growth, renewal, aging and disease. ” It is important then to find a way to generalize the concept of the lineage map to the population level. In this paper we propose a solution by introducing lineage variability maps. These involve the variances of, and covariances between, every position in the tree. A key idea is that, to interpret lineage patterns, it is not only the phenotypic values at each lineal position that are important, but also the phenotypic associations between different lineal positions. By developing a generalized spectral analysis for binary trees, we show how to estimate variability maps for a lineage of any depth using measurements from only a few, unordered pedigrees. For complete data, our approach is a non-parametric one, involving first and second moments of the data but assuming no distribution function. We could thus, alternatively, refer to these maps as second order lineage maps. The rest of the paper is organized as follows. In the Methods section we describe essential aspects of the data used in this paper, the analysis framework and labeling convention, and how to estimate all pairwise associations by employing constraints on symmetry and sparsity. In the Results section we interpret the covariance matrix in terms of graphical models, called lineage variability maps, and show how these can be used to understand fate restriction and expression throughout the lineage. In the Discussion section we examine the implications and future prospects of this analysis.
  All experiments using mice were performed in accordance with the Animal Experimentation Ethics Committee of the Peter MacCallum Cancer Centre (Approval E427) and mice were sacrificed by anaesthetic overdose. Data from 3 types of lineages are analyzed in this paper. Experimental data from T cells provide an example of a lineage with significant variability and no obvious structure. Previously-published data from C. elegans are the example of a lineage with complicated but highly-reproducible structure. Finally, a simulated, stationary branching process provides the benchmark of a featureless, variable lineage and to test the accuracy of the inference procedure. In more detail: Sample pedigrees from these 3 lineage types are shown in Fig 1 while the expression of each phenotype as a function of generation is shown in Fig 2. In this study, lineage data are regarded as repeated measurements on pedigrees arising from individual founder cells, each selected at random from a population of similar cells. We restrict our attention to modeling a single trait from pedigrees subject to the same conditions. A sample consisting of multiple replicate pedigrees can then be represented by a two-factor array (Yij), where i has n levels corresponding to the number of pedigrees and j has p levels corresponding to the number of lineal positions within a pedigree. With no meaningful distinctions among pedigrees (they are all of the same cell type and subject to the same conditions) we assume they are independent and identically-distributed replicates. The data can thus be represented by a matrix Y with n replicates (rows) and p variables (columns). Each of the p dimensions corresponds to a lineal position. We use a binary number to label each position so that, for example, the first 3 generations are labeled as founder (1), daughters (10,11), and granddaughters (100,101,110,111), where each label thus encodes the lineal position. We will also need to define a nomenclature for generations and subtrees. Generations, g, refer to the depth in the tree where the founder cell is defined to be at generation g = 1. Subtrees are defined by two indices, (ℓ, τ), where ℓ refers to the longitudinal coordinate and τ to the transverse coordinate of the root node (see Fig 3). In our convention, the subtree at ℓ = 1 is the entire tree. As we will show, subtrees will be associated with sources of variation, and we need to define a ‘subtree’ at ℓ = 0 that sits outside the lineage and represents variation among lineages. This concept does not exist for a lineage map but is essential for a lineage variability map since different pedigrees may not be the same. Often in lineage measurements there are many more lineal positions (p) than there are pedigrees (n). Thus p ≳ n, with the disparity getting exponentially worse with the number of generations studied. Performing reliable inference when p/n > 1 is an open research question [61]. Best results are achieved when prior knowledge of the problem can be incorporated. In the next section we describe increasingly more sophisticated steps to reduce the effective dimensionality of the inference calculation, first by exploiting known symmetry properties and then by using observed sparsity properties. Our goal is to identify a scheme where the number of replicates required is independent of the number of generations studied. This will allow variability maps over many generations to be constructed using data from only a few pedigrees. The essential idea for this analysis is to measure second-order variation throughout the lineage by estimating the variance of, and covariance between, every lineal position. This population covariance matrix Σ for the lineage involves no assumption about the underlying distribution. It involves just the first and second order moments of the data.
  We can visualize Σ ^ G and Σ ^ Ω using graphical models to produce different ‘maps’ of the variation in the lineage. We call these lineage variability maps. To depict Σ ^ G we use undirected graphs, since lineal positions within a generation have no ordering. We call the result a lineage correlation map. To depict Σ ^ Ω we can use directed graphs, since the natural variables in an irreducible subspace are ordered in a sequence. Thus the spectral transformation enables the undirected graph to be converted into a directed one. This graph, which we call a dynamic lineage map, compactly represents the dynamics of the bifurcated expression pattern in each subtree. Lineage variability maps describe the pattern of phenotypic associations throughout the lineage. However, as with lineage maps, our interest is often in using them to infer where fate is specified. In the introduction, we described how this involves identifying the most recent common ancestor of cells with shared fate. For a clonal pattern, where a cell fate is exclusive to a single subtree, we would infer that fate was specified at (or near) a single lineal position—the root of that subtree. For a non-clonal pattern, where cell fate is expressed in multiple subtrees, we would infer that fate was specified at multiple lineal positions. In C. elegans these inferences can be made visually [5]. Here we show how, by knowing the lineage variability map Σ, we can make these inferences statistically, overcoming the problem of how to identify subtrees with shared phenotypes when lineages are variable. Before we begin, we must define what we mean by cell fate. In this study we define cell fate to be the measured phenotype of a cell at the latest generation studied, G. This practical definition allows us to analyze cell fate whether or not the phenotype in the last generation is actually a terminal fate. Also, by defining cell fate as the phenotype itself rather than as a cell type (determined by that phenotype), we can use the phenotypic measurements as is, without having to cluster or threshold them. Such discretization procedures can be difficult to define when phenotypes exist on a continuum of differentiation, as is often the case [79]. Having defined fate, we turn now to explaining its variability in terms of aspects of the lineage. We first partition the variability among the subtrees, or sources of variation. This quantifies how much of a cell’s fate is restricted by, or specified by, each subtree. We then examine the correlation of a cell’s fate with the phenotypes of its ancestors. This identifies the earlier generations over which a phenotypic fate has been stably expressed. Together these two measures, of fate restriction and fate expression, make up what we call fate profiles.
  The lineage map, which has been instrumental in the discovery of fate specification mechanisms in simple organisms, was born from the study of invariant lineages and is not a particularly useful concept for understanding the more ubiquitous case of variable lineages. To address this limitation, we have introduced lineage variability maps, which provide a way to describe lineages at the population level. Whereas the lineage map is a description of the fixed pattern of phenotypes across a pedigree, the lineage variability map describes the pattern of phenotypic associations across a pedigree. This map of phenotypic associations, Σ, provides quantitative answers to essential scientific questions such as those about cell potency, fate restriction, and the sources of variation in a lineage. We have constructed lineage variability maps from a sample of highly-variable pedigrees from CD8+ T-lymphocytes up to five generations. These show that most of the variation in cell fate, defined here to be average cell size at generation 5, is explained by the choice of naive cell. Yet, despite the pivotal role played by this founder in restricting cell fate, its phenotype is not predictive of fate: though a naive cell may specify that its descendants be large, it may not be large itself (at least on average). Although we expect to apply our technique primarily to variable pedigrees which are difficult to interpret by visualization alone, we can also apply it to invariant lineages to check our results. In fact, by constructing lineage variability maps from sample wild-type pedigrees from C. elegans marked for pharyngeal expression, we successfully recovered essential information in the known lineage map, identifying global features such as the small degree of inter-pedigree variation characteristic of a totipotent zygote, and the several-generation delay between fate specification and expression. Yet our lineage variability maps capture important finer detail as well. Consider the peak in fate restriction at ℓ = 2 observed in Fig 11b. This arises from the strong bifurcation of fate traced back to the division of both P1 and of AB, progenitors located at ℓ = 2 (see Section S1. 6 in S1 Appendix for the labeling of lineal positions in the worm). That only a single daughter from P1 and from AB exhibit pharyngeal fate results in the spike in fate restriction that we observe. Interestingly, this phenomenon, of pharyngeal fate ensuing from two cousins at generation 3 (ABa and EMS) but not from their sisters, is a phenomenon that has been investigated in detail [80]. Such work laid the foundation for further studies leading to the understanding of the molecular and cellular mechanisms for specification of pharyngeal tissue [81]. This demonstrates how, even though we may be ignorant of the ordering of the lineage, we can still detect a fate bifurcation event of biological relevance that had previously required knowledge of this ordering. In other words, although we must assume lineage relationships are symmetric, this does not prevent us from detecting the effects of asymmetric lineage patterns from the ‘boost’ they give to the variance in particular subtrees. Recent technological innovations have introduced a variety of methods for recording lineage data, involving both advanced imaging [19,46–48] and genetic barcoding [20,49,51–57] techniques. With the statistical lineage mapping and fate profiling methods described in this manuscript, it should be possible to quantify fundamental properties of these lineages, such as the potency of progenitors, whether heterogeneity is clonal, and at what depth such heterogeneity appears. Just as the visual identification of fate bifurcations in the worm lineage map enabled the location of fate specification events to be discovered, the capacity to perform systematic screens to rapidly identify the important stages of fate restriction should contribute to a deeper understanding of the mechanisms of fate specification in more complex, more variable systems.
summary: Multicellular organisms develop from a single fertilized egg by sequential cell divisions. The progeny from these divisions adopt different traits that are transmitted and modified through many generations. By tracking how cell traits change with each successive cell division throughout the family, or lineage, tree, it has been possible to understand where and how these modifications are controlled at the single-cell level. This helps address questions about, for example, the developmental origin of tissues, the sources of differentiation in immune cells, or the relationship between primary tumors and metastases. Such lineages often show large variability, with apparently similar founder cells giving rise to different patterns of descendants. In addition, questions about the range of accessible cell types at different stages of the lineage tree are actually questions about lineage variability. To characterize this variation, and thus understand the lineage at the population level, we introduce lineage variability maps. Using data from worm and mammalian cell lineages we show how these maps provide quantifiable answers to questions about any developing lineage, such as the potency of progenitor cells and the restriction of cell fate at different stages of the tree.
keywords: blood cells
  infographics
  invertebrates
  fourier analysis
  medicine and health sciences
  immune cells
  caenorhabditis
  immunology
  random variables
  covariance
  animals
  cell differentiation
  animal models
  developmental biology
  caenorhabditis elegans
  model organisms
  mathematics
  experimental organism systems
  discrete mathematics
  combinatorics
  research and analysis methods
  computer and information sciences
  white blood cells
  animal cells
  animal studies
  mathematical and statistical techniques
  t cells
  probability theory
  eukaryota
  permutation
  cell biology
  data visualization
  phenotypes
  graphs
  genetics
  nematoda
  biology and life sciences
  cellular types
  physical sciences
  organisms
title: Maps of variability in cell lineage trees