A composite likelihood approach to the analysis of longitudinal clonal data on multitype cellular systems under an age-dependent branching process

doi:10.1093/biostatistics/kxq050

Biostatistics. 2011 January; 12(1): 173–191.

Published online 2010 August 23. doi: 10.1093/biostatistics/kxq050

PMCID: PMC3006127

A composite likelihood approach to the analysis of longitudinal clonal data on multitype cellular systems under an age-dependent branching process

Rui Chen and Ollivier Hyrien^*

Department of Biostatistics and Computational Biology, University of Rochester Medical Center, 601 Elmwood Avenue, Rochester, NY 14642, USA ; Email: ollivier_hyrien/at/urmc.rochester.edu

Mark Noble and Margot Mayer-Pröschel

Department of Biomedical Genetics, University of Rochester Medical Center, 601 Elmwood Avenue, Rochester, NY 14642, USA

^*To whom correspondence should be addressed.

Received August 26, 2008; Revised June 3, 2010; Accepted June 29, 2010.

This article has been cited by other articles in PMC.

Abstract

A recurrent statistical problem in cell biology is to draw inference about cell kinetics from observations collected at discrete time points. We investigate this problem when multiple cell clones are observed longitudinally over time. The theory of age-dependent branching processes provides an appealing framework for the quantitative analysis of such data. Likelihood inference being difficult in this context, we propose an alternative composite likelihood approach, where the estimation function is defined from the marginal or conditional distributions of the number of cells of each observable cell type. These distributions have generally no closed-form expressions but they can be approximated using simulations. We construct a bias-corrected version of the estimating function, which also offers computational advantages. Two algorithms are discussed to compute parameter estimates. Large sample properties of the estimator are presented. The performance of the proposed method in finite samples is investigated in simulation studies. An application to the analysis of the generation of oligodendrocytes from oligodendrocyte type-2 astrocyte progenitor cells cultured in vitro reveals the effect of neurothrophin-3 on these cells. Our work demonstrates also that the proposed approach outperforms the existing ones.

Keywords: Bias correction, Cell differentiation, Composite likelihood, Discrete data, Monte Carlo, Neurotrophin-3, Oligodendrocytes, Precursor cell, Stochastic model

1. INTRODUCTION

In stem cell biology, there exists considerable interest in studying signals that may modulate or alter the processes that regulate the formation of tissues during development or repair. Understanding such processes has important clinical implications as it may offer potential means to restore or maintain tissue and organ functions. One of the challenges in conducting detailed analysis of the effects of these signals, however, is the general lack of quantitative approaches that allow analysis of different aspects of progenitor or stem cell function.

The present paper is focused on oligodendrocyte type-2 astrocyte (O-2A) progenitor cells, which is one of the best studied precursor cell populations (e.g. Noble and others, 2004). These progenitor cells may either divide into 2 new progenitor cells or differentiate into terminally differentiated oligodendrocytes, which produce the myelin sheaths that enwrap axons in the central nervous system. The differentiation of progenitor cells into oligodendrocytes needs to be properly regulated to maintain normal function of the central nervous system. In addition to the function of O-2A progenitor cells in vivo, their ability to grow and differentiate in vitro has established these cells as an important tool in modern cell biology that allows the study of the most basic cellular functions (i.e. self-renewal by division, differentiation, and death) at the clonal level in tissue culture (e.g. Ibarrola and others (1996); Smith and others, 2000). The development of a typical clone of these cells cultured in vitro started from a single O-2A progenitor cell is displayed in Figure 1. It may be altered in multiple ways following exposure to agents, such as cytokines, drugs, or toxicants, and it is of great interest to biologists to determine and quantify these changes (Dietrich and others, 2006).

Fig. 1.

Oligodendrocyte development under in vitro conditions. When cultured under appropriate growing conditions, any O-2A progenitor cell will either divide into 2 new O-2A progenitor cells or it will differentiate into a single terminally differentiated oligodendrocyte. (more ...)

Over the past decade, a great deal of attention has been paid to the development of a statistical framework for gaining a quantitative understanding of the biological processes that govern the division of O-2A progenitor cells and their differentiation into oligodendrocytes (Yakovlev, Boucher, and others, 1998, Yakovlev, Mayer-Pröschel, and others, 1998, Yakovlev and others, 2000; von Collani and others, 1999; Boucher and others, 1999, Boucher and others, 2001; Zorin and others, 2000; Hyrien and others, 2005a, Hyrien and others, 2005b, Hyrien and others, 2006; Hyrien, 2007). These publications used multitype age-dependent branching processes as a means to model the temporal development of cell clones. The resulting methods provided quantitative insights into the generation of oligodendrocytes from cultured O-2A progenitor cells. For instance, Yakovlev, Boucher, and others (1998) observed that the generation of oligodendrocytes was regulated by a combination of cell-intrinsic factors and environmental signals that may modulate the probability of differentiation of O-2A progenitor cells. Other studies based on similar quantitative approaches (Yakovlev, Boucher, and others, 1998, Yakovlev, Mayer-Pröschel, and others, 1998, Yakovlev and others, 2000; von Collani and others, 1999; Boucher and others, 1999, Boucher and others, 2001; Zorin and others, 2000; Hyrien and others, 2006) confirmed that thyroid hormone induces the differentiation of O-2A progenitor cells into oligodendrocytes (Barres and others, 1994; Ibarrola and others (1996); Ahlgren and others, 1997; Rodríguez-Peña, 1999), and that the decision of O-2A progenitor cells to either divide or differentiate is made early in the cell cycle, if not before, rather than at the end (Hyrien and others, 2005a, Hyrien and others, 2006, Hyrien and others, 2010).

Most of these previous studies were concerned either with clonal experiments in which every cellular clone was observed only once in the course of the experiment thereby yielding independent observations (Yakovlev, Boucher, and others, 1998, Yakovlev, Mayer-Pröschel, and others, 1998, Yakovlev and others, 2000; von Collani and others, 1999; Boucher and others, 1999, Boucher and others, 2001; Zorin and others, 2000; Hyrien and others, 2005a, Hyrien and others, 2005b, Hyrien and others, 2006, Hyrien and others, 2010; Hyrien, 2007) or with time-lapse experiments (Hyrien and others, 2006), where the complete history of each clone was fully recorded. While time-lapse experiments provide complete information about the temporal development of each clone, clonal experiments are far less time-consuming. This is the reason why clonal analyses remain the experiment of choice to investigate the ability of agents to impact the regulation of precursor cell functions. There exists a third experimental setting in which each cell clone is observed longitudinally over time and which offers a compromise between the above 2 experiments. In this experimental setting, the family trees are not completely observed, but several snapshots of the temporal development of each individual clone are obtained at multiple time points thereby providing more information than independent clonal experiments. We shall refer to these alternate experiments as longitudinal clonal experiments (as opposed to “independent” clonal experiments where each cell clone is examined only once).

Longitudinal clonal data are not independent across time points, and their analysis presents new challenges. Some of the estimation methods that have been proposed for independent clonal data could also apply to the longitudinal setting. Since maximum likelihood estimation remains generally unattainable with age-dependent branching processes, simulation-based approaches have been adopted by some authors as a means to construct viable estimators. For instance, Hyrien and others (2005a), Hyrien and others (2005b) and Hyrien (2007) considered the utility of a simulated pseudolikelihood approach. This method relies solely on the expectation and variance–covariance functions of the process. It is attractive because it provides consistent and asymptotically Gaussian estimators, while remaining easy to implement, even for longitudinal clonal data. The disadvantage of this estimator is that it does not enjoy optimality properties, and, for instance, it is not as efficient as the method of maximum likelihood (Hyrien, 2007). The simulated maximum likelihood estimator proposed by Zorin and others, 2000 could be more efficient than the simulated pseudo maximum likelihood estimator but its implementation faces serious limitations that are amplified when cell clones are observed repeatedly over time. Specifically, it is very time-consuming, and it suffers from mismatches, as described by Zorin and others (2000). We shall discuss these limitations in greater details in Section 4 when motivating the construction of the proposed estimator.

This article proposes an alternative approach that provides a compromise between the statistical efficiency of the maximum likelihood estimator and the computational advantage of the simulated maximum pseudolikelihood estimator. The proposed method is based on a composite likelihood approach (sometimes also referred to as method of pseudolikelihood). The use of composite likelihood traces back to Besag (1974). It has been investigated in other settings by several authors, including Azzalini (1983), Lindsay (1988), Heagerty and Lele (1998), Cox and Reid (1987), Chandler and Bate (2007), Varin and Vidoni (2005), Hyrien and Zand (2008), to name a few. Composite likelihoods can be defined as estimating functions formed by combining together likelihood objects that remain tractable for the problem at hands. Our work shows that composite likelihood estimators provide a viable solution to the analysis of longitudinal clonal data.

The proposed method relies on a complex stochastic process. Simpler statistical approaches could be invoked to describe the time course of cell counts and assess whether specific treatment conditions alter the kinetics of the cell population. The level of sophistication of the proposed method is required for gaining a quantitative insight into the processes of division and differentiation of precursor cells because the events of interest (division and differentiation) are not observable during clonal experiments.

The remainder of this article is organized as follows. Section 2 introduces the data structure. Section 3 presents a branching process model of the generation of oligodendrocytes from cultured O-2A progenitor cells. Section 4 describes simulated composite likelihood estimation as applied to longitudinal clonal data. The finite-sample properties of the proposed method are investigated via simulations in Section 5. Section 6 presents an application to the analysis of the proliferation of O-2A progenitor cells and their differentiation into oligodendrocytes exposed to neurotrophin-3 (NT-3) in tissue culture.

2. LONGITUDINAL CLONAL DATA

2.1. Clonal data on the generation of oligodendrocytes in vitro

A typical longitudinal clonal experiment conducted on oligodendrocytes begins by plating initiator O-2A progenitor cells in culture dishes at a density that is low enough so these cells (and the pools of subsequent progenies) will not interact with each other. Over time, these progenitor cells divide into O-2A progenitor cells and/or differentiate into oligodendrocytes to give rise to individual cell clones. These clones may contain O-2A progenitor cells or oligodendrocytes, or, as is more typically the case, a mixture of both cell types. O-2A progenitor cells and oligodendrocytes are morphologically distinguishable, which enables experimentalists to count separately the number of O-2A progenitor cells and the number of oligodendrocytes contained in any given clone through visual examination using a microscope. The location of each clone in the dish is identified for subsequent follow up, and the composition of each clone is examined repeatedly over time so one can observe how the numbers of O-2A progenitor cells and the number of oligodendrocytes change over time.

2.2. An example of longitudinal clonal data

We performed a longitudinal clonal experiment, as described above, to investigate the impact of NT-3 on the processes of division and differentiation of O-2A progenitor cells. We purified O-2A progenitor cells isolated from the corpus callosum tissue of 7-day old rats as described previously (Ibarrola and others, 1996) and plated cells at a clonal density in defined medium supplemented with NT-3 at 20 ng/ml. No platelet-derived growth factor (PDGF) was added to the culture medium at any time. In this experiment, n = 40 clones were followed for up to 6 days. Every clone was generated by a single O-2A progenitor cell. The composition of each clone was examined daily, so t_i = (1,2,3,4,5,6) in this particular experiment (with time being expressed in days). The number of O-2A progenitor cells and the number of oligodendrocytes were reported at each time point for every clone. Half of these clones were cultured in the presence of NT-3, and the other half was cultured without NT-3.

Figure 2 shows the histograms for the number of progenitor cells (left panels) and for the number of oligodendrocytes (right panels) from day 0 to day 6 (from top to bottom). At the start of the experiment (day 0), all clones contained exactly one O-2A progenitor cell and zero oligodendrocyte. The composition of each clone evolved over time according to whether O-2A progenitor cells and their subsequent progenies divided or differentiated. In this particular experiment, all clones growing without NT-3 contained between 0 and 6 O-2A progenitor cells and between 0 and 4 oligodendrocytes. In the presence of NT-3, the number of O-2A progenitor cells ranged between 0 and 9, and the number of oligodendrocytes between 0 and 6.

Fig. 2.

Histograms for the number of O-2A progenitor cells (PCs) and for the number of oligodendrocytes along with the fitted model in the absence and in the presence of NT-3.

The primary objectives of our experiment were to assess and investigate the effects of NT-3 on the proliferation of O2-A progenitor cells and their differentiation into oligodendrocytes. A visual comparison of the histograms of the number of O2-A progenitor cells and of the number of oligodendrocytes clearly suggests that the numbers of O-2A progenitor cells were substantially smaller in clones treated with NT-3 starting from day 2. In contrast, the numbers of oligodendrocytes appeared similar among treated and untreated clones at all time points. The marginal distributions of the number of O-2A progenitor cells (respectively, oligodendrocytes) at any time point in clones treated with and in clones treated without NT-3 could be compared using (e.g.) Wilcoxon rank sum statistics. In order to assess the overall effect of NT-3, irrespective of time, and avoid issues associated with multiple testing, the resulting p-values could be combined using Fisher's combination method. Since the individual p-values are dependent across time points, an overall p-value could be computed using a permutation testing approach, where clones (not just individual observations) are randomly assigned to one group or the other so the dependencies among observations are properly accounted for when performing the test.

We assessed the overall effect of NT-3 on the numbers of O-2A progenitor cells and on the numbers of oligodendrocytes separately using this approach. Our test suggested that the numbers of O-2A progenitor cells were significantly different depending upon whether they had been cultured with or without NT-3 (p = 0.02), but it did not detect any significant difference among the numbers of oligodendrocytes (p = 0.98). Thus, this analysis would suggest that NT-3 had a significant impact on the number of O-2A progenitor cells but not on the number of oligodendrocytes, at least for the first 6 days of culture.

The conclusion of these analyses is of interest in itself because it reveals the existence of a potential effect of NT-3 on the regulation of the processes of division and of differentiation of O-2A progenitor cells. From a biological standpoint, however, it remains limited in scope because it neither offers any mechanistic interpretation on how these processes might have actually been altered by NT-3 nor quantifies the effects of NT-3 on cellular functions. A number of biological reasons could be invoked to explain the observed effects of NT-3. For instance, O-2A progenitor cells exposed to NT-3 might have differentiated more frequently in the presence of NT-3, causing the number of O-2A progenitor cells to decrease prematurely. An alternative explanation is that O-2A progenitor cells divided more slowly and/or differentiated faster when exposed to NT-3. A combination of these 2 scenarios is also quite plausible. In addition, experimentalists would also be interested in evaluating the impact of such changes on the number of oligodendrocytes produced by O-2A progenitor cells over an extended period of time. In order to gain such a quantitative insight into cell proliferation kinetics, we propose to model the dynamics of cell clones using an age-dependent branching process.

3. A BRANCHING PROCESS MODEL OF OLIGODENDROCYTE GENERATION

We model the temporal development of the population of O-2A progenitor cells and of terminally differentiated oligodendrocytes using a multitype age-dependent branching process defined from the following set of assumptions:

(A1)The process begins at t = 0 with a single cell of type 1 (an initiator O-2A progenitor cell).
(A2)For every k = 1,2,…, upon completion of its lifespan, every type-k cell (an O-2A progenitor cell of generation k) either divides into 2 new cells of age 0 and type k + 1 with probability π_k or it differentiates into a single type-0 cell (an oligodendrocyte) with probability 1 − π_k. Following Boucher and others (2001) and Hyrien and others (2005a), Hyrien and others (2005b), we assumed that the probability of division is given by π_k = min{1,a + bc^k}, k ≥ 1, where a, b, and c are unknown positive constants with, for instance, a denoting the limiting probability of division of O-2A progenitor cells as the number of division increases when c < 1.
(A3)Time to division of initiator O-2A progenitor cells: The time to division of any type-1 cell is represented by a nonnegative random variable (r.v.) T₁ assumed to follow a gamma distribution with cumulative distribution function (c.d.f.) G₁(t;ξ₁), where ξ₁ = (m₁,σ₁), and where m₁ and σ₁ denote the expectation and standard deviation of T₁.
(A4)Time to division of O-2A progenitor cells of generation greater than or equal to 2: The time to division of any O-2A progenitor cells of generation k, k ≥ 2, is described by a nonnegative r.v. T₂ assumed to follow a gamma distribution with c.d.f. G₂(t;ξ₂), where ξ₂ = (m₂,σ₂), and where m₂ and σ₂ denote the expectation and standard deviation of T₂.
(A5)Time to differentiation of O-2A progenitor cells: The time to division of any O-2A progenitor cells of generation k, k ≥ 1, is described by a nonnegative r.v. T₀ assumed to follow a gamma distribution with c.d.f. G₀(t;ξ₀), where ξ₀ = (m₀,σ₀), and where m₀ and σ₀ denote the expectation and standard deviation of T₀. The distribution of the time to differentiation is thus independent of the generation number.
(A6)Independence assumption: Cells evolve independently from each other.

The above process extends a model of the generation of oligodendrocytes previously proposed by Hyrien and others (2005a) by allowing the distribution of the time to division of O-2A progenitor cells of first generation to be dissimilar from that of O-2A progenitor cells of subsequent generations. This extension was motivated by an analysis of time-lapse data on oligodendrocytes (Hyrien and others, 2006), which indicated that first generation O-2A progenitor cells may have a shorter time to division than cells of subsequent generations.

Let θ = (m₀,σ₀,m₁,σ₁,m₂,σ₂,a,b,c)^′ denote the complete set of unknown model parameters. For every k = 0,1,…, let Z_k(t) denote the number of type-k cells at time t ≥ 0, with Z₀(t) corresponding to the number of oligodendrocytes, and, for every k ≥ 1, Z_k(t) representing the number of O-2A progenitor cells of generation k. It is not possible to determine the generation of O-2A progenitor cells during clonal experiments, and our observations consist of the number of O-2A progenitor cells Y₁(t) = ∑_{k = 1}^∞Z_k(t) and of the number of oligodendrocytes Y₂(t) = Z₀(t). The observable process is therefore the vector Y(t) = {Y₁(t),Y₂(t)}. The distribution of the process Y(t) does not generally admit a closed-form expression, which is why we and others have resorted to simulations to approximate this distribution.

4. STATISTICAL INFERENCE

4.1. A composite likelihood estimator

Let Y_ij1 and Y_ij2 denote the number of O-2A progenitor cells and the number of oligodendrocytes counted in the ith clone at time t_ij, where j = 1,…,n_i and i = 1,…,n. Write Y_ij = (Y_ij1,Y_ij2) and Y_i = (Y_i1,…,Y_{in_i}) for the ordered time points t_i = (t_i1,…,t_{in_i}). The log-likelihood function for the sample Y₁,…,Y_n takes the form

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx1_ht.jpg

where p_ij(θ) = _θ{Y(t_ij) = Y_ij|Y(t_{i,j − 1}) = Y_{i,j − 1},…,Y(t_i,1) = Y_i,1} denote the conditional distribution function of the r.v. Y(t_ij), given the trajectory of the process Y(t) is observed at all sampling times prior to t_ij.

In the present setting, the p_ij(θ)s have usually no closed-form expressions, making (exact) maximum likelihood inference impossible. Zorin and others (2000) proposed a simulated maximum likelihood approach when a single observation is available per clone (n_i = 1). The probabilities p_i1(θ)s are approximated by simulating independent clones from the branching process, and by computing the empirical proportions of simulated clones whose composition at time t_i1 matches that of the ith experimental clone. Denoting the corresponding approximation by An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx3_ht.jpg

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx3_ht.jpg

, L(θ) is replaced by An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx2_ht.jpg

, and a simulated maximum likelihood estimate is computed using a multidimensional version of the Kiefer–Wolfowitz algorithm (Kiefer and Wolfowitz, 1952; Blum, 1954).

There are 2 limitations to such an approach. The first one is that An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx3_ht.jpg

will be identically zero if the composition of the ith clone as observed during the actual experiment cannot be matched by any of the simulated clones. Consequently, An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx4_ht.jpg

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx4_ht.jpg

is undefined, making it difficult to compute a simulated maximum likelihood estimate. This problem was reported by Zorin and others (2000) in their analysis of oligodendrocytes generation. The second limitation has to do with the fact that the simulated log-likelihood function An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx4_ht.jpg

is a biased estimator of L(θ). The Kiefer–Wolfowitz algorithm considered by these authors is designed to converge to the parameter vector that maximizes the map An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx5_ht.jpg

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx5_ht.jpg

, where the expectation is taken with respect to the simulations. Since the maximum of g(θ) will be generally different than the maximum likelihood estimate, the simulated maximum likelihood estimator will be biased. The severity of the bias depends, of course, upon the noise of the simulated log-likelihood that is attributable to the simulations and upon the gradient of L(θ) with respect to θ. One can argue that the bias will become small if the number of simulations S is large enough. Computing time, however, may become an issue with branching processes, especially when the process is supercritical (i.e. loosely speaking, when each cell produces on average more than one cell), which causes cell clones to grow large; see Hyrien and others (2010) for additional comments on this practical issue. Thus, the strategy that consists of letting S increase is not always applicable in this case. Furthermore, there exists a simple way to approximate L(θ) that is able to mitigate both limitations simultaneously (see Section 4.2).

The simulated maximum likelihood estimator proposed by Zorin and others (2000) can be extended in an obvious way to the situation where clones are observed longitudinally over time (n_i > 1). However, the method becomes less appealing than when n_i = 1 because now simulated clones need to match the trajectories of experimental clones observed at several time points, not just one. Since the probability of mismatching increases quickly with n_i, the number of simulations would have to increase with n_i thereby making the approach even more time-consuming.

As an alternative to the method of maximum likelihood and simulated maximum likelihood, we propose to construct estimators of θ using the estimating function

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx6_ht.jpg

where q_ij(θ) = _θ{Y(t_ij) = Y_ij|Y(t_{i,j − k}) = Y_{i,j − k},k [set membership]

I_ij}, and where I_ij [subset or is implied by]

{1,…,j − 1} for every i and j. The probability q_ij(θ) differs from p_ij(θ) in that the conditioning event is restricted to a subset of the past observations. By selecting I_ij appropriately, we can thus reduce the curse of dimensionality encountered with the simulated log-likelihood. These probabilities q_ij(θ) do not generally have explicit expressions either, but we shall leave aside this aspect of the method until Section 4.2. The estimating function P(θ) is a composite likelihood, and the maximum composite likelihood estimator maximizes P(θ). We shall denote it by An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx7_ht.jpg

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx7_ht.jpg

The conditioning sets I_ij need not be identical for all i and j nor do they have to be collections of consecutive integers, although they will likely be so in most applications. Thus, I_ij = {1,3,5} would be allowed for instance. They are assumed to be independent of θ, however. Depending upon the choice of I_ij, we obtain different estimating functions, the complexity of which (from a computational standpoint) increases as I_ij “approaches” {1,…,j − 1}. The following cases are the most noteworthy:

If I_ij = {1,…,j − 1}, P(θ) = L(θ), and is a maximum likelihood estimator.
If I_ij = , no conditioning (other than upon the state of the process at the start of the experiment) is done. This is the simplest case. Hyrien and Zand (2008) used it in the context of CFSE-labeling data.
If I_ij = {1}, all probabilities are conditioned upon the most recent observation. The method leads to a maximum likelihood estimator when the process is Markovian and no efficiency is lost in this particular case. For non-Markovian models, this property does not hold true since the conditional distribution of Y(t_ij) given {Y(t_{i,j − k}),k = 1,…,j − 1}, will still depend on Y(t_i,1),…,Y(t_{i,j − 2}).
By taking subsets of the form I_ij = {1,…,k_ij}, for some k_ij{1,…,j − 1}, we can modulate how far in the past the probability q_ij(θ) is conditioned. Depending upon the memory of the process, P(θ) will provide a substitute to L(θ) that will yield more or less efficient estimators. For “near” Markovian processes, the proposed approach may not sacrifice much efficiency when using k_ij = 1 or 2 (e.g.) while greatly simplifying simulation-based inference.

In practice, the choice of I_ij will be driven by practical considerations (such as reducing the number of mismatches). The maximum composite likelihood estimator possesses desired asymptotic properties stated below.

PROPOSITION 1

Under classical regularity conditions, the maximum composite likelihood estimator is consistent as n→∞. Furthermore, denoting the true value of θ by θ^*, we have that

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx8_ht.jpg

where

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx9_ht.jpg

and

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx10_ht.jpg

The proof of Proposition 1 rests on standard arguments used to establish asymptotic properties of M-estimators (for instance, see van der Vaart, 1998) and of multinomial likelihood estimators (Andersen, 1980). In particular, the consistency relies on the fact that the composite likelihood function satisfies the inequality _θ^*{P(θ)} ≤ _θ^*{P(θ^*)} for every θ, with the equality holding true for θ = θ^* only if a condition on the identifiability of the model parameters applies. The inequality is easily established by noticing that each individual terms logq_ij(θ) is a log-likelihood object and therefore satisfies the classical inequality _θ^*{logq_ij(θ)} ≤ _θ^*{logq_ij(θ^*)}. The matrices A(θ) and B(θ) will generally be different from each other, and, in most cases, the composite likelihood estimator will not be fully efficient. Lindsay (1988) commented further on the properties of composite likelihood estimators.

4.2. Practical considerations

The probabilities q_ij(θ)s will generally have no closed-form expressions and estimation of θ cannot proceed through direct maximization of P(θ). We describe below 2 algorithms designed to compute a composite likelihood estimate for θ. Both algorithms use bias-corrected simulation-based approximations to logq_ij(θ) that are computed in the following way.

Let _s(θ,U_s), s = 1,…,S, be S independent clones simulated under the assumed model for the parameter value θ, where U_s denotes a set of r.v.'s uniformly distributed on the interval [0,1]. These r.v.'s are used as seeds to generate the lifespans and the progeny vectors of all cells in the sth simulated clone. Write U^(S) = (U₁,…,U_S). Let Y_ij^s(θ,U_s) = {Y_ij1^s(θ,U_s),…,Y_ijD^s(θ,U_s)} be a vector of cell counts for the D morphologically distinguishable cell types for the sth simulated clone at time t_ij. To reduce the notational burden, we shall omit the notation U^(S) whenever unnecessary, and, for instance, write Y_ij^s(θ) in place of Y_ij^s(θ,U_s). Let N_ij(θ;I_ij) = N_ij(θ,U^(S);I_ij) denote the number of simulated clones whose composition matches that of the ith experimental clone at the time points {t_{i,j − k},∀k [set membership]

I_ij}. Formally, we have

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx11_ht.jpg

where 1_{·} is the indicator function. The conditional probability q_ij(θ) can be approximated by

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx12_ht.jpg

and a natural substitution in P(θ) produces the simulated composite likelihood function

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx13_ht.jpg

This simulated composite likelihood, however, suffers from the same limitations as those previously mentioned for the simulated log-likelihood function; that is, it is a biased estimator of P(θ), and it is not defined when N_ij(θ;I_ij) = 0. We mitigate both limitations simultaneously by considering a bias-corrected version of P(θ) that is built as follows.

Conditional on {N_ij(θ;I_ij) = s_ij} and on {Y_{i,j − k}^s(θ) = Y_{i,j − k},∀k [set membership]

I_ij,∀s = 1,…,s_ij}, we have that N_ij(θ;I_ij [union or logical sum]

{0})~Bin(s_ij,q_ij(θ)), where Bin(s,q) denotes the binomial distribution with parameters s and q. Following Pettigrew and others (1986), a bias-corrected estimator of logq_ij(θ) takes the form

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx14_ht.jpg

The bias of An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx15_ht.jpg

is of smaller order—compared to that of An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx16_ht.jpg

—when s_ijq_ij(θ)→∞. When s_ijq_ij(θ)→0, however, the bias of either approximation increases proportionally to s_ij^{− 2}q_ij(θ)^{− 2}. Therefore, even with the biased-corrected estimator of logq_ij(θ), it will remain important to reduce as much as possible the number of mismatches so the simulation-based approximation to the objective function will not be too biased. The performance of the estimator will otherwise be affected. The composite likelihood approach offers flexibility in choosing an estimating function that reduces the number of mismatches or even eliminates them. Note that An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx15_ht.jpg

remains defined when there are mismatches.

This motivated us to consider the simulated composite likelihood defined as

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx17_ht.jpg

There exist multiple ways to compute an estimator of θ based on P_bc^S(θ). Two methods are discussed below.

Method 1: sample path.

In order to estimate θ, a first approach consists in maximizing the simulated composite likelihood P_bc^S(θ,U^(S)) for a given vector U^(S). In doing so, the simulated composite likelihood function is not subject to random fluctuations during the optimization. The resultant estimator can be proven to converge to An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx7_ht.jpg

as S→∞ (see Hyrien, 2007, for a proof of a similar result). Due to the discrete nature of branching processes, the function P_bc^S(θ,U^(S)) is discontinuous, and optimizations algorithms that do not rest upon differentiability of the objective function (such as the simplex method of Nelder and Mead) can be used.

Method 2: stochastic approximation.

Alternatively, θ can be estimated using a multidimensional version of stochastic procedure of Kiefer–Wolfowitz (Kiefer and Wolfowitz, 1952; Blum, 1954). In this alternate schema, an estimate of θ is computed through the sequence of r.v.'s {θ^(k)}_{k ≥ 1} defined as θ^{(k + 1)} = θ^(k) + a_kΔ_k, where Δ_k = (Δ_k,1,…,Δ_k,q)^′, with

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx18_ht.jpg

where U_{k, −}^(S) and U_{k, +}^(S), k = 1,2,…, are i.i.d. copies of U^(S), where e_l, l = 1,…,q, are q×1 vectors with all elements equal to 0, except for a 1 in the lth place, and where (a_k)_{k ≥ 1} and (b_k)_{k ≥ 1} are 2 given sequences ofnonrandom positive numbers such that

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx19_ht.jpg

Usual choices for (a_k)_{k ≥ 1} and (b_k)_{k ≥ 1} include the following: a_k = a₀/k^α and b_k = b₀/k^β, for some given positive constants a₀, b₀, α, and β. In our analysis of clonal data on O-2A progenitor cells, we used a₀ = 1, b₀ = 1, α = 1, and β = 1/3. Other choices are possible. If P_bc^(S)(θ) was an unbiased estimator for P(θ), θ^(k) would converge to An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx7_ht.jpg

as k→∞. This is not formally the case, but because P_bc^S(θ) is first-order bias-corrected, θ^(k) is expected converge to a value close to An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx7_ht.jpg

The computation of composite likelihood estimates may take up to more than a day, depending on the number of simulations S, the number of model parameters, and the sample size. The number of mismatches increases substantially with the number of conditioning observations, and eliminating the mismatches may require increasing the number of simulations out of practical limits. The numbers of simulations being equal, the stochastic approximation procedure is expected to provide more accurate estimates than the sample path approach since the noise arising from the simulations is averaged out as the sequence {θ^(k)}_{k ≥ 1} progresses. This procedure is more time-consuming because it requires evaluation of the objective function at 2q points at every step k, while the sample path approach needs only one such evaluation at every step. In practice, it is probably best to begin with the sample path approach and complete estimation using stochastic approximation, with a₀ and b₀ appropriately chosen to avoid large jumps of the sequence {θ^(k)}_{k ≥ 1}.

5. SIMULATION STUDIES

We performed simulation studies to evaluate the performance of the proposed simulated composite likelihood estimator in samples of finite sizes. We also compared it with the simulated pseudo maximum likelihood estimator proposed in past publications (Hyrien and others, 2005a, Hyrien and others, 2005b; Hyrien, 2007). In order to keep the computing time reasonable (each estimates required between 1 h of computing time with S = 10000 simulated clones, and up to 6 h with S = 50000 simulated clones), the study was conducted using 2 simple branching processes that are special cases of our model of the generation of oligodendrocytes. The conclusions of these studies are therefore not meant to be exhaustive but they corroborate what can be intuitively expected. The 2 processes that we considered were defined as follows:

Model 1 (a binary splitting process): This model assumes that every cell divides into 2 new cells upon completion of its lifespan, which follows a gamma distribution with mean m and variance σ². When m = σ, the process reduces to the (Markov) linear birth process.
Model 2 (a 2-type process): This model considers 2 cell types: every type-1 cell divides into 2 new type-1 cells with probability p or it turns into a single type-0 cell with probability 1 − p. The lifespan of any type-1 cell follows a gamma distribution with mean m and variance σ². Type-0 cells have infinite lifespans.

Longitudinal clonal data were simulated from these 2 models, and cell clones were observed every day from day 1 to day 6 (same experimental design for both models). We considered various sample sizes: n = 20, 50, 100; and used S = 10000 and S = 20000 simulated clones to compute the simulated composite likelihoods for Models 1 and 2, respectively. For each simulated data set, we computed the simulated composite likelihood under various conditioning sets, all of the form I_ij = or I_ij = {1,…,j∧j₀}, for j₀ = 1,…,5, where j∧j₀ denotes the smallest value of j and j₀. The particular case I_ij = {1,…,j∧5} corresponds to a simulated maximum likelihood estimator. The pseudo maximum likelihood estimator (Hyrien and others, 2005a, Hyrien and others, 2005b; Hyrien, 2007), computed for comparison, is defined as the parameter vector that maximizes the pseudolikelihood function

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx20_ht.jpg

where m_i(θ) = _θ(Y_i) and V_i(θ) = Var_θ(Y_i), and where detV denotes the determinant of the matrix V. When these moments have no closed-form explicit expressions, they can be replaced by approximations computed either from simulations (Hyrien and others, 2005a, Hyrien and others, 2005b; Hyrien, 2007) or using the method of saddlepoint (Hyrien and others, 2010). The pseudo maximum likelihood estimator is consistent and asymptotically Gaussian (Hyrien, 2007). A summary of the simulation results is displayed in Tables 1 and and22.

Table 1.

Results of the simulation study based on a binary splitting process (see text for explanation). CL_k = simulated composite likelihood with I_ij = {1, 2…, j ∧ k}; ML = simulated maximum likelihood PL = simulated pseudolikelihood; s.e. = (more ...)

Table 2.

Results of the simulation study based on a 2-type process (see text for explanation). CL_k = simulated composite likelihood with I_ij = {1, 2…, J ∧ k}; ML = simulated maximum likelihood; PL = simulated pseudolikelihood; s.e. = standard error (more ...)

The biases of the estimators tended to decrease as the sample size increased. The standard errors depended on the method of estimation being used. Overall the pseudo maximum likelihood estimator appeared to be the least efficient of all tested estimators. This was particularly true for the smallest sample sizes. The simulated maximum likelihood estimator was not particularly more efficient than the simulated composite likelihood estimators. For instance, in the case of Model 2, the simulations indicated that the marginal composite likelihood estimator and the composite likelihood estimators defined by conditioning upon the most recent one or two observations had smaller standard errors. This was particularly true for the largest sample sizes. The most likely explanation of the fact that some estimators performed worse when we conditioned upon more observations is that the simulation-based approximations to the conditional probabilities that make up the expression for the estimating functions are more noisy in the latter case (due to the simulations) than when we condition on fewer observations.

We also found that the bias decreases generally quickly with increasing sample size but not always. For instance, the estimator of the probability of division (Model 2), say An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx21_ht.jpg

An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx21_ht.jpg

, remained substantially biased even with n = 100 observed clones. When conditioning upon the most recent observation only, the bias of An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx21_ht.jpg

was about 0.06 with a sample size of n = 20 clones and close to 0.05 with a sample size of n = 50 clones. The bias became larger as the composite likelihood estimator was conditioning upon an increasing number of past observations. We performed additional simulations to assess whether this bias could be attributed (in part at least) to a value of S that was too small. Thus, we considered S = 50000 simulated clones (instead of S = 20000). As a result of time constraints, we restricted our investigations to the above 2 cases and found that the bias of An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx21_ht.jpg

decreased to 0.03 both when n = 20 and when n = 50, and thus was reduced by almost a factor of 2. Also, the standard error of An external file that holds a picture, illustration, etc.
Object name is biostskxq050fx21_ht.jpg

decreased to 0.06 and 0.04, respectively. The estimators of m and σ improved also in this last run of simulations, with standard errors smaller than those reported in Table 2.

6. AN APPLICATION TO OLIGODENDROCYTES GENERATION EXPOSED TO NT-3

Oligodendrocytes produce the myelin sheaths that enwraps axons in the central nervous system. They are involved in signal propagation along the nerves and are particularly essential for normal impulse conduction. Some demyelinating diseases, such as multiple sclerosis, have been shown to be associated with injury or dysfunction of the oligodendrocyte population (Blakemore, 2008; Korn, 2008). The biology of oligodendrocytes and their lineage is well documented in the literature (e.g. Franklin and Ffrench-Constant, 2008). Oligodendrocytes are assumed to be nonproliferating, terminally differentiated cells, and are generated from their immediate precursor cells, the O-2A progenitor cells. These O-2A progenitor cells have the ability to self-renew by dividing, to differentiate into oligodendrocytes, or to die. Due to the critical role of oligodendrocytes, the proper regulation of the differentiation of O-2A progenitor cells into oligodendrocytes is therefore particularly important for the normal functioning of the central nervous system. O-2A progenitor cells have been shown to respond to the presence of various environmental signals or factors that may modulate their ability to self-renew and differentiate, and one of the primary objectives of the research performed on the oligodendrocyte lineage is to identify factors that could restore the normal balance in the generation of oligodendrocytes.

NT-3 is one member of a family of growth factors known as neurotrophins. Althaus and others (2008) gave an overview on neurotrophins and their effects on oligodendroglial cells. It has been demonstrated that, in rodents and humans at least, NT-3 may be synthesized by oligodendrocytes. Although NT-3 has been shown to be involved in the growth, differentiation, and survival of neurons (Zhou and Rush, 1996; Kalb, 2005), its role (and that of other neurotrophins as well) on the regulation of nonneuronal cells of the central nervous system is still unclear.

In an attempt to better understand the effect of NT-3 on the regulation of the processes of division and differentiation of O-2A progenitor cells, we conducted the longitudinal clonal experiments presented in Section 2.2. In particular, the objectives of this experiment were to assess the hypothesis that the presence of NT-3 alone (i.e. without coexposure to PDGF) in pure cultures has no effect on the proliferation and differentiation of O-2A progenitor cells, as suggested by Barres and others, 1994, and to gain some quantitative insight into how NT-3 may alter division and differentiation of O-2A progenitor cells.

Histograms for the number of O-2A progenitor cells and for the number of oligodendrocytes per clone with and without NT-3 suggested a striking difference between these 2 experimental conditions in respect to the number of O-2A progenitor cells. When cells were cultured without NT-3, the number of O-2A progenitor cells appeared to increase over time, suggesting that these cells maintained some capacity for self-renewal after 6 days of culture even in the absence of PDGF. By day 6, more than 60% of the clones still contained at least one O-2A progenitor cell, and the clones contained on average 1.9 O-2A progenitor cells and 2.1 oligodendrocytes. In contrast, in clones grown in the presence of NT-3, O-2A progenitor cells underwent division first, just as when they were cultured without NT-3, and then appeared to differentiate massively. By day 6, less than 20% of the clones still contained at least one O-2A progenitor cell. The mean numbers of oligodendrocytes per clone (~eq1.9) was not much different than that without NT-3, but the number of O-2A progenitor cells was considerably lower (~eq0.8).

We investigate further the effects of NT-3 on the processes of division and differentiation of O-2A progenitor cells to determine whether NT-3 affected the decision of O-2A progenitor cells to divide or differentiate the timing of cellular division or the timing of cellular differentiation, or all of the above. Changes in cell fate would be seen through either an increase or a decrease of the probability of division of O-2A progenitor cells, whereas a modification of the time it takes for O-2A progenitor cells to divide or to differentiate would be reflected by changes in the parameters of the distribution of the time to division and/or of the distribution of the time to differentiation of these cells.

We fitted our branching process model using the proposed simulated composite likelihood approach. The ultimate parameter estimates were obtained using 15 000 simulated clones to approximate the distribution of the number of cells of each type per clone under the above model. We started from multiple initial values and let optimization run for more than a day each time. We used the conditioning sets I_ij = ,{1},{1,j∧2},…,{1,2,…,j∧5}). In the last case, we obtain a simulated maximum likelihood estimator. The model was fitted separately to the data obtained with and without NT-3 because none of the model parameters were common to both experimental conditions.

In our experiments, the follow-up period was long enough so all initiator precursor cells had either divided or differentiated by the time the experiment ended. By determining the composition of a clone, it is possible to decipher whether the corresponding initiator cell divided or differentiated (in the latter case, the clone always contains a single oligodendrocyte). Thus, let N_div and N_diff = n − N_div denote the numbers of initiator cells having divided and differentiated by the end of the experiment, respectively. For every 0 ≤ j₁ < j₂, let N_div(j₁,j₂) denote the number of initiator cells having divided between day j₁ and j₂. Define the vector N_div = {N_div(t − 1,t);t = 1,…,6}. Conditional on N_div = N, N_div follows a multinomial distribution with parameters N and {p₁(m₁,σ₁),…,p₆(m₁,σ₁)}, where p_t(m,σ) = G₁(t;m,σ) − G₁(t − 1;m,σ) denotes the probability that any initiator cell that will ultimately divide, completes its cycle between day t − 1 and t. The estimates for m₁ and σ₁ can therefore be computed by maximizing the multinomial likelihood function associated with the vector of counts N_div. The mean and variance of the time to differentiation of O-2A progenitor cells can be estimated in a similar fashion for clones starting with an initiator cell that ended up differentiating. We used these maximum likelihood estimates as starting values in our composite likelihood analyses.

Table 3 reports parameter estimates resulting from these analyses. Also reported are the total numbers of mismatches for the final parameter estimates. In our analyses of cell clones cultured without NT-3, the number of mismatches remained low (0-1) until we started to condition upon the most recent 3 (or more) observations. For cell clones exposed to NT-3, the number of mismatches increased substantially when we started to condition upon the most recent observations. In either case, the largest number of mismatches we encountered ranged between 6 and 9, which is relatively high compared to the number of cell clones observed for each experimental condition (n = 20). Alternatively, we could have increased the number of simulations to reduce the number of mismatches, but each analysis required between 1 and 2 days. We therefore discuss results obtained with the marginal composite likelihood (no conditioning) for cells exposed to NT-3, and with the composite likelihood conditioning upon the most recent 2 observations for control (unexposed) cells. Standard errors for the corresponding parameter estimates were obtained using a parametric bootstrap approach. We compared parameter estimates across treatment groups using a Wald test. The fitted marginal distribution for the number of O-2A progenitor cells and for the number of oligodendrocytes per clone is displayed in Figure 2 as solid lines. The estimated probabilities of division of O-2A progenitor cells are presented in Figure 3 (panels A and B) as a function of the cell generation for each culture condition. These analyses pointed out the following conclusions:

(a)In each experiment, the probability of division of O-2A progenitor cells was found to decrease when the number of divisions increased. This is consistent with the results of previously analyzed data sets (Yakovlev, Boucher, and others, 1998, Yakovlev, Mayer-Pröschel, and others, 1998, Yakovlev and others, 2000; von Collani and others, 1999; Boucher and others, 1999, Boucher and others, 2001; Zorin and others, 2000; Hyrien and others, 2005a, Hyrien and others, 2005b, Hyrien and others, 2006).
(b)Exposure to NT-3 decreased the probability of division of O-2A progenitor cells. For instance, cells of generation 3 divided with a probability estimated as 0.35 (vs. 0.53) when cells were exposed (vs. unexposed) to NT-3. This may explain that the numbers of O-2A progenitor cells were much smaller in the presence of NT-3. The difference between the probabilities of division across treatment groups was almost significant (p = 0.07).
(c)In both culture conditions, the mean time to division was shorter for initiator O-2A progenitor cells than for O-2A progenitor cells of subsequent generations (40 vs. 50 h with NT3, and 31.5 vs. 44 h without NT3). This is consistent with results from previous studies based on time-lapse experiments (Hyrien and others, 2006).
(d)Exposure to NT-3 increased the mean time to division of O-2A progenitor cells (e.g. 50 vs. 44 h for cells of generation greater than or equal to 2) and decreased their mean time to differentiation (36 vs. 44 h). When comparing the distribution of the time to division across treatment groups using the Wald test, we found that the difference was significant (p < 0.0001). The effect of NT-3 on the time to differentiation was almost significant (p = 0.05).

Table 3.

Estimates of the model parameters (means and standard deviations expressed in hours). See text for explanations. CL_k = simulated composite likelihood with I_ij = {1, 2…, j ∧ k}; ML = simulated maximum likelihood; s.e. = standard errors

Fig. 3.

Panels A–B: Probability of division of O-2A progenitor cells versus the number of divisions estimated from the fitted model in the absence and in the presence of NT-3. Panels C–F: Mean number of O-2A progenitor cells and mean number of (more ...)

The overall conclusion of our analysis is that NT-3 may modulate the proliferation of O-2A progenitor cells and their differentiation into oligodendrocytes, and that it does not act only in concert with PDGF. In particular, we found that NT-3 induces differentiation of O-2A progenitor cells into oligodendrocytes by decreasing the probability of division, by increasing the mean time to division, and by decreasing the mean time to differentiation. The results of our analysis are tentative and should be confirmed in other experiments.

The effect of NT-3 may be perceived as subtle when expressed as a change of the probability of division, of the distribution of the time to division, and of the distribution of the time to differentiation. To better appreciate the impact of such changes on the composition of cell clones over time, we simulated the fitted branching process over 20 days and calculated the mean number of O-2A progenitor cells and the mean number of oligodendrocytes cultured with and without NT-3 as a function of time (Figure 3). While the number of O-2A progenitor cells declines under either culture conditions (as expected from the theory of branching processes since lim_k→∞p_k < 0.5), the number of oligodendrocytes differed substantially. For instance, by day 20, the number of oligodendrocytes was 2-fold larger in untreated cells. This seems to contradict our conclusion that NT-3 induces differentiation, but it does not actually. Indeed the best strategy to increase the number of oligodendrocytes over a medium-term period is to first increase the pool of progenitor cells by self-renewal, and then have some of these progenitor cells differentiate into oligodendrocytes. This seems to be what happened when cells were not exposed to NT-3 compared to when they were.

The mean numbers of oligodendrocytes in clones exposed and unexposed to NT-3 drifted apart from each other after day 6 only. This would explain why our permutation test did not detect any significant difference between day 0 and day 6 in counts of oligodendrocytes.

7. DISCUSSION

The proposed simulated composite likelihood method provides a good alternative to existing methods for the analysis of clonal data using age-dependent branching processes. It possesses good properties both asymptotically and in samples of finite sizes. It appeared to outperform the pseudo maximum likelihood estimator proposed in past publications (Hyrien and others, 2005a, Hyrien and others, 2005b; Hyrien, 2007), and it mitigates the computational issues encountered with the simulated likelihood approach (Zorin and others, 2000). Another advantage of the method is that it is adaptable to the complexity of the situation. It allows the selection of an estimating function that reduces the number of mismatches. The proposed method was framed within the context of a particular multitype age-dependent branching process of the generation of oligodendrocytes in vitro. The structure of this model could be modified, if needed, using the more general process described by Hyrien and others (2005a).

FUNDING

National Institutes of Health (1R01 CA134839-01, 2R01 NS039511, 1R01AI069351, R01 HD39702, 2 P30 ES001247, T32 ES007271).

Acknowledgments

The authors acknowledge the reviewers and an associate editor for their careful reading of the manuscript and for their suggestions. Conflict of Interest: None declared.

References

Ahlgren SC, Wallace H, Bishop J, Neophytou C, Raffa MC. Effects of thyroid hormone on embryonic oligodendrocyte precursor cell development in vivo and. in vitro. Molecular and Cellular Neuroscience. 1997;9:420–432. [PubMed]
Althaus HH, Klöppner S, Klopfleisch S, Schmitz M. Oligodendroglial cells and neurotrophins: a polyphonic cantata in major and minor. Journal of Molecular Neuroscience. 2008;35:65–79. [PubMed]
Andersen EB. Discrete Statistical Models with Social Science Applications. New York: North-Holland; 1980.
Azzalini A. Maximum likelihood estimation of order m for stationary stochastic processes. Biometrika. 1983;70:381–387.
Barres BA, Raff MC, Gaese F, Bartke I, Dechant G, Barde YA. A crucial role for neurotrophin-3 in oligodendrocyte development. Nature. 1994;367:371–375. [PubMed]
Besag J. Spatial interaction and statistical analysis of lattice systems. Journal of the Royal Statistical Society, Series B. 1974;25:192–236.
Blakemore WF. Regeneration and repair in multiple sclerosis: the view of experimental pathology. Journal of the Neurological Sciences. 2008;265:1–4. [PubMed]
Blum JR. Multidimensional stochastic approximation methods. The Annals of Mathematical Statistics. 1954;23:737–744.
Boucher K, Yakovlev AY, Mayer-Pröschel M, Noble M. A stochastic model of temporarily regulated generation of oligodendrocytes in vitro. Mathematical Biosciences. 1999;159:47–78. [PubMed]
Boucher K, Zorin A, Yakovlev AY, Mayer-Pröschel M, Noble M. An alternative stochastic model of generation of oligodendrocytes in cell culture. Journal of Mathematical Biology. 2001;43:22–36. [PubMed]
Chandler RE, Bate S. Inference for clustered data using the independence loglikelihood. Biometrika. 2007;94:167–183.
Cox DR, Reid N. Parameter orthogonality and approximate conditional inference. Journal of the Royal Statistical Society, Series B. 1987;49:1–39.
Dietrich J, Han R, Yang Y, Mayer-Pröschel M, Noble M. CNS progenitor cells and oligodendrocytes are targets of chemotherapeutic agents in vitro and in vivo. Journal of Biology. 2006;5:22. [PMC free article] [PubMed]
Franklin RJ, Ffrench-Constant C. Remyelination in the CNS: from biology to therapy. Nature Reviews Neuroscience. 2008;9:839–855. [PubMed]
Heagerty PJ, Lele SR. A composite likelihood approach to binary spatial data. Journal of the American Statistical Association. 1998;93:1099–1111.
Hyrien O. A pseudo maximum likelihood estimator for discretely observed multi-type Bellman-Harris branching processes. Journal of the Statistical Planning and Inference. 2007;137:1375–1388.
Hyrien O, Ambeskovic I, Mayer-Pröschel M, Noble M, Yakovlev A. Stochastic modeling of oligodendrocyte generation in cell culture: model validation with time-lapse data. Theoretical Biology and Medical Modelling. 2006;3:21. [PMC free article] [PubMed]
Hyrien O, Chen R, Mayer-Pröschel M, Noble M. Saddlepoint approximations to the moments of multitype age-dependent branching processes, with applications. Biometrics. 2010;66:567–577. [PMC free article] [PubMed]
Hyrien O, Mayer-Pröschel M, Noble M, Yakovlev A. A stochastic model to analyze clonal data on multi-type cell populations. Biometrics. 2005a;61:199–207. [PubMed]
Hyrien O, Mayer-Pröschel M, Noble M, Yakovlev A. Estimating the life-span of oligodendrocytes from clonal data on their development in cell culture. Mathematical Biosciences. 2005b;193:255–274. [PubMed]
Hyrien O, Zand MS. A mixture model with dependent observations for the analysis of CFSE-labeling experiments. Journal of the American Statistical Association. 2008;103:222–239.
Ibarrola N, Mayer-Pröschel M, Rodriguez-Pena A, Noble M. Evidence for the existence of at least two timing mechanisms that contribute to oligodendrocyte generation in vitro. Developmental Biology. 1996;180:1–21. [PubMed]
Kalb R. The protean actions of neurotrophins and their receptors on the life and death of neurons. Trends in Neurosciences. 2005;28:5–11. [PubMed]
Kiefer J, Wolfowitz J. Stochastic estimation of the maximum of a regression function. The Annals of Mathematical Statistics. 1952;23:462–466.
Korn T. Pathophysiology of multiple sclerosis. Journal of Neurology. 2008;255(Suppl 6):2–6. [PubMed]
Lindsay BG. Composite likelihood methods. Contemporary Mathematics. 1988;24:221–239.
Noble M, Pröschel C, Mayer-Pröschel M. Getting a GR(i)P on oligodendrocyte development. Developmental Biology. 2004;265:33–52. [PubMed]
Pettigrew HM, Gart JJ, Thomas DG. The bias and higher cumulants of the logarithm of a binomial variate. Biometrika. 1986;73:425–435.
Rodríguez-Peña A. Oligodendrocyte development and thyroid hormone. Journal of Neurobiology. 1999;40:497–512. [PubMed]
Smith J, Ladi E, Mayer-Pröschel M, Noble M. Redox state is a central modulator of the balance between self-renewal and differentiation in a dividing glial precursor cell. Proceedings of the National Academy of Sciences of the United States of America. 2000;18:10032–10037. [PubMed]
van der Vaart AW. Asymptotic Statistics. Cambridge: Cambridge University Press; 1998.
Varin C, Vidoni P. A note on composite likelihood inference and model selection. Biometrika. 2005;92:519–528.
von Collani E, Tsodikov A, Yakovlev A, Mayer-Pröschel M, Noble M. A random walk model of oligodendrocyte generation in vitro and associated estimation problems. Mathematical Biosciences. 1999;159:189. [PubMed]
Yakovlev AY, Boucher K, Mayer-Pröschel M, Noble M. Quantitative insight into proliferation and differentiation of oligodendrocyte-type 2 astrocyte progenitor cells in vitro. Proceedings of the National Academy of Sciences of the United States of America. 1998;95:14164–14167. [PubMed]
Yakovlev AY, Mayer-Pröschel M, Noble M. A stochastic model of brain cell differentiation in tissue culture. Journal of Mathematical Biology. 1998;37:49–60. [PubMed]
Yakovlev AY, Mayer-Pröschel M, Noble M. Stochastic formulations of a clock model for temporally regulated generation of oligodendrocytes in vitro. Mathematical and Computer Biology. 2000;32:125–137.
Zhou XF, Rush RA. Functional roles of neurotrophin 3 in the developing and mature sympathetic nervous system. Molecular Neurobiology. 1996;13:185–197. [PubMed]
Zorin AA, Yakovlev AY, Mayer-Pröschel M, Noble M. Estimation problems associated with stochastic modeling of proliferation and differentiation of O-2A progenitor cells in vitro. Mathematical Biosciences. 2000;167:109–121. [PubMed]

Articles from Biostatistics (Oxford, England) are provided here courtesy of
Oxford University Press