A case study on the use of scale separation-based analytic propagators for parameter inference in stochastic gene regulation

About: The article was published on 2015-01-01 and is currently open access. It has received 9 citations till now. The article focuses on the topics: Inference.

Summary

1. INTRODUCTION

• Gene expression is a complex and highly regulated multi-step process that is responsible for the timely synthesis of proteins necessary for cellular function.
• A simple, “two-stage” model for stochastic gene expression consists of a constitutively active gene from which an mRNA molecule can be transcribed, and protein, the production of which depends on the instantaneous abundance of mRNA (see Fig. 1A).
• In recent work [14], the procedure developed in [2] was extended to capture departure from the assumption of perfect scale separation: the ratio of degradation rates of protein and mRNA, denoted ε, was taken to be small and positive instead of zero, as was the case in [2].
• While fluorescence microscopy yields only time-series for the intensity, these can nonetheless be converted into absolute protein numbers if a calibration factor of molecules per unit intensity can be estimated, see e.g. [19].
• For comparison, both propagators are also contrasted with an approximate solution of the CME that is computed using a finite state projection.

2.1. Two-stage Gene Expression Model

• The authors model gene expression as a two-stage process, whereby DNA is transcribed to mRNA, which is then translated into protein (see Fig. 1A).
• (1) Here,m and n denote mRNA and protein copy numbers, respectively, a is the non-dimensional transcription rate and b is the non-dimensional translation rate, while the degradation rates of mRNA and protein are given by γ and 1, respectively (cf. Fig. 1A).
• Finally, τ denotes a suitably nondimensionalized time variable.
• It follows that for ε sufficiently small, the dynamics of Eq. (1) will vary on two distinct time-scales: the long-term behavior of the system is naturally described on the “slow” τ -scale, while the “fast” transients evolve according to the rescaled time t := τε .

2.2.2. Uniform (First-Order) Propagator

• Here, ε denotes the per- turbation parameter, as before, while t is the fast time variable.
• The authors remark that the transition probability Pn|n0(τ, ε) contributes on the slow τ -scale in Eq. (3), while the t-dependent contribution in Eq. (3) accounts for the transient dynamics on the fast time-scale.

2.3. Special Cases of the Hypergeometric Functions

• Care must be taken when evaluating the hypergeometric function 2F1(a, b, c, z).
• The following special cases are of use [20].

2.4. Stochastic Simulation

• Stochastic simulations were performed using the StochKit 2.0 [21] simulation framework and the standard stochastic simulation algorithm [3], with a non-dimensionalized transcription rate a = 20 and a non-dimensionalized translation rate b = 2.5, corresponding to “regime I” as defined in [14].
• Each value of γ was simulated 20 times, and the resulting trajectories were used for computing the probability landscapes of the rescaled model parameters a and b.
• Protein quantities were observed without measurement noise at intervals of 0.1 time units.
• All simulation runs assumed zero molecules of mRNA and protein initially, i.e., m0 = 0 = n0.

2.5. Implementation

• Both the zeroth-order propagator Pn|n0 , Eq. (2), and the uniform propagator Pn|n0 , Eq. (3), were implemented in C++ with a Matlab mex-file interface.
• Special functions were evaluated using the GNU scientific library [22], the Hyp_2F1 function implementation of the Gauss hypergeometric function [23], and the Algorithm 910multiprecision special function library [24].
• While the difference of such numbers is potentially below a double precision machine error of approximately 10−13, they are nonetheless essential in the correct computation of the transition probabilities.
• Their C++ implementation is still inaccurate in some extreme cases, typically for very large protein numbers n, due to numerical differences which are sometimes as small as 10−370 in Eq. (7), but which unfortunately cannot be neglected as they are inflated by the remaining terms in the expression.
• Such inaccuracies are infrequent, though, and generally occur during transitions for which the uniform propagator yields nonphysical values; thus, they do not substantially affect their analysis, or the conclusions obtained in this study.

3. RESULTS AND DISCUSSION

• To assess the applicability of the zeroth-order propagator Pn|n0(τ, 0), Eq. (2), and the uniform propagator Pn|n0(τ, t, ε), Eq. (3), for parameter inference in the two-stage gene expression model, the authors simulated time-series with a specific parameter pair (a∗, b∗).
• Then, the authors computed the likelihood of the observed data set on the basis of the two propagators for a range of values for the parameters a and b.
• For simplicity, the authors assumed the scale separation γ between mRNA and protein lifetimes to be known (see Methods for definitions).

3.2. Parameter Inference

• The term ni represents the number of proteins at measurement time ti.
• Thus, the authors compute the logarithm of the probability of each transition, from ni−1 protein molecules at time ti−1 to ni molecules at time ti, in the sequence of observed measurements (see Fig. 1B, inset).
• In order to estimate the parameters a and b from simulated protein time-courses, Eq. (13) has to be evaluated very frequently.
• The authors thus developed a numerically stable expression for the uniform propagator Pn|n0 (see Section 2 2.2), and they used an efficient implementation in C++ for both propagators that results in reasonable runtimes; see Section 2 2.5 for details.

3.3. Comparison of Propagator Accuracy and Efficiency

• For each pair (a, b), the authors computed the log-likelihood L(a, b), thus obtaining a likelihood landscape that should ideally have its maximum, the maximum likelihood estimator (MLE), at the true parameter values (a∗, b∗).
• The authors immediately encountered the obstacle that the uniform propagator Pn|n0 yields negative transition probabilities, or even probabilities larger than one, for some choices of (a, b).
• The resulting non-physicality is discussed in [14], and is due to the fact that Pn|n0 is derived from an asymptotic approximation; nonetheless, it is problematic when computing the overall log-likelihood, as the definition in Eq. (13) becomes meaningless.
• In Fig. 3B, a typical protein time-course with γ = 10 is shown (top), along with the log-likelihood obtained from the uniform propagator Pn|n0 , Eq. (3), for the true parameter values , and for the MLE (cyan).
• From that plot, it is obvious that large regions of the transition space yield non-physical values, shown in gray.

4. CONCLUSION

• The authors have investigated the utility of a propagator-based approach for approximating the transition probabilities in a simple two-stage gene expression model by attempting parameter inference from protein time-series.
• The simulations were initialized with zero molecules of both mRNA and protein; this simplification, as compared to a typical biological setting, does not affect the subsequent analysis.
• Alternatively, one could use Markov Chain Monte Carlo (MCMC) techniques to sample directly from the posterior in order to obtain the log-likelihood landscape [26].
• The variance of the noise then constitutes an additional unknown parameter σ which would have to be inferred.

Figures

FIG. 2: Averaged simulated log-likelihood landscapes for (a, b) = (20, 2.5) and γ = 100. Landscapes for single time-courses (left) are shown with contour lines drawn at intervals of 1 unit; contours for landscapes obtained from 20 time-courses (right) are drawn at intervals of 10 units. The averaged log-likelihood landscapes generated using a single time-course for the zeroth-order propagator, Eq. (2), (A), the uniform propagator, Eq. (3), (B), and a finite state projection approximation (C) for a single time-course each display a bias of the MLE with respect to the true model parameters (a∗, b∗). Notably, the landscape of the uniform propagator (B) shows distortions arising from non-physical transition probabilities for some parameter pairs (a, b). As the number of trajectories is increased to 20, the MLE converges to (a∗, b∗) for each of the zeroth-order propagator (D), uniform propagator (E), and the finite state projection (F). The averaged log-likelihood resulting from the finite state projection seems to be most tightly-peaked around the true parameter values (a∗, b∗).

FIG. 3: A. The averaged log-likelihood landscape of (a, b) for γ = 10 of the uniform propagator Pn|n0 shows prominent distortions in the contours caused by frequent non-computable transitions. The MLE (cyan) exhibits an obvious bias with respect to the true model parameters (a∗, b∗) (black). B. Inspection of a single time-course (shown on top), evaluated at (a∗, b∗) and at the MLE (bottom), reveals more non-computable transitions (indicated with white boxes below) for the MLE than for the true parameters; however, for those points that can be computed, the MLE probability is higher than for the true parameters, leading to a higher averaged probability and thus to a biased estimate of the parameters (a, b). C. Transition probability in regime I, with (a, b) = (20, 2.5). The transition marked with arrows in (B), from 55 to 57 molecules, results in a negative transition probability. D. The transition matrix for the uniform propagator in regime I from n(t) to n(t+∆t) proteins reveals a large region of non-computable transitions, shown in gray. E. The computability score C(a, b) shows that the MLE is biased towards the region with the lowest computability, for which most transitions are omitted from the averaged log-likelihood score L̄(a, b). F. By contrast, the transition matrix is fully computable in regime II, with (a, b) = (0.5, 100), corresponding to the region of bursty protein synthesis, i.e., to translational bursting.

FIG. 1: A. The two-stage model for gene expression, which captures stochastic birth and death of mRNA and protein, with non-dimensionalized parameters a for transcription, b for translation, and γ for mRNA degradation. B. Time-courses were simulated using the stochastic simulation algorithm, shown here for a = 20, b = 2.5, and γ = 10 or γ = 100. Probabilities can be computed for each protein transition from the analytical two-stage propagators given in Eqs. (2) and (9) (inset, probability distributions shown in gray). C. Analytical propagators can be used to compute the probability of observing a particular number of protein molecules at arbitrary points in time, conditional on the initial conditions. The prediction from the uniform propagator, Eq. (3) (blue background), shows good qualitative agreement with stochastic simulation (gray lines), as illustrated for γ = 10 here.

Full text

Edinburgh Research Explorer
A case study on the use of scale separation-based analytic
propagators for parameter inference in stochastic gene
regulation
Citation for published version:
Feigelman, J, Popovic, N & Marr, C 2015, 'A case study on the use of scale separation-based analytic
propagators for parameter inference in stochastic gene regulation', Journal of Coupled Systems and
Multiscale Dynamics, vol. 3, no. 2, pp. 164-173. https://doi.org/10.1166/jcsmd.2015.1074
Digital Object Identifier (DOI):
10.1166/jcsmd.2015.1074
Link:
Link to publication record in Edinburgh Research Explorer
Document Version:
Peer reviewed version
Published In:
Journal of Coupled Systems and Multiscale Dynamics
General rights
Copyright for the publications made accessible via the Edinburgh Research Explorer is retained by the author(s)
and / or other copyright owners and it is a condition of accessing these publications that users recognise and
abide by the legal requirements associated with these rights.
Take down policy
The University of Edinburgh has made every reasonable effort to ensure that Edinburgh Research Explorer
content complies with UK legislation. If you believe that the public display of this file breaches copyright please
contact openaccess@ed.ac.uk providing details, and we will remove access to the work immediately and
investigate your claim.
Download date: 09. Aug. 2022

Journal of Coupled Systems and Multiscale Dynamics
A case study on the use of scale separation-based analytical propagators for
parameter inference in models of stochastic gene regulation
Justin Feigelman
1,2
, Nikola Popovi
´
c
3
, and Carsten Marr
1
1
Helmholtz Zentrum München - German Research Center for Environmental Health,
Institute of Computational Biology, Ingolstädter Landstraße 1, 85764 Neuherberg Germany
2
Technische Universität München, Center for Mathematics,
Chair of Mathematical Modeling of Biological Systems,
Boltzmannstraße 3, 85748 Garching, Germany and
3
University of Edinburgh, School of Mathematics and Maxwell Institute for Mathematical Sciences,
James Clerk Maxwell Building, King’s Buildings,
Peter Guthrie Tait Road, Edinburgh EH9 3FD, United Kingdom
(Dated: October 2, 2015)
Advances in long-term fluorescent time-lapse microscopy have made it possible to study the expression of
individual genes in single cells. In a typical setting, the intensity of one or more fluorescently-labeled proteins
is measured at regular time intervals. Such time-courses are inherently noisy due to both measurement noise
and intrinsic stochasticity of the underlying gene expression regulation. Fitting stochastic models to time-series
data remains a difficult task, partly because analytical and tractable expressions for the transition probabilities
cannot usually be derived in closed form.
In the present work, we employ a recently developed approach that is based on geometric singular perturba-
tion theory, as applied to the chemical master equation of a simple two-stage gene expression model, to compute
parameter likelihoods using synthetic protein time-series. We study the identifiability of model parameters in
this simple setting, and compare the performance of the perturbative (uniform) propagator to a previously pub-
lished, idealized (zeroth-order) propagator that assumes perfect time-scale separation between degradation of
mRNA and protein. We find that both propagators are useful for parameter inference when the scale separation is
sufficiently large. However, with decreasing separation, the uniform propagator sometimes yields non-physical
negative transition probabilities which render parameter inference difficult. Finally, we discuss the utility of
both propagators, and possible extensions thereof, for inference. For computational efficiency, the propagators
were implemented in C++ and embedded in Matlab; the code is available upon request.
Keywords: Multi-scale gene expression dynamics. Propagator approximation. Parameter inference. Geometric
singular perturbation.
1. INTRODUCTION
Gene expression is a complex and highly regu-
lated multi-step process that is responsible for the
timely synthesis of proteins necessary for cellular
function. At the molecular level, gene expression is
inherently stochastic due to random binding events
of transcription factors and the transcriptional ma-
chinery, which ultimately leads to mRNA transcrip-
tion with probabilities depending on the concentra-
tion of the reaction educts. Protein synthesis re-
quires a chance encounter of mRNA with ribosomes,
and mRNA or protein degradation an encounter with
the degradation machinery. Thus, models for gene
expression have to capture the stochasticity at both
mRNA and protein levels.
A simple, “two-stage” model for stochastic gene
expression consists of a constitutively active gene
from which an mRNA molecule can be transcribed,
and protein, the production of which depends on the
instantaneous abundance of mRNA (see Fig. 1A).
Both mRNA and protein are subjected to stochas-
tic degradation. Such a qualitative model can
be described mathematically as a two-dimensional
Markov jump process in the copy numbers of mRNA
and protein, with reaction probabilities that are func-
tions of the current state only (hence the Markov
property), and suitably chosen kinetic constants [1,
2].
While the two-stage model is easily simulated us-
ing stochastic simulation algorithms such as Gille-
spie’s algorithm [3], it is nonetheless a difficult task
to derive analytical expressions for the evolution of
mRNA and protein copy numbers with time. The
Markov process itself obeys the chemical master
equation (CME), an infinite-dimensional system of
ordinary differential equations, for which no exact
(closed-form) solutions are known in general. Nu-
merous approaches exist for the approximate solu-
tion of the CME, such as the linear noise approxi-
mation [4], a second-order Taylor series expansion
in the system size of the reaction volume; moment
equations and variants thereof [5, 6], which capture
an arbitrary number of statistical moments of the
stochastic process; finite state projection [7], a trun-
cation of the state-space of possible copy number
combinations, and many others (for an overview, see
1

Journal of Coupled Systems and Multiscale Dynamics
[8]). We further note that this particular model has
been studied using a variety of analytical and com-
putational techniques, see e.g. [9–12] or [13] for a
review of related modelling approaches.
An alternative analytical approach was developed
by Shahrezaei and Swain [2], wherein it is assumed
that mRNA molecules decay much faster than pro-
tein, a realistic assumption in many prokaryotic
cells. In the limit of a perfect scale separation in
which the decay of mRNA is instantaneous, the
CME underlying the two-stage model can be solved
analytically by the introduction of a generating func-
tion. The latter then obeys a first-order linear par-
tial differential equation, the solution of which can
be obtained via the method of characteristics. The
resulting analytical expression for the general time-
dependent joint probability density of mRNA and
protein, called the propagator of the system, is of
great utility for understanding its dynamics in time.
However, it is not valid when the assumption of
scale separation is violated, as is commonly the case
for eukaryotic cells. In recent work [14], the proce-
dure developed in [2] was extended to capture de-
parture from the assumption of perfect scale sepa-
ration: the ratio of degradation rates of protein and
mRNA, denoted ε, was taken to be small and pos-
itive instead of zero, as was the case in [2]. The
presence of the (singular) perturbation parameter ε
allows for the application of asymptotic techniques,
such as geometric singular perturbation theory [15]
and matched asymptotic expansions [16].
In the present case study, we explore the util-
ity of this newly developed perturbative approach
for propagator-based parameter inference in systems
with varying degrees of scale separation. Specif-
ically, our goal is to estimate molecular parame-
ters in the model from observations of protein abun-
dance only. Trajectories are simulated via Gille-
spie’s stochastic simulation algorithm in a param-
eter regime in which mRNA and protein are pro-
duced continuously, i.e., not in translational bursts.
The protein time-courses are sampled at regular time
intervals, thus mimicking a typical time-lapse flu-
orescence microscopy setup [17, 18]. While fluo-
rescence microscopy yields only time-series for the
intensity, these can nonetheless be converted into
absolute protein numbers if a calibration factor of
molecules per unit intensity can be estimated, see
e.g. [19]. We note that mRNA time-courses are not
observed, and that they are hence not used for pa-
rameter inference.
The zeroth-order propagator obtained by setting
ε = 0 [2] is then compared to a first-order propaga-
tor (in ε > 0) that is uniformly valid both on short
and on long time-scales [14], in terms of the abil-
ity of each to capture the correct parameters – i.e.,
the kinetic constants in the underlying chemical re-
action network – in the two-stage model for gene
expression. For comparison, both propagators are
also contrasted with an approximate solution of the
CME that is computed using a finite state projection.
A number of simplifying assumptions are made; no-
tably, we ignore impeding factors such as measure-
ment noise, uncertainty in the conversion from flu-
orescence intensity to protein numbers or low sam-
pling frequency of fluorescent signal. Rather, our
focus in this case study is on assessing the general
efficiency and accuracy of the propagator-based ap-
proach for parameter inference.
2. METHODS
2.1. Two-stage Gene Expression Model
We model gene expression as a two-stage process,
whereby DNA is transcribed to mRNA, which is
then translated into protein (see Fig. 1A). Denoting
the probability of observing m molecules of mRNA
and n molecules of protein in the system at time τ
by P
m,n
(τ), we find that the latter evolves according
to the non-dimensionalized CME [2, 4]
∂P
m,n
∂τ
= a(P
m−1,n
− P
m,n
)
+ γbm(P
m,n−1
− P
m,n
)
+ γ[(m + 1)P
m+1,n
− mP
m,n
]
+ [(n + 1)P
m,n+1
− nP
m,n
]. (1)
Here, m and n denote mRNA and protein copy num-
bers, respectively, a is the non-dimensional tran-
scription rate and b is the non-dimensional trans-
lation rate, while the degradation rates of mRNA
and protein are given by γ and 1, respectively
(cf. Fig. 1A). Finally, τ denotes a suitably non-
dimensionalized time variable.
As in [2, 14], we define the perturbation parame-
ter ε = γ
−1
here. It follows that for ε sufficiently
small, the dynamics of Eq. (1) will vary on two
distinct time-scales: the long-term behavior of the
system is naturally described on the “slow” τ -scale,
while the “fast” transients evolve according to the
rescaled time t :=
τ
ε
.
2.2. Propagator Expressions
In this section, we collect a number of analytical
results that underly the present case study; details
can be found in [2, 14].
2

Journal of Coupled Systems and Multiscale Dynamics
2.2.1. Zeroth-Order Propagator
The zeroth-order propagator for the two-stage
gene expression model (Fig. 1A) represents an ap-
proximation to the CME, Eq. (1), under the assump-
tion of infinitely fast mRNA degradation. Math-
ematically speaking, it is obtained in the singular
limit of γ → ∞, i.e., of ε → 0. Following [2],
we have
P
n|n
0
(τ, 0) = (1 − e
−τ
)
n
0

1 + be
−τ
1 + b

a

b
1 + b

n
n
X
k=0

(−1)
k
k!(n − k)!
Γ(a + n − k)Γ(n
0
+ 1)
Γ(a)Γ(n
0
− k + 1)
×
h
1 + b
b(1 − e
τ
)
i
k
2
F
1

− n + k, −a, 1 − a − n + k,
1 + b
e
τ
+ b


(2)
for the zeroth-order marginal probability P
n|n
0
(τ, 0)
of observing n protein molecules after time τ , given
m
0
= 0 molecules of mRNA and n
0
molecules of
protein initially. Here,
2
F
1
(a, b, c, z) is the Gauss
hypergeometric function [20]. We remark that, by
construction, P
n|n
0
(τ, 0) neglects any contributions
from the fast t-scale, as the decay of mRNA is in-
stantaneous to leading order in ε.
2.2.2. Uniform (First-Order) Propagator
The uniform propagator, denoted P
n|n
0
(τ, t, ε),
was derived as in [14]. Here, ε denotes the per-
turbation parameter, as before, while t is the fast
time variable. We emphasize that P
n|n
0
describes
the probability of transitioning from n
0
protein
molecules initially to n molecules at time τ = εt,
uniformly on the two time-scales. After some alge-
braic rearrangement, we find
P
n|n
0
(τ, t, ε) = P
n|n
0
(τ, ε)
+ εa

b
1 + b

n−n
0
1
(1 + b)
2
× [n − n
0
− b − (1 + b)t] +
εa
Γ(n − n
0
+ 2)
(bt)
n−n
0
t
×

1
F
1
(n − n
0
+ 1, n − n
0
+ 2, −(1 + b)t)t

1 −
n − n
0
− b
1 + b

+
n − n
0
+ 1
1 + b
e
−(1+b)t

(3)
to first order in ε; here,
1
F
1
(a, b, z) is the Kum-
mer function of the first kind (or confluent hyper-
geometric function) [20]. We remark that the transi-
tion probability P
n|n
0
(τ, ε) contributes on the slow
τ-scale in Eq. (3), while the t-dependent contribu-
tion in Eq. (3) accounts for the transient dynamics
on the fast time-scale.
Specifically, P
n|n
0
(τ, ε) denotes the marginal
probability, up to and including O(ε)-terms, of ob-
serving n protein molecules after time τ given m
0
=
0 molecules of mRNA and n
0
molecules of protein
initially:
P
n|n
0
(τ, ε) =
∞
X
m=0
P
m,n|0,n
0
(τ, ε) (4)
As shown in [14], the probability of encountering
more than 1 molecule of mRNA at time τ is negligi-
ble to first order in ε; thus, Eq. (4) reduces to
P
n|n
0
(τ, ε) = P
0,n|0,n
0
(τ, ε) + P
1,n|0,n
0
(τ, ε).
(5)
3

Journal of Coupled Systems and Multiscale Dynamics
After some algebraic simplification, the two tran-
sition probabilities P
0,n|0,n
0
and P
1,n|0,n
0
in the
above relation are found to be
P
0,n|0,n
0
(τ, ε) = (1 − e
−τ
)
n
0

b
1 + b

n

1 + be
−τ
1 + b

a
×
n
X
k=0
1
(n − k)B(a, n − k)
2
F
1

− n + k, −a, 1 − a − n + k,
1+b
e
τ
+b

×

g(n
0
, k) −
ε
2
a
(1 + b)
2
(k + 1) ×
h
2
F
1

− k, −n
0
, −1 − k,
1+b
b(1−e
τ
)

+

1 + b
e
τ
+ b

k+2
e
2τ
2
F
1

− k, −n
0
, −1 − k,
e
τ
+b
b(1−e
τ
)

i

, with
g(n
0
, k) =
(
0 for k > n
0
(−1)
k

n
0
k

(1+b)
b(1−e
τ
)

k
for k ≤ n
0
;
(6)
P
1,n|0,n
0
(τ, ε) = aε

b
b + 1

n
1
b + 1
(1 − e
−τ
)
n
0

1 + be
−τ
1 + b

a
×
n
X
k=0

1
(n − k)B(a, n − k)
2
F
1

k − n, −a, −a + k − n + 1,
b+1
e
τ
+b

×

h(n
0
, k) + (−1)
n
0
h
be
τ
+ 1
b(1 − e
τ
)
i
n
0

, with
h(n
0
, k) =
(
0 for k ≥ n
0

n
0
k+1

b+1
b(1−e
τ
)

k+1
2
F
1

1, k − n
0
+ 1, k + 2,
b+1
b(1−e
τ
)

for k < n
0
.
(7)
Here, B(a, b) :=
Γ(a)Γ(b)
Γ(a+b)
is the Beta function, with
the proviso that
1
(n−k)B(a,n−k)
= 1 when n = k.
Finally, Eq. (3) can be simplified by substituting
1
F
1
(n − n
0
+ 1; n − n
0
+ 2; −(1 + b)t) = [(1 + b)t]
−(n−n
0
+1)
Γ(n − n
0
+ 2)
− (n − n
0
+ 1)Γ

n − n
0
+ 1, (1 + b)t

(8)
to achieve the computationally more tractable for-
mulation
P
n|n
0
(τ, t, ε) = P
n|n
0
(τ, ε)
+ εa

b
1 + b

n−n
0
1
(1 + b)
2
[n − n
0
− b − (1 + b)t] + εat

−

b
1 + b

n−n
0
1
(1 + b)t
×

b + n
0
− n
1 + b
− t


1 − Q(n − n
0
+ 1, (1 + b)t)

+
(bt)
(n−n
0
)
1 + b
e
−(1+b)t
Γ(n − n
0
+ 1)

. (9)
4