domain: info math

text:
Vincent Danos

Tobias Heindel

Ilias Garnier

Jakob Grue Simonsen

Continuous-time Markov chains as transformers of unbounded observables

This paper provides broad sufficient conditions for the computability of time-dependent averages of stochastic processes of the form f (Xt) where Xt is a continuous-time Markov chain (CTMC), and f is a real-valued function (aka an observable). We consider chains with values in a countable state space S, and possibly unbounded f s. Observables are seen as generalised predicates on S and chains are interpreted as transformers of such generalised predicates, mapping each observable f to a new observable Ptf defined as (Ptf )(x) = Ex(f (Xt)), which represents the mean value of f at time time t as a function of the initial state x. We obtain three results. First, the well-definedness of this operator interpretation is obtained for a large class of chains and observables by restricting Pt to judiciously chosen rescalings of the basic Banach space C0(S) of S-indexed sequences which vanish at infinity. We prove, under appropriate assumptions, that the restricted family Pt forms a strongly continuous operator semigroup (equivalently the time evolution map t → Pt is continuous w.r.t. the usual topology on bounded operators). The computability of the time evolution map follows by generic arguments of constructive analysis. A key point here is that the assumptions are flexible enough to accommodate unbounded observables, and we give explicit examples of such using stochastic Petri nets and stochastic string rewriting. Thirdly, we show that if the rate matrix (aka the q-matrix) of the CTMC is locally algebraic on a subspace containing f , the time evolution of projections t → (Ptf )(x) is PTIME computable for each x. These results provide a functional analytic alternative to Monte Carlo simulation as test bed for mean-field approximations, moment closure, and similar techniques that are fast, but lack absolute error guarantees.

Introduction

The study of properties of stochastic processes is a currently highly active field of research in computer science, but has been avidly studied in mathematics even prior to Kolmogorov's development of the axiomatic approach to probability. For the particular case of continuous-time Markov chains (CTMCs), the crucial property of interest is how it acts on a function from its underlying state space to the reals -a "generalised predicate" in analogy to the usual notion of predicate transformers in computer science [Koz83]. Here, we are interested in the problem of computing time-dependent expectations of real-valued functions of the state space (so-called observables) on a continuous-time Markov chain X t . Ultimately, one wants provably correct, practical algorithms that given some continuous-time Markov chain X t and some observable f computes the transient conditional mean t → (x → E x (f (X t ))), i.e., the average of f conditioned on the initial state as a function of time. However, it is not a priori clear that this function is even computable let alone practically so; hence, it must be established that there is an algorithm that computes the transient conditional mean to arbitrary desired precision. As the basic objects are typically (computable) real numbers, we will employ the type 2 theory of effectivity (TTE) for computability.

It is useful to recall the nature of the objects that we manipulate and want to compute in the finite case. A continuous-time Markov chain (CTMC) on a finite state space S is entirely captured by its q-matrix, which is an S × Sindexed real matrix in which every row sums to zero and all negative entries lie on the diagonal. Any finite q-matrix Q induces a semigroup t → P t = e tQ which describes the time evolution of the CTMC. For any π on S (viewed as a row vector), the map t → πP t corresponds to the probability of being at any given state at time t for the initial probabilistic state π. Therefore, if X 0 is distributed according to π, the associated stochastic process X t is distributed according to πP t . Dually, one can interpret a CTMC as a transformer of observables, using the discrete time setting [Koz83] as an analogy. For any f : S → R (seen as a column vector), P t f = (x → y∈S p t,xy f (y)) is the vector of conditional means 3 of f at time t as a function of the initial state. The function we seek to compute is precisely t → P t f . Formally, P t f is the (column) vector of conditional means (i.e., (P t f )(x) = E x (f (X t )) for x ∈ S) where X t is the CTMC with transition function P t . For the finite state case, t → P t f = e tQ f is the unique solution to the initial value problem

d dt u t = Qu t u 0 = f (1)
and this can be solved by any method for computing matrix exponentials or numerical solvers for finite ODEs. In the countably infinite state case, observables are possibly unbounded-this makes computing transient conditional means fundamentally harder than in the finite case and calls for more sophisticated mathematics, as we shall describe in § 3.2.

For motivation, we give two paradigmatic examples of transient means. The classic example of a stochastic process of the form f (X t ), i.e., a pair of a CTMC X t and an observable f , where X t models a set of chemical reactions with the molecule count of a certain chemical species as observable: states of X t are multisets over a finite set of species; then, we are interested in computing the evolution of the mean count of a certain species.

3 By definition, Ex(f (Xt)) = ω∈Ω f (Xt(ω)) dpx for some suitable probability space (Ω, p), where px is the probability measure p on Ω conditioned by X0 = x. Performing a change of variables, Ex(f (Xt)) = S f dX * t (px), where X * t (px) is the image measure of px through Xt. Since pt,x-= X * t (px)(-), we deduce Ptf (x) = Ex(f (Xt)).

Another natural example native to computer science is a stochastic interpretation of any string rewriting system as a CTMC X t . An obvious class of observables for string rewriting are functions f that count the occurrence of a certain word as sub-string in each state of the CTMC X t ; note that this is different from counting "molecules", as there is always only a singe word! For example, consider the string rewriting system with the single rule a aba and initial state a: the mean occurrence count of the letter a grows as the exponential function e t while the mean occurrence count of the word aa is zero at all times; adding the rule ba ab does not change the mean a-count but renders the mean count of the word aa non-trivial.

Note that these two classes of models are only meant as simple examples and the results of this paper are not geared towards any particular modelling language for CTMCs; the results apply equally well to Kappa models [DFF + 10] or stochastic graph transformation. 4Our main result is a sufficient condition for E x (f (X t )) to be computable, powerful enough to encompass many interesting unbounded observables. We proceed by construction of a suitable computable Banach space where to restrict our operators P t , in such a way that they form a strongly continuous semigroup (SCSG). It follows that the initial value problems corresponding to observables in the domain of the generator of the semigroup admit computable solutions. The construction of our Banach space rests on functional analytic techniques recently developed in Ref. [Spi12,Spi15], which we combine with the work on SCSGs by Weihrauch & Zhong [WZ07].

We do not assume that the reader has some acquaintance with continuoustime Markov chains on a countable state space and recall the basic concepts of transition functions and q-matrices (Def. 5 and 6) in the preliminaries. Having said that, and both in order to save space and to keep the logic of the paper clearly apparent, most of the mathematical material used is postponed to a series of appendices.

Two motivating examples of CTMCs with observables

We illustrate our constructions with: (i) chemical reaction networks (CRN), aka stochastic Petri nets, and (ii) stochastic string rewriting as a simple example of (rule-based) modelling. In both cases, the construction of the q-matrix implied by a model is readily done, and so is the definition of a natural set of unbounded observables with clear relevance to the dynamics of a model: word occurrence counts for stochastic string rewriting (Def. 2) and multiset inclusions for Petri nets (Def. 3).

Stochastic string rewriting and word occurrences

Stochastic string rewriting can be thought of as never ending, fair competition between all redexes of rules, "racing" for reduction; the formal definition is as follows, in perfect analogy to Ref.

[HLM06] which covers the case of graphs.

Definition 1 (Stochastic string rewriting). For each rule ρ = l r ∈ Σ + × Σ + we define the q-matrix of ρ, denoted by Q ρ , as the q-matrix Q ρ = (q ρ uv ) u,v∈Σ + on the state space of words Σ + with off-diagonal entries

q ρ uv = (w, w ) ∈ Σ * × Σ * u = wlw , v = wrw
for each pair of words u, v ∈ Σ + such that u = v, and diagonal entries,

q ρ uu = -v =u q ρ uv for all u ∈ Σ + . For a finite set of rules R ⊆ Σ + × Σ + , we define Q R = ρ∈R Q ρ ,

and with additional choices of rate constants

k : R Q + , we define Q R,k = ρ∈R k ρ Q ρ .
For a given rule set R, each entry q R uv of the q-matrix corresponds to the propensity to rewrite: it is just the number of ways in which u can be rewritten to v. We shall usually work without rate constants for the sake of readability. Note that the use of Σ + for the left and right hand side of rules is convenient to get string rewriting as a special case of graph transformation in a straightforward manner.

The occurrence counting function of a word as sub-string in the state of the CTMC of R is as follows.

Definition 2 (Word counting functions). Let w ∈ Σ + be a word. The wcounting function, denoted by w :

Σ + → R ≥0 , maps each word x ∈ Σ + to w (x) = |{(u, v) ∈ Σ * × Σ * | x = uwv}|.

Stochastic Petri nets and sub-multiset occurrences

We recall the definition of stochastic Petri nets and occurrence counting of a multisets. Note that for the purposes of the present paper, places and species are synonymous.

Definition 3 (Multisets and multiset occurrences).

A multiset over a finite set P of places is a function x : P → N that maps each place to the number of tokens in that place. Given a multiset, x ∈ N P , the x-occurrence counting function x : N P → N is defined by

x (y) = y! (y-x)! x ≤ y 0 otherwise
where z! = p∈P z(p)! is the multiset factorial for all z ∈ N P .

Definition 4 (Stochastic Petri net). Let P be a finite set of places. A stochastic Petri net over P is a set

T ⊆ N P × R >0 × N P
where N P is the set of multisets over P, which are called markings of the net; elements of the set T are called transitions. The q-matrix Q l,k,r on the set of markings for a transition (l, k, r) ≡ l → k r ∈ T has off-diagonal entries

q l,k,r xy = k • l (x) l (x) > 0, y = x -l + r 0 otherwise
where addition and subtraction is extended pointwise to

N P . The q-matrix of T is Q T = (l,k,r)∈T Q l,k,r .

Preliminaries

For the remainder of the paper, we fix an at most countable set S as state space.

Transition functions and q-matrices

We first recall the basic definitions of transition functions and q-matrices. We make the usual assumptions [And91] one needs to work comfortably: namely that q-matrices are stable and conservative and that transition functions are standard and also minimal as described at the end of § 3.1 (and detailed in Appendix A). With these assumptions in place, transition functions and q-matrices determine each other, and one can freely work with one or the other as is most convenient. z∈S p t,xz ≤ 1 for all x ∈ S and t ∈ R ≥0 . Thus, each row of a transition function corresponds to a sub-probability measure, and transition functions converge entry-wise to the identity matrix at time zero. Unless stated otherwise, all transition functions in this paper are standard.

Taking entry-wise derivatives at time 0 gives a q-matrix.

Definition 6 (q-matrix). A q-matrix on S is an S × S-matrix Q = (q xy ) x,y∈S with real entries q xy such that q xy ≥ 0 (if x = y), q xx ≤ 0, and z∈S q xz = 0 for all x, y ∈ S.

Conversely, for each q-matrix, there exists a unique entry-wise minimal transition function that solves Equation (2) [And91, Theorem 2.2],

d dt P t = QP t , P 0 = I (2)
which is called the transition function of Q. From now on, we assume that all transition functions are minimal solutions to Equation (2) for some q-matrix Q (cf. Definition 15, see also [Nor98,p. 69]).

The Abstract Cauchy problem for P t f

Abstract Cauchy problems (ACPs) in Banach spaces [Ein52], are the classic generalisation of finite-dimensional initial value problems. (See Appendices C and D for details.) Concretely, we want to obtain P t f as unique solution u t of the following generalisation of our earlier equation (1):

d dt u t = Qu t (t ≥ 0) u 0 = f (acp)
where f is an observable and Q is a linear operator which plays the role of the q-matrix.

There are a few points worth noting. First, the topological vector space of all observables R S cannot be equipped with a suitable complete norm. Therefore, one has to look for a Banach subspace B ⊂ R S wherein to interpret the above equation. There are several ways to do this, and they are not equally interesting. Second, as P t f = u t is the desired solution, and P 0 = I, it follows that d dt P t f | t=0 = Qf . If this derivative does not exist, Qf is simply not defined. In fact, as is clear from the examples, we can only expect Q to be partially defined and an unbounded operator.5 On the positive side, if P t is a strongly continuous semigroup (Def. 20) on B, meaning lim h 0 P h f = f for all f in B, and we take Q as its generator defined on a subspace of B (Def. 21) by the above formula, then P t f is the unique solution of (acp), provided f is in the domain of Q [EN00, Proposition II.6.2].

Even better, in this case, not only does (acp) have P t f as unique solution, but we get an explicit approximation scheme:

P t f = lim n→∞ e tAn f ( 3
)
where θ is a constant of the SCSG such that nI -Q is invertible for n > θ and

A n = nQ(nI -Q) -

Banach space wanted!

Table 1 gives an overview of initial value problems for transient distributions (first row) and transient conditional means (second row). Transient distributions are summable sequences, and transition functions form SCSGs [Reu57] and therefore allow for a well-posed corresponding ACP. But the classic example of a Banach space to reason about conditional means [RR72] is the space C 0 (S) of functions vanishing at infinity, i.e., functions f : S → R such that for all > 0, the set {x ∈ S | f (x) ≥ } is finite, equipped with the supremum norm (Def. 17). The corresponding processes are called Feller transition functions [And91, § 1.5] and verify a principle of finite velocity of information flow (for all t, y, p t,xy vanishes as x goes to infinity).

Spieksma's theorem

A solution is provided by a result of Spieksma [Spi12, Theorem 6.3], which constructs a Banach space for a given q-matrix Q and an observable f of interest. The formulation will be in terms of so-called drift functions.

Definition 7 (Drift function).

Let Q be a q-matrix on S, and let c ∈ R.

A function W : S → R >0 is called a c-drift function for Q if for all x ∈ S (QW )(x) := y∈S q xy W (y) ≤ cW (x).
We will say that W is a drift function for Q if there exists c ∈ R such that it is a c-drift function for Q. One can show that P t W ≤ e ct W in this case. Thus, drift functions control their own growth under the transition function.

We also need to work with weighted variants of C 0 (S).

Definition 8 (Weighted C 0 (S)-spaces). Let S be a set and let W : S → R >0 be a positive real-valued function, referred to as a weight. The Banach space C 0 (S, W ) consists of functions f : S → R such that f /W vanishes at infinity, where

( f /W )(x) = f (x) /W (x). The norm _ W on such functions f is f W = sup x∈S | f (x) /W (x)|.
As C 0 (S, W ) is isometric to C 0 (S) it is indeed a Banach space. It is also a closed subspace of L ∞ (S, W ), the set of functions such that f /W is bounded. We will use later the fact that: Lemma 1. Finite linear combinations of indicator functions,7 form a dense subset of C 0 (S) Indeed, suppose S = N for convenient notations, and define g n = f truncated after n. It is easy to see that g n converges to f in the _ W norm. In fact,

g n -f W = sup x>n |f (x)|/W (x) → 0 iff f /W vanishes at infinity. Now the idea is to find positive drift functions V, W for Q such that V ∈ C 0 (S, W ), i.e.
, such that V /W vanishes at infinity. Thus, both functions can grow at most exponentially in mean and V is negligible compared to W at infinity. Intuitively, functions on the order of V are as good as functions vanishing at infinity, in the case of Feller processes [RR72], i.e., CTMCs such that transition functions induce SCSGs on C 0 (S).

Theorem 1. Let P t be a transition function on S with q-matrix Q, let V, W : S → R >0 be drift functions for Q such that V ∈ C 0 (S, W ). Then:

1. P t induces an SCSG on C 0 (S, W ). 2. For all f ∈ C 0 (S, W ), P t f is given by Equation (3) in C 0 (S, W ) where Q is the generator of P t .
The first part of the theorem is proved in [Spi12, Theorem 6.3]; the second part follows from the general theory of ACPs. Note that f does not need to be in the domain of Q, in which case we only obtain a 'mild' solution to the ACP, i.e., a solution to its integral form which might not be everywhere differentiable. If on the other hand, f is in the domain of Q, then P t f is continuously differentiable.

We have now covered the mathematical ground needed to make sense of conditional means of not necessarily bounded observables. This, however, does not immediately yield an algorithm for computing transient means, as it does not yet ensure that transient means are computable even in presence of a complete specification of the underlying system. Even transient conditional probabilities can fail to be computable [AFR11]! Before we proceed to computability questions, let us return to our two classes of examples.

Applications: string rewriting and Petri nets

We now give examples of drift functions for stochastic string rewriting and Petri nets. The former case is well-behaved since the mean letter count grows at most exponentially. The case of Petri nets will be more subtle and we shall give an example of an explosive Petri net such that we can nevertheless reason about conditional means of unbounded observables.

For string rewriting, we have canonical drift functions.

Lemma 2 (Powers of length are drift functions). Let R ⊆ Σ + × Σ + be a finite string rewriting system and let n ∈ N + be a positive natural number. There exists a constant

c n ∈ R >0 such that |_| n : Σ + → R ≥0 is a c n -drift function.
Proof. The proof is deferred to Appendix E. Now, we can apply Spieksma's method to get a Banach space for reasoning about conditional means and moments of word counting functions.

Corollary 1 (Stochastic string rewriting).

Let R be a finite string rewriting system, let n ∈ N \ {0}, and let |_| : Σ + N be the word length function. The transition function

P t of q-matrix Q R is a SCSG on C 0 (Σ + , |_| n ).
Thus, all higher conditional moments of word counting functions can be accommodated in a suitable Banach space. The case of Petri nets is more subtle, since, in general, the (weighted) token count is not a drift function.

Example 1. Consider the Petri net with single transition 2A → 1 3A and with one place A. The token count A is not a drift function. In fact, the corresponding CTMC is explosive (by [Spi15, Theorem 2.1]).

Our final example is an extension of the previous explosive CTMC with a new species whose count can nevertheless be treated using Theorem 1.

Example 2 (Unobserved explosion). Consider the Petri net with transitions

{2A → 1 3A, B → 1 2B}.
The underlying CTMC is explosive, and we cannot apply Theorem 1 to compute the transient conditional mean of the A-count for the exact same reason as in Example 1. However, we can do so for the B-count, using the weight function W = B 2 and observable f = B . Putting V = f allows to apply Spieksma's recipe (ruling out states with B-count 0 for convenience). The conditional mean E 2A+B ( B (X t )) can be best understood by adding a coffin state, on which both the A-and B-count are zero and in which the Markov chain resides after (the first and only) explosion.

Computability

We follow the school of type-2 theory of effectivity, as reviewed in Appendix G following [BHW08]. A real number x is computable iff there is a Turing machine that on input d ∈ N (the desired precision), outputs a rational number r with |r -x| < 2 -d . Next, a function g : R → R is computable if there is a Turing machine that, for each x ∈ R, takes an arbitrary Cauchy sequence with limit x as input and generates a Cauchy sequence that converges to g(x)-where convergence has to be sufficiently rapid, e.g., by using the dyadic representation of the reals.

Computability extends naturally to any Banach space B other than R. We only need a recursively enumerable dense subset on which the norm, addition and scalar multiplication are computable, thus making B a computable Banach space (Definition 30); usually, the dense subset is induced by a basis of a dense subspace. For weighted C 0 -spaces (with computable weight functions) and their duals, we fix an arbitrary enumeration of all rational linear combinations of indicator functions 1 x ; for the Banach space of bounded linear operators on weighted C 0 -spaces we use the standard construction for continuous function spaces [WZ07, Lemma 3.1]. A SCSG P t is computable if the function t → P t from the reals to the Banach space of bounded linear operators is computable. The computable SCSGs correspond to those obtained from CTMCs through Theorem 1.

For simplicity, motivated by the observation that we do not lose any of the intended applications to rule-based modeling, we restrict to row-and columnfinite q-matrices with rational entries in our main result.

Theorem 2 (Computability of CTMCs as observation transformers).

Let Q be a q-matrix on S, let V, W : S → R ≥0 be positive drift functions for Q such that V ∈ C 0 (S, W ). If the q-matrix Q is row-and column-finite, consists of rational entries, and is computable as a function Q : S 2 → Q, and -W : S → Q is computable, the following hold.

1. The SCSG P t , i.e., the function t → P t , is computable.

The evolution of conditional means

(t, f ) → P t f is computable. 3. The evolution of means (π, t, f ) → πP t f is computable for all π : S → R such that x∈S W (x)|π(x)| < ∞.
Proof. We will apply a result by Weihrauch & Zhong on the computability of SC-SGs [WZ07, Theorem 5.4]. Applying this result requires some extra information:

1. the SCSG P t must be bounded in norm by e θt for some positive constant θ; 2. we must have a recursive enumeration of a dense subset of the graph of the infinitesimal generator of the SCSG P t .

We first show that the constant θ, featuring in Theorem 1, i.e., the witness that W is a θ-drift function for Q, satisfies P t ≤ e θt (using the first part of the proof of Theorem 6.3 of Ref. [Spi12]). Next, we obtain a recursive enumeration of a dense subset A ⊆ Q of the domain of the generator Q by applying Q to all rational linear combinations of indicator functions 1 x . Note that for the latter, we use that indicator functions belong to the domain of the generator and

Q1 x = Q1 x (see Proposition 2).
Now, by Theorem 5.4.2 of Ref.

[WZ07], we obtain the first two computability results, as θ, A, 1 is a so-called piece of type IG-information [WZ07, p. 513]. Finally, the third point amounts to showing computability of the duality pairing

_, _ : L 1 (S, W ) × C 0 (S, W ) → R where L 1 (S, W ) is (isomorphic to) the dual of C 0 (S, W ) (Lemma 3).
This theorem immediately gives computability of the CTMCs and (conditional) means for stochastic string rewriting and Petri nets discussed in Sect. 4. Note that the theorem does not assume that V itself is computable; its role is to establish that the transition functions is a SCSG, but V plays no role in the actual computation of the solution. Note also that the algorithms that compute the functions t → P t , (t, f ) → P t f , and (π, t, f ) → πP t f push the responsibility to give arbitrarily good approximations of the respective input parameters π, t and f to the user. This however is no problem for any of our examples or rulebased models in general: t is typically rational, f is computable and even to the natural numbers, and π is often finitely supported or a Gaussian.

Remark 1 (Domain of the generator). The difficulty in working with SCSGs is to find a tractable description of the domain of their generator. Such a description is important for two reasons: first, as said before the theorem, for computability, we need a recursively enumerable dense subset of the graph of the domain of the generator; second, as the unique solution to the ACP is also differentiable, it is useful for numerical integration.

Computability ensures existence of algorithms computing transient means, but yields no guarantees of the efficiency of such algorithms. We now proceed to a special case that (i) encompasses a number of well-known examples, including context-free string rewriting, and (ii) leads to PTIME computability, by reducing the problem of transient conditional means to solving finite linear ODEs.

The finite dimensional case and PTIME via ODEs

We now turn to the special case where we can restrict to finite dimensional subspaces B ⊆ R S . 8 The prime example will be word counting functions and context-free string rewriting systems. Hyperedge replacement systems [DKH97], the context-free systems of graph transformation, can be handled mutatis mutandis. The main result is PTIME computability of conditional means. For convenience, we extend the usage of the term locally algebraic as follows.

Definition 9 (Locally algebraic). We call a q-matrix Q on S locally algebraic for an observable f ∈ R S if the set {Q n f | n ∈ N}, containing all multiple applications of the q-matrix Q to the observable f , is linearly dependent, i.e., if there exists a number N ∈ N such that the application Q N f of the N -th power is a linear combination

N -1 i=0 α i Q i f of lower powers of Q applied to f .
Using local algebraicity of a q-matrix Q for an observable f , one can generate a finite ODE with one variable for each conditional mean

E[Q n f (X t ) | X 0 = x]
(as detailed in the proof of Theorem 3); then, recent results from computable analysis [PG16] entail PTIME complexity.

Theorem 3 (PTIME complexity of conditional means).

Let Q be a qmatrix on S, let x ∈ S, let f : R → S be a function such that f (x) is a computable real and

Q N f = N -1 i=0 α i Q i f for some N ∈ N and computable coefficients α i .
The time evolution of the conditional mean P t f (x), i.e., the function t → P t f (x), is computable in polynomial time.

Proof. Consider the N -dimensional ODE with one variable E n for each n ∈ {0, . . . , N -1} with time derivative

d dt E n (t) = E n+1 (t) if n < N -1 N -1 i=0 α i E i (t) if n = N -1 and initial condition E n (0) = Q n f (x). Solving this ODE is in PTIME [PG16]
(even over all of R ≥0 , using the "length" of the solution curve as implicit input).

Finally,

E i (t) = P t Q i f (x) is the unique solution.
Note that the linear ODE that we construct has a companion matrix as evolution operator, which allows to use special techniques for matrix exponentiation [TR03,BCO + 07].

Proposition 1 (Local algebraicity of context-free string rewriting). Let R be a string rewriting system, let w ∈ Σ + , let m ∈ N. The q-matrix Q R of the string rewriting system R is locally algebraic for the m-th power of w-occurrence

counting w m if R ⊆ Σ × Σ + .
Proof. For every product of word counting functions m i=1 wi , applying the q-matrix Q R to this product yields the observable Q R m i=1 wi . Using previous work on graph transformation [DHHZS15], restricted to acyclic, finite, edge labelled graphs that have a unique maximal path (with at least one edge), the observable Q R m i=1 wi is a linear combination k j=1 α j kj l=1 w j,l of word counting functions w j,l (with all

k j ≤ m). Moreover, if R is context-free (R ⊆ Σ × Σ + ), we have kj l=1 |w j,l | ≤ m i=1 |w i |
for all j ∈ {1, . . . , k}. Thus, we stay in a subspace that is spanned by a finite number of products of word counting functions.

Corollary 2. For context-free string rewriting, conditional means and moments of word occurrence counts are computable in polynomial time.

We conclude with a remark on lower bounds for the complexity.

Remark 2. The complexity of computing transient means, even for context-free string rewriting, is at least as hard as computing the exponential function. This becomes clear if we consider the rule a → aa, and the observable of a-counts a . Now, the time evolution of the a -mean conditioned on the initial state to be a, i.e., the function t → E a ( a (X t )), is exactly the exponential function e t . Tight lower complexity bounds for the exponential function are a longstanding open problem [Ahr99].

The main result is computability of transient (conditional) means of Markov chains X t "observed" by a function f , i.e., stochastic processes of the form f (X t ). For this, we have described conditions under which a CTMC, specified by its q-matrix, induces a continuous-time transformer P t that acts on observation functions. In analogy to predicate transformer semantics for programs, this should be called observation transformer semantics for CTMCs; formally, P t is a strongly continuous semigroup on a suitable function space. Finally, motivated by important examples of context-free systems -be it the well-known class from Chomsky's hierarchy or the popular preferential attachment process (covered in previous work [DHHZS15]) -we have considered the special case of locally finite q-matrices. For this special case, we obtain a first complexity result, namely PTIME computability of transient conditional means.

The obvious next step is to implement our theoretical results since one cannot expect that the general algorithms of Weihrauch & Zhong [WZ07] perform well for every SCSG on a computable Banach space. For example, the Gauss-Jordan algorithm for infinite Matrices [Par12] should already be more practicable for inverting the operator nI -Q from Equation (3) compared to the brute force approach used by Weihrauch & Zhong [WZ07]. Computability ensures existence of algorithms for computing transient means, but yields no guarantees of the efficiency of such algorithms.

Even if it should turn out that efficient algorithms are a pipe dream -after all, transient probabilities p t,xy are a special case of transient conditional meanswe expect that already implementations that are slow but to arbitrary desired precision will be useful for gauging the quality of approximations of the "meanfield" of a Markov process, especially in the area of social networks [Gle13], but possibly also for chemical systems [SSG15]. Theoretically, they are a valid alternative to Monte Carlo simulation, or even preferable. A From Transition functions to q-matrices and back The first two properties in Definition 10 make transition functions monoids of real-valued S × S-matrices equipped with matrix multiplication; the last property requires that rows induce sub-probability measures. 10 Thus, an honest transition function P t consists of proper probability matrices, i.e., the entries of every row sum up to 1. Unless explicitly allowed to be non-standard, we shall assume transition functions to be standard, i.e., entry-wise continuous at time 0.

9 The reader familiar with CTMCs might object that in shifting focus from Markov chains to pairs of initial distributions and transition functions, we lose information. However, "[a]ll the probabilistic information about the process, insofar as it concerns only countably many time instants, is contained in the transition function and initial distribution. One could almost say that the transition function is the Markov chain." [And91, p. 4]. 10 A sub-probability measure on an at most countable set S (equipped with the σalgebra of all subsets) can be safely identified with a function to the non-negative reals f : S → R ≥0 that sums up to at most 1, i.e., such that x∈S f (x) ≤ 1 holds.

Example 3 (Yule-Furry process). The Yule-Furry process of parameter β > 0, the archetypal example of a branching process, can be thought of as a simplistic model of population growth; for n 0 ≤ n, its transition probabilities have the following closed form.

p t,n0 n = n -1 n -n 0 (e -βt ) n0 (1 -e -βt ) n-n0
Using the metaphor of population growth, if we start with a population of size 1, the probability at time t to have grown to exactly size n ≥ 1 is e -βt (1 -e -βt ) n-1 .

Typically, transition functions of CTMCs do not have simple closed expressions for their entries; the Yule-Furry process is one of the few exceptions. In contrast, taking (time) derivatives of transition functions (at time zero) typically yields natural functions, intuitively corresponding to a flow of probability mass; for a Yule process of parameter β with initial population n 0 , we think of n 0 independent individuals that can replicate at any given moment in time, and thus the outflow of probability mass is -βn 0 .

Transition functions are always entry-wise differentiable [Kol51,Aus55] and the derivatives at time zero p 0,xy form the so-called q-matrix of the transition function [And91,p. 14]. Provided that all states are stable, i.e., p 0,xx = -∞ for all x ∈ S, q-matrices of transition functions can safely be identified with stable q-matrices according to the following definition [And91, Theorem 2.2].

Definition 11 (q-Matrix [And91, p. 64]). A q-matrix is an S × S-matrix Q = (q xy ) x,y∈S with entries real numbers or -∞ such that

∞ > q xy ≥ 0 (if x = y),
-∞ ≤ q xx ≤ 0, and z∈S q xz ≤ 0 hold for all x, y ∈ S; it is stable (resp. conservative) if q x := -q xx < ∞ (resp.

z∈S q xz = 0) holds for all x ∈ S. In this paper, we assume q-matrices to be stable and conservative. For each q-matrix Q = (q xy ) x,y∈S on an at most countable state space S, the Kolmogorov backward equation

d dt P t = QP t , P 0 = I (4)
has a unique minimal non-negative solution P t = (p t,xy ) x,y∈S (see, e.g., Theorem 2.2 of Ref.

[And91]), which is a possibly dishonest transition function.

Definition 12 (Transition functions of q-matrices). For a stable q-matrix Q on S, the transition function of q-matrix Q, is the unique minimal transition function P t solving the Kolmogorov backward equation, i.e., any other transition function Pt = ( pt,xy ) x,y∈S that solves Equation (4) satisfies pt,xy ≥ p t,xy for all (t, x, y) ∈ R ≥0 × S × S.

In case of a finite state space, all the above is trivial: the transition function P t is just the matrix exponential e Qt , and similarly, if the entries of the q-matrix are bounded. In general, the hard part is to find for a given q-matrix a corresponding transition function, while the opposite direction, passing from a transition function to its q-matrix, is always via differentiation [Kol51,Aus55]; however, with the usual assumptions on CTMCs (as in Norris's standard textbook [Nor98]), qmatrices are stable and conservative, which allows to pass back and forth between transition functions and q-matrices without complications [And91, Theorem 2.2].

B Markov chains

For completeness sake, we recall the classic definitions of continuous-time Markov chains in terms of a stochastic process satisfying the Markov property. However, note that all proofs rely on the functional analytic perspective described in Appendix A, which does not hinge on any non-trivial measure theory. The following definitions are based on [And91].

Definition 13 (Continuous-time Markov chain (CTMC)).

A stochastic process {X t } t∈R ≥0 , defined on a probability space (Ω, F, Pr), with values in a countable set S (to be called the state space of the process), is called a continuoustime Markov chain (CTMC) if for any finite set

0 ≤ t 1 ≤ t 2 ≤ • • • ≤ t n ≤ t n+1
of "times", and corresponding set i 1 , i 2 , . . . , i n-1 , i, j of states in S we have

Pr X tn+1 = j X tn = i, X tn-i = i n-1 , . . . X t1 = i 1 = Pr X tn+1 = j X tn = i whenever Pr X tn = i, X tn-i = i n-1 , . . . X t1 = i 1 > 0.
If for all 0 ≤ s ≤ t and all i, j ∈ S the equation

Pr X t = j X s = i = Pr X t-s = j X 0 = i
holds, the stochastic process {X t } t∈R ≥0 is homogeneous.

As described by Anderson [And91,p. 4], all probabilistic information about a CTMC is given by their transition function.

Definition 14 (Transition function of a CTMC).

Let {X t } t∈R ≥0 be a homogeneous CTMC. The t-dependent matrix P t = (p t,ij ) i,j∈S with

p t,ij = Pr X t = j | X 0 = i for all i, j ∈ S and t ∈ R ≥0 is the transition function of the CTMC.
We only consider homogeneous CTMCs that are stable, standard11 and minimal (cf. Assumption 2.1 of Ref. [Spi12]).

Definition 15 (Stable, standard, minimal). Let {X t } t∈R ≥0 be a homogeneous CTMC. It is stable if the q-matrix of its transition function is stable, it is standard if its transition function is standard, and it is minimal if its transition function is the minimal transition function that solves the Kolmogorov (backward) equation (4).

C Linear operators on Banach spaces

In the present paper, we only consider vector spaces over the field R. We recall some core definitions (see also [EN00, Appendices A and B]).

Definition 16 (Banach space). A Banach space is a complete normed vector space.

We now give examples of Banach sequence spaces, i.e., Banach spaces consisting of functions f : S → R from a fixed set S, which will always be the state space of some CTMC is this paper.

Definition 17 (C 0 -spaces). A function f : S → R vanishes at infinity if for all > 0, f (x) < holds for almost all x ∈ S, i.e., the set {x ∈ S | f (x) ≥ } is finite. The Banach space C 0 (S) is the Banach space that consists of all functions f : S → R that vanish at infinity and its norm _ is given by f

= sup x∈S |f (x)|.
The duals of C 0 -spaces exists and are isomorphic to L 1 -spaces; however, C 0 -spaces are in general not isomorphic to their double dual (which is isomorphic to the corresponding L ∞ -space).

Definition 18 (L 1 -spaces). A function

f : S → R is summable if the sum x∈S |f (x)
| converges (for some order of summation). The Banach space L 1 (S) consists of summable functions and its norm _ is given by f

= x∈S |f (x)|.
Definition 19 (Linear operators, bounded linear operators). Let B be a Banach space. A linear operator on B is a pair (A, D(A)) where D(A) ⊆ B is a linear subspace, called domain, and A is a linear map A : D(A) → B. A linear operator on B is bounded if there exists a constant K ∈ R ≥0 such that A(u) ≤ K u holds for all u ∈ D(A). We denote by L(B) the set of all linear operators (A, D(A)) on B that are bounded and defined on the whole space, i.e., satisfy D(A) = B.

D Abstract Cauchy problems of SCSGs

Every transition function P t = (p t,xy ) x,y∈S (of some CTMC X t with state space S) induces a partial function on the vector space R S that maps a function f : S → R to the function P t f := (x → y∈S p t,xy f (y)) (of the same "type" S → R) whenever the respective sums converge absolutely for all states x ∈ S, and it is undefined on f otherwise. This function is indeed always a linear operator on the vector space of bounded functions [And91, Sect. 1, Lemma 4.4.]. However, we are interested in a stronger property (Theorem 1).

Definition 20 (Strongly continuous semigroup).

Let B be a Banach space, and let L(B) be the set of bounded linear operators on B with B as domain of definition. A strongly continuous semi-group (SCSG) on B is a family {P t } t∈R ≥0 of bounded linear operators P t ∈ L(B) satisfying 1. P 0 = I B , the identity on B; 2. P t+s = P t P s for all s, t ∈ R ≥0 ; and 3. lim h 0 P h f -f = 0 holds for all f ∈ B. 12In the finite state case, each transition functions P t can be obtained from their q-matrix Q -simply by taking the matrix exponential P t = e tQ . In general, a SCSG is only the limit of a sequence of matrix exponentials of bounded operators, the so-called Yosida-approximants (see, e.g. Ref.

[EN00]). Part of the intricacy is explained by the fact that the time derivative d dt P t f at 0, i.e., the limits lim h 0 1 /h(P h f -f ), do not exists for all elements f of the Banach space. This brings us to a new object, which is typically related but in general not identical to the q-matrix (cf. [Reu57, Theorem 2]).

Definition 21 (Infinitesimal generator). The infinitesimal generator Q of a strongly continuous semigroup P t on a Banach space B is defined by

Qf = lim h 0 1 /h(P h f -f ) for all f ∈ dom(Q) = f ∈ B The limit lim h 0 1 /h(P h f -f ) exists. .
Infinitesimal generators determine the semigroup uniquely [EN00, Theorem 1.4]. To understand why, note first that it is sufficient to describe each P t on a dense subset of the Banach space; the domain of definition of the infinitesimal generator is simply the canonical such subset. Finally, for each f in the domain of the domain of the generator, the function t → P t f is characterised as the unique differentiable solution to the abstract Cauchy problem [Ein52] in Equation (5).

d dt u t = Qu t (t ≥ 0) (5) u 0 = f

E Drift functions for observables

Powers of length are drift functions (Lemma 2)

Let R ⊆ Σ + × Σ + be a finite string rewriting system and let n ∈ N + be a positive natural number. There exists a constant c n ∈ R >0 such that

|_| n : Σ + → R ≥0 is a c n -drift function.
Proof. Let K = max l r∈R |r|; it is a (rough) bound for the maximal length "growth" of rules in R. Define c n = (2K) n |R|. For all words v ∈ Σ + and all rules l r ∈ R,

(Q l r |n|)(v) = w∈Σ + q l r vw |w| n (by definition of Q l r |n|) = w∈Σ + \{v} q l r vw |w| n -   w∈Σ + \{v} q l r vw   q l r v |v| n (by conservativity) = w∈Σ + \{v} q l r vw (|w| n -|v| n ) ≤ w∈Σ + \{v} q l r vw ((|v| + K) n -|v| n ) (by definition of K) ≤ ((|v| + K) n -|v| n ) w∈Σ + \{v} q l r vw ≤ |v|((|v| + K) n -|v| n ) (≤ |v| redexes of l r in v) = |v| n k=1 n k |v| n-k K k (by Binomial Theorem and cancelling -|v| n ) ≤ |v| n k=1 n k |v| n-k K k ≤ |v| n k=1 n k |v| n-k K k |v| k-1 K n-k = |v| n k=1 n k |v| n-1 K n ≤ 2 n K n |v| n
The desired follows since Q R = l r∈R Q l r .

F Spieksma, corrigendum

The statement of Theorem 6.3 of Ref.

[Spi12] contains an error (corrected in our version) that can be spotted if one tries to fill in the details that are left out at the end of Spieksma's proof in op. cit. A counterexample is given in Appendix I.

The following proposition gives a useful and correct weakening of her statement, on a non-empty subset of the domain of the generator of the relevant SCSGs.

Proposition 2. Let P t be a transition function on S with q-matrix Q such that there exists a positive function W : S → R >0 and a constant θ such that QW ≤ θW , (i.e., (QW )(x) ≤ θW (x) for all x ∈ S). Let the family {P t } t∈R ≥0 be a SCSG on C 0 (S, W ), witnessed by V ∈ C 0 (S, W ) that is positive and a drift function for Q. For all f ∈ C 0 (S, W ) that satisfy f V < ∞ and Qf ∈ D(Q), we have

Qf = Qf = lim h 0 1 /h(P t f -f ), (6)
i.e., the latter limit exists in C 0 (S, W ), and in particular f ∈ D(Q).

G Computable analysis

We collect the main definitions and terminology from [BHW08], as required in [WZ07].

G.1 Naming schemes for sets

The central point of computable analysis in the tradition of the type-2 theory of effectivity (aka tte) are so-called representations of elements of sets (that later will always come equipped with a topology).

Definition 22 (Notations, representations, naming systems). Let X be a set, and let Σ be finite alphabet. A notation of X is a surjective partial function ν :⊆ Σ * → X, and, similarly, a representation is a surjective partial function δ :⊆ Σ ω → X; finally, a naming system for X is a notation or a representation of X. Given an element x ∈ X and a naming system of X, a γ-name is an element of γ -1 ({x}) (which thus is a finite or infinite string over Σ that is mapped to x).

G.2 Computable reals

There are several equivalent definitions of computable real numbers [BHW08, Theorem 3.2].

Definition 23 (Computable real [BHW08, Definition 3.1]). A real number

x ∈ R is computable if {q ∈ Q | q < x}
is a decidable set (of rational numbers).

Definition 24 (Computable sequence of reals [BHW08, Definition 3.8]).

A sequence (x n ) n∈N of real numbers is called computable if there exists a computable sequence (q k ) k∈N of rational numbers such that |x n -q n,i | < 2 -i for all n, i ∈ N.

Definition 25 ([BHW08, Definition 3.11.]).

1. If (r n ) n∈N is a convergent sequence of real numbers with limit x, and m :

N → N is a function such that i ≥ m(n) implies |x -r i | < 2 -n (for all i, n ∈ N),
we call the function m a modulus of convergence of the sequence (r n ) n∈N . 2. A sequence (r n ) n∈N of real numbers converges computably if it has a computable modulus of convergence.

G.3 Computable functions on infinite sets

A function f :⊆ R n → R is computable if there is an oracle Turing machine that, given any k ∈ N, may ask for arbitrarily good rational approximations of the input x ∈ dom(f );13 after finitely many steps, it writes a rational number q on the output tape with |f (x) -q| < 2 -k [BHW08, Definition 4.1.]. However, for the purposes of computability, instead of using oracle Turing machines, it suffices to consider a simpler model, namely Turing machines with infinite I/O, i.e., Turing machines equipped with one infinite input and one infinite output tape that are read only and one-way write only, respectively. This computation model allows for a natural notion of computability on infinite strings.

Definition 26 (Computable stream functions [BHW08, Definition 4.7]).

Let Σ be an alphabet. A function F :⊆ Σ ω → Σ ω is computable if there exists a Turing machine that, given a stream p ∈ dom(F ) on the input tape, writes the output stream F (p) on the one-way output tape, and that, given a stream p ∈ Σ ω dom(F ) does not write infinitely many symbols on the output tape (and instead moves to the halting state at some point).

If X is a space, a (computable) representation of X is a (computable) map δ X :⊆ Σ ω → X.

A particular representation that we shall invariably use for R with the usual (Euclidean metric) is ρ, given as follows: Let, for each n

∈ N, Q n = {m2 --n : m ∈ Z}, and let Q D = n Q n be the set of dyadic numbers. Let ν D : N → Q D be any computable bijective numbering such that if •, •, • : N 3 → N is a computable pairing function, we have ν D ( i, j, k ) = (i-j)2 -k . The representation ρ : {0, 1} * → R is now defined by dom(ρ) = {p ∈ {0, 1} * : ∀k.(ν D (p(k) ∈ Q k ∧|ν D (p(k))-ν D (p(k +1)| < -(k+1) )))} and ρ(p) = lim n ν D (p(n)
Note that ρ is merely a representation that uses rapidly converging Cauchy sequences and that the choice of {0, 1} * is immaterial: we could have chosen a larger alphabet, or even a unary alphabet if need be.

The above definitions carry over with ease to computation between arbitrary sets w.r.t. suitable representations.

Definition 27 (Computable function [BHW08, Definition 4.7]). Let Σ be an alphabet, and let

f :⊆ X → Y be a function; moreover, let δ X :⊆ Σ ω → X and δ Y :⊆ Σ ω → Y be representations. A function F :⊆ Σ ω → Σ ω is called a (δ X , δ Y )-realizer of f if δ Y F (p) = f δ X (p) holds for all p ∈ dom(f δ X ), i.e., such that for any δ X -name of some x ∈ dom(f ), the value F (p) is a δ Y -name of f (x). The function f is (δ X , δ Y )-computable if a computable (δ X , δ Y )-realizer of f exists.

G.4 Computable metric spaces

Definition 28 (Computable metric space [BHW08, Definition 7

.1]). A triple B = (X, d, α) is called a computable metric space, if 1. d : X × X → R is a metric on X, 2. α : N → X is a sequence such that the set {α(n) | n ∈ N} is dense in X, 3. d • (α × α) : N 2 → R is a computable map (with R equipped with the repre- sentation ρ).
Observe in the above definition that countable metric spaces are separable spaces and they come with a choice of a dense subset via an enumeration.

H Computability of duality pairing

Weighted C 0 -spaces C 0 (S, W ) and their duals L 1 (S, W ) with recursively enumerable S will be equipped with a recursive enumeration α of all finite rational linear combinations of indicator functions 1 x of states x ∈ S (observe that the set of such finite linear combinations is countable and can obviously be enumerated by an algorithm). As the indicator functions themselves are computable, and the map W → Q + is computable, it follows easily that both spaces are computable Banach spaces. Moreover, we equip the product space L 1 (S, W ) × C 0 (S, W ) with the sum of the norms of the underlying spaces, i.e., (π, f

) = π * W + f W where f W = sup x∈S |f (x)| /W (x) and π * W = x∈S W (x)|π(x)|.
Finally, we choose a recursive enumeration of all linear combinations of the set x∈S {(1 x , 0), (0, 1 x )} to make the product space a computable Banach space.

Lemma 3 (Computability of duality pairing).

Let S be a recursively enumerable set, and let W : S → Q + be a computable function. The dual pairing , : L 1 (S, W ) × C 0 (S, W ) → R, mapping (π, f ) to x∈S π(x)f (x), is computable.

Proof. For any two elements (π, f ), (π , f ) ∈ L 1 (S, W )×C 0 (S, W ) of the product space, if π is finitely supported, i.e., the set of states K = {x ∈ S | π (x) = 0} is finite, we have the following estimate for the distance of π , f and π, f in R.

Next, we describe an algorithm that for each desired precision = 2 -d as input, first, makes three successive requests for approximations (π , f ) of (π, f ) of ascending precision based on the estimate from (8) and Equation ( 9), followed by computation of a rational approximation r of π, f , before finally, giving the result r as output by writing up to d+1 digits of precision in binary representation.

1. Request (π 1 , f 1 ) as approximation of (π, f ) up to precision 1. Store F = 2 f 1 W + 1; now we have F > f W as a save over-approximation of the norm f W of f . 2. Request (π 2 , f 2 ) as approximation of (π, f ) up to precision min{ /3F , 1}). This ensures f 2 W π -π 2 * W < /3 (cf. ( 8)). Moreover, with ψ x = π(x) -π 2 (x),

x∈S |ψ x f 2 (x)| = x∈S W (x)|ψ x | |f 2 (x)|/W (x) ≤ f 2 W x∈S |W (x)ψ x | = f 2 W π -π 2 * W
< /3 [START_REF] Vincent Danos | Abstracting the differential semantics of rule-based models: Exact and automated model reduction[END_REF] 3. Compute the support of π 2 , i.e., K 2 = {x ∈ S | π 2 (x) = 0} using the representation of (π 2 , f 2 ). Compute M = max x∈K2 W -1 (x). 4. Request (π 3 , f 3 ) up to precision min{ /3M2 π2 * W , /3F , 1}. thus, the second approximation of π 2 to keep the support K 2 , and the best approximation of f 3 . 6. Return enough digits of r in the dyadic representation, d + 1 after the dyadic point.

To show that the algorithm is correct, we derive as in ( 8)

|r -π, f | = | π 2 , f 3 -π, f | ≤ x∈K |π 2 (x)f 3 (x) -π(x)f (x)| + f W π -π 2 * W
and then, as in 9, we obtain

≤ x∈K π 2 (x)f 3 (x)- π(x)(f 3 (x) -f (x)) + ψ x f 3 (x) + π 2 (x)f 3 (x) + f W π -π 2 * W = x∈K |-π(x)(f 3 (x) -f (x)) -ψ x f 3 (x)| + f W π -π 2 * W ≤ x∈K |π(x)(f 3 (x) -f (x))| + |ψ x f 3 (x)| + f W π -π 2 * W ≤ x∈K |π(x)(f 3 (x) -f (x))| + |ψ x f 3 (x)| + f W π -π 2 * W
and then as in ( 10) and ( 11)

≤ /3M2 π2 * W M 2 π 2 * W + f 3 W π -π 2 * W + f W π -π 2 * W ≤ /3 + f 3 W /3F + f W /3F < .

I Counterexample [Spi12, Theorem 6.3]

Example 4 (Counterexample to Theorem 6.3 of Ref. [Spi12]). Consider the positive natural numbers N \ {0} as state space of a CTMC with q-matrix Q on S = N \ {0} given by

q xy =      x 2 y = x -1 -x 2 y = x = 1 0 y = x = 1 (12)
and consider the function W (x) = x 2 , which is positive as x > 0. Now, W is a drift function since, for x > 1, we have

QW (x) = y∈S q xy W (y) = x y=x-1
q xy W (y) = q x(x-1) W (x -1) -W (x))

= q x(x-1) (x -1) 2 -x 2 )

< 0 ≤ W (x)
and QW (1) = 0 ≤ W (1). Thus, Q is a drift function.

Next, the identity function V (x) = x on S belongs to C 0 (S, W ) and is positive. Moreover, QV (x) = q x(x-1) V (x -1) -V x)

= q x(x-1) x -1 -x)

= -q x(x-1)

= -x 2

and QV (1) = 0, whence

QV (x) = -(x) 2 • 1 S\{1} (x).
This implies that V is also a drift function and hence the minimal transition function P t with q-matrix Q is a SCSG on C 0 (S, W ), applying the first part of Theorem 6.3 of Ref. [Spi12]. Finally, we have V V = 1 < ∞ an also QV W = -W W = 1 < ∞. Thus, by the second part of Theorem 6.3 of Ref.

[Spi12], we conclude that V belongs to the domain D(Q) of the generator Q. This means by definition that the limit

|

  π , f -π, f | = | x∈S π (x)f (x) -π(x)f (x)| (7) ≤ x∈S |π (x)f (x) -π(x)f (x)| = x∈K |π (x)f (x) -π(x)f (xof f W ) ≤ x∈K |π (x)f (x) -π(x)f (x)| + f W x∈S\K |W (x)(π(x) -π (x) e x = f (x) -f (x) and φ x = π(x) -π (x), π(x)f (x) = π(x)(f (x) + e x ) = π(x)e x + π(x)f (x) = π(x)e x + (φ x + π (x))f (x)= π(x)e x + φ x f (x) + π (x)f (x).

  x∈K2 |π(x)(f (x) -f 3 (x))| ≤ f -

lim h 0 P

 0 h V -V h exists in C 0 (S, W ), thus V ∈ D(Q), and QV ∈ C 0 (S, W ). Now QV = QV in analogy to Sect. 1.4 of Ref [And91] which leads to QV ∈ C 0 (S, W ), which is false: QV = -(•) 2 • 1 S\{1} and -(•) 2 • 1 S\{1} / ∈ C 0 (S, (•) 2 )!

Definition 5 (Standard transition function [And91, p. 5f.]). A transition function on S is a family {P t } t∈R ≥0 of S × S-matrices P t

	= (p t,xy ) x,y∈S with
	non-negative, real entries p t,xy such that
	1. lim t 0 p t,xx = 1 for all x ∈ S;
	2. lim t 0 p t,xy = 0 for all x, y ∈ S such that y = x;
	3. P t+s = P t P s = ( z∈S p s,xz p t,zy ) x,y∈S for all s, t ∈ R ≥0 ; and
	4.

Table 1 .

 1 Transition functions acting on Banach spaces: state of the art

	1 is called a Yosida approximant. 6 The Yosida approximants
	are the cornerstone of the generation theorems [EN00, Corollary 3.6] that allow
	one to pass from the generator Q to the corresponding SCSG. The constant θ
	also bounds the growth of the SCSG in norm, that is to say P t ≤ M e θt for
	some M . This should already make clear that Equation (3) is crucial to obtain
	error bounds for results on the computability of SCSGs. In fact it is the starting
	point of the proof of the main result on the computability of SCSGs [WZ07,
	Theorem 5.4.2,p. 521].
	It remains to see whether we can find a Banach space to build ACPs which
	accomodate interesting (specifically unbounded) observables.

In fact, string rewriting is simply the restriction of stochastic graph transformation [HLM06] to directed, connected, acyclic, edge labelled graphs with in and out degree of all nodes bounded by one, i.e., to graphs consisting of a unique maximal path.

Even when Qf is defined, one has to check Qf = Qf , that is to say: 1 /h(P h f -f ) converges to Qf in the Banach space norm. But this will turn out to be easy compared to finding sufficient conditions for Qf to be defined.

In case f does not belong to the domain of the infinitesimal generator Q, (3) will generate a 'mild' solution that is a solution to the integral form of (acp) which might not be differentiable everywhere.

The indicator function 1x is defined as usual as 1x(y) = δxy.

Note that Anderson only considers standard transition functions [And91, p. 6].

This last condition, in presence of the first two, is equivalent to continuity of the function t → Pt, from R ≥0 to L(B), if we consider L(B) with the strong operator topology.

Thus, the oracle Turing machine may ask finitely many questions of the kind "Give me a vector p ∈ Q n of rational numbers with d(x, p) < 1/2 i ," where the exponent i may depend on the answers to the previous questions.


Emmanuel Ponsot

Patrick Susini

Sabine Meunier

Loudness Processing of Time-Varying Sounds: Recent advances in psychophysics and challenges for future research

Keywords: Loudness, Time-Varying sounds, Psychophysical methods -I-INCE Classification of Subjects Number(s): 63.1 (See

come   L'archive ouverte pluridisciplinaire

INTRODUCTION

For sound sequences of several seconds varying through time, the loudness of a stimulus cannot be considered as a static perceptual variable. As soon as listeners have access to the dynamics of their sensory representations, i.e. they are able to judge their momentary perception at different times, one needs to distinguish momentary sensations from what can be experienced and evaluated at the end of a sound, here defined as global loudness. Global loudness corresponds to the "overall loudness" perceived after the end of the sound, when the whole momentary loudness pattern (i.e., the vector made of "momentary" loudness snapshots over time) is merged into one scalar value by the human brain. This retrospective subjective merging is an operation of interest for the industry and the media who need to assess as accurately as possible how loud a sequence of several seconds or minutes will be perceived (e.g., the sound produced by the passing-by of an airplane, or an advertisement to be broadcasted on the radio).

Psychoacoustical experiments investigating dynamic loudness perception showed that global loudness does not correspond to a simple average of momentary loudness, but is rather strongly influenced by the loudest events (e.g. [START_REF] Kuwano | Continuous judgment of level-fluctuating sounds and the relationship between overall loudness and instantaneous loudness[END_REF][START_REF] Susini | Global and continuous loudness estimation of time-varying levels[END_REF][START_REF] Susini | Loudness asymmetries for tones with increasing and decreasing levels using continuous and global ratings[END_REF]. Current indicators of global loudness are based on this finding. In the media, the overall loudness of a program is taken as the integration of its momentary loudness values (predicted by simplified auditory models) that exceed a certain threshold (e.g., see the R-128 EBU-Recommendation in 4). Other indicators employed in the industry, such as N 5 or L Aeq also rely on the assumption that global loudness of a time-varying sound does not depend on the shape of its loudness pattern, i.e. that it can be predicted by its loudness distribution only. Even in the context of psychoacoustical research, Zwicker [START_REF] Fastl | Psychoacoustics: facts and models[END_REF] or Glasberg and Moore (6) suggested to use the peak of the "short-term loudness" or the "long-term loudness" time-series predicted by their models to estimate the global loudness of time-varying sounds. However, recent studies conducted even with very basic sounds have pointed out limitations to such assumption, showing for example that even with 1-kHz stimuli ramping up or down in level, global loudness is not simply based on the maximum of their short-term or long-term loudness (e.g. 7).

Why current models fail to predict global loudness of time-varying stimuli? How could this limit be overcome? The present paper attempts to provide answers to these questions. In Section 2, we discuss why, in our view, current loudness models cannot yield to accurate predictions of loudness for long time-varying sounds, as they were not designed for that purpose. In Section 3, we intend to expose the main outcomes of recent psychophysical studies on loudness of time-varying stimuli that shed light onto mechanisms not yet identified and thus not considered in current models neither. These mechanisms were exposed by the use of two broad classes of psychophysical approaches, known as "molar" and "molecular" methods [START_REF] Green | Consistency of auditory detection judgments[END_REF]. We explain why a wise combination of these approaches can be valuable for anyone interested in characterizing a whole complex psychophysical processing such as the one considered here.

CURRENT LOUDNESS MODELS WERE NOT CONCEIVED FOR GLOBAL LOUDNESS PREDICTIONS

Historically, psychoacoustics has mostly been dedicated to the study of fine timescale sensory and perceptual processing (e.g., temporal integration, backward and forward masking, detection of intensity increments/decrements), generally by conducting experiments involving very brief and simple stimuli (e.g., a few hundred milliseconds). This was the necessary step to reach a good characterization of the basic components of auditory processing, such as outer and middle ear filtering, inner ear coding, or more central stages of processing where the information incoming from both ears are merged together (e.g., binaural loudness summation). This led to quite sophisticate auditory models such as loudness models that, for the most advanced versions, are said to be suitable for temporally varying sounds (e.g., [START_REF] Fastl | Psychoacoustics: facts and models[END_REF][START_REF] Glasberg | A model of loudness applicable to time-varying sounds[END_REF]. We could schematize this approach by saying that the study of loudness has mainly been focused on what occurs in the Stage A (see Fig. 1) of the perceptual processing. The experiments were often designed so as to minimize the effects of the higher level processes potentially present in the later stages (e.g., memory processes, heuristic rules), notably through the design of experiments employing short and simple stimuli. In these cases, loudness could simply be considered as a unitary variable (i.e., the concepts of momentary loudness and global loudness were confounded).

However, in certain cases of psychophysical judgments, our perceptual system operates on a larger timescale that inescapably involves the processes present in Stage B. This is the case of global loudness processing of time-varying sounds lasting several seconds: our perceptual system needs to integrate and to store what occurred during the whole stimulus presentation, thus potentially involving "higher-level" processing stages such as memory and decisional mechanisms. We believe that this is the main reason why current loudness models cannot yield to accurate predictions of global loudness for long time-varying stimuli: they were not conceived to reveal and thus to implement the processes occurring in Stage B. In order to overcome this limitation, future auditory models that seek to predict global perceptual variables would have to incorporate higher-level processing stages identified through experiments specifically designed for that purpose. In the next section, we present several recent studies on loudness that attempted to identify and characterize these latter processes by using two classes of psychophysical approaches. The stimuli employed in these studies were 1-kHz tones or noises varying in intensity over several seconds; they were either fluctuating around a given flat intensity profile or with rising / falling intensity profiles, the latter constituting a first-order approximation of intensity dynamics encountered in natural (musical or environmental) sounds.

RECENT EXPERIMENTAL STUDIES ON LOUDNESS OF TIME-VARYING SOUNDS: TWO CLASSES OF PSYCHOPHYSICAL APPROACHES

Generalities concerning molar and molecular approaches

The terms "molar" and "molecular" psychophysics were initially introduced by Green (8) to differentiate between two classes of psychophysical approaches. The main opposition could be formulated as follows (a tentative illustration is proposed in Figure 2).

In molar psychophysics, a limited set of stimuli is typically used at the input of the perceptual system. Each stimulus is presented many times until a significant number of responses are reached, so that performance metrics (e.g. means of responses) can be computed in each condition. In the end, only these metrics, based on an average of individual responses, are considered in the analyses. In molecular psychophysics, one main stimulus is used to create an infinity of very similar but always different stimuli, this by adding some noise on the physical dimension(s) of interest (in Fig. 2, this is illustrated for the intensity profile of a rising-intensity stimulus). Observers then evaluate these noisy stimuli in hundreds or thousands of trials. In the end, their responses are reverse-correlated with noise properties to identify the mechanisms that underpin the way they processed the stimuli. As compared to molar psychophysics, the analysis undertaken in molecular psychophysics (i.e. reverse-correlation) relies on the noise introduced in each stimulus and the observer's response on each trial.

Main results from recent loudness studies using the two classes of approaches

Molar psychophysics for characterizing processes in a wide range

Molar psychophysics is probably the most direct tool for psychophysicians. Indeed, many questions can very often be addressed experimentally using a restricted number of conditions. A dependent variable is measured precisely in each condition (e.g., using adaptive tracking procedures) and an analysis of variance is conducted on the data to see whether the dependent variable measured in these conditions differ significantly or not. The main advantages of this approach are that very different conditions can be examined and that there is virtually no restriction in the design of such experiments. If one wants to assess the loudness difference between the loudness of a 1-second stimulus and the loudness of a 10-s stimulus (an order of magnitude apart), an appropriate experimental design can easily be proposed.

Molar psychophysics was thus adopted in a couple of recent psychophysical studies investigating global loudness processing of rising and falling-intensity tones with very different characteristics. Notably, specific temporal integration and memory processes were exposed by significant stretching of stimuli amplitude envelopes: durations of the ramp stimuli were varied from 1 to 16 seconds, the slope of the stimuli and their dynamics (range of intensity variation) were also manipulated in wide ranges [START_REF] Ponsot | The role of duration in global loudness evaluation of rising and falling-intensity sounds[END_REF]. A molar approach was also used to investigate global loudness processing of rising and falling ramps at very different intensity-regions, apart from more than 20 dB (7, 10).

Molecular psychophysics for characterizing processes at a given operating point

Even though the molecular psychophysical approach exists for a long time and was widely used in the neuroscience and physiology community to characterize sensory systems, the method has only been recently (for a decade or so) transposed to the study of loudness processing. It appeared to be particularly relevant for investigating how humans integrate the different temporal portions of a time-varying stimulus, because the dimension of interest being the intensity profile, some noise can be introduced on it directly to examine temporal weighting strategies in the global loudness judgments produced by the listeners. Weighting strategies have thus been exposed with time-varying sounds using both "flat-level" profiles stimuli (e.g., fluctuating stimuli) and "non-flat" profiles stimuli (e.g., rising and falling-intensity stimuli).

When the sounds are varying around overall "flat" profiles, there are evidences that the underlying temporal weighting processes employed by listeners to judge global loudness are not uniform. The first events and, to a lesser extent, the last portions of the stimuli receive stronger weights compared to the others; "primacy" and "recency" effects have been proposed to account for these findings, respectively (e.g. 11). Most previous studies demonstrated these effects with broadband or narrowband noises lasting about 1 s; we recently showed that it was also true for noises of longer durations (with 3-s white noises fluctuating in level; see [START_REF] Ponsot | Temporal weighting of loudness: Comparison between two different psychophysical tasks[END_REF]. When the overall level profile of a time-varying sound is not flat, other mechanisms seem to be engaged. Oberfeld and Plank [START_REF] Oberfeld | The temporal weighting of loudness: Effects of the level profile[END_REF] showed that when a gradual increase of intensity was imposed on the beginning of the stimuli, the primacy effect was transferred to the first next highest-level portions. Using stimuli with rising and falling-intensity profiles, Ponsot et al. [START_REF] Ponsot | Temporal loudness weights for sounds with increasing and decreasing intensity profiles[END_REF] found that "primacy" and "recency" effects disappeared and that the loudest portions of the stimuli then dominated global loudness judgments; a phenomenon referred to as "level-dominance" [START_REF] Lutfi | Molecular analysis of the effect of relative tone level on multitone pattern discrimination[END_REF]. Second-order analyses of these data were recently conducted, showing nonlinear components in the integration process of the loudest portion that had not been considered so far (see 16, Appendix F).

How combining molar and molecular methods could be helpful to characterize the whole system

On the one hand, molar psychophysics appears to be an indispensable tool for quantifying how the perceptual system responds to stimuli whose characteristics vary in a wide range. Temporal integration and memory mechanisms involved in the case of rising and falling-intensity sounds would have never been exposed without the use of stimuli of very different durations. On the other hand, level-dominance effects and the shape of the temporal integration window used by the listeners could directly be accessed to by the molecular psychophysical approach and nonlinear components of global loudness processing were identified through a molecular analysis too.

CONCLUSION AND PERSPECTIVES

In this paper, we present several recent psychophysical studies that focused on the identification of the mechanisms underlying global loudness processing of time-varying sounds lasting several seconds, which adopted either molar or molecular approaches. We show that specific aspects of this issue can more easily be tackled with one approach or the other, depending on the range and the precision of the characterization required. When one starts to investigate a complex psychophysical issue that operate at different timescales, we believe that it would worth identifying upstream the processes that can be exposed by molar or molecular approaches because the two do not make use of the same conditions and stimuli.

All the processes reported above invalidate the assumption that the global loudness of a time-varying sound does not depend on the shape of its momentary loudness pattern; thus, current models cannot account for these effects. Further attention is deserved to refine our understanding of the whole underlying perceptual machinery involved in global loudness processing. In particular, memory mechanisms involved in the integration process were highlighted but not quantified precisely yet. The respective contributions of "primacy" and "recency" effects exposed by molecular psychophysics seem to depend on the stimuli characteristics, such as their overall duration and both the size and duration of the introduced level perturbations. Such aspects have to be further explored, as for instance the fact that "primacy" effects appear to be downsized when the duration increases or the fact that the type of judgment used by listeners to evaluate global loudness alters the strength of these edge-emphasis weightings [START_REF] Ponsot | Temporal weighting of loudness: Comparison between two different psychophysical tasks[END_REF]. Efforts should also be given to examine the role of nonlinear components in the whole integration process.

Figure 1 -

 1 Figure 1 -(A) General schematic view of human perceptual processing operating at different timescales; early processes in Stage A provide continuous representation of the sensory input while latter processes of Stage B are solicited for integrating past information into an overall judgment. (B) This processing chain is rendered explicit in the case of global loudness processing. Since, the mechanisms in Stage B can potentially influence those involved in Stage A, this figure must be taken as a schematic illustration only.

Figure 2 -

 2 Figure 2 -The two main classes of psychophysical approaches illustrated in the case of the investigation of how humans process loudness of rising intensity-stimuli.

ponsot@ircam.fr

ACKNOWLEDGEMENTS

The studies conducted by the authors which are reported in this paper were realized in the framework of the LoudNat project funded by the French National Research Agency (ANR-11-BS09-016-01).