Antoine Lejay
email: antoine.lejay@univ-lorraine.fr

A Monte Carlo estimation of the mean residence time in cells surrounded by thin layers

Keywords: semi-permeable membrane, elastic Brownian motion, snapping out Brownian motion, thin layer problem, mean-residence time, Monte Carlo methods, homogenization, first eigenvalue estimation

teaching and research institutions in France or abroad, or from public or private research centers.

Introduction

Understanding the dynamics of a diffusive process in a medium with interfaces is of great importance both for modeling and simulation. An interface could be a permeable or a semi-permeable barrier. It arises as a limit of a sharp change in the properties of the underlying material. In geophysics, the inclusion of rocks, fissures, ..., leads to such interfaces, just to cite one among the many possible domains of applications (see, e.g., [START_REF] Bossy | Probabilistic interpretation and random walk on spheres algorithms for the Poisson-Boltzmann equation in molecular dynamics[END_REF][START_REF] Lejay | Simulating diffusion processes in discontinuous media: a numerical scheme with constant time steps[END_REF]).

We consider here a medium with a "thin layer" of low diffusivity seen in the limit as a semipermeable barrier the particle has difficulty to pass through (see [START_REF] Sánchez-Palencia | Nonhomogeneous media and vibration theory[END_REF]Chap. 13] or [START_REF] Perrussel | Asymptotic expansion of steady-state potential in a high contrast medium with a thin resistive layer[END_REF] for a proof relying on homogenization techniques, and [10, Th. 2] for a probabilistic proof in one dimension). The membrane surrounding living cells is an example of such a layer. A brain imaging technique such as diffusion Magnetic Resonance Imaging (dMRI) records the mean square displacement of particles of water [START_REF] Kärger | Principles and applications of self-diffusion measurements by nuclear magnetic resonance[END_REF]. Reconstructing the brain activity is done by estimating the diffusion coefficient through an inverse problem. Understanding the diffusive behavior of particles in living tissue, where thin layers are frequently present, is then essential for applying Monte Carlo methods [START_REF] Nilsson | The role of tissue microstructure and water exchange in biophysical modelling of diffusion in white matter[END_REF].

In brain imaging, the mean residence time (MRT) in the cells is an important macroscopic parameter related to the rate of convergence toward equilibrium. It could serve as an input for simplified, homogenized models such as double porosity models [START_REF] Coatléven | A macroscopic model including membrane exchange for diffusion MRI[END_REF][START_REF] Fieremans | Monte Carlo study of a twocompartment exchange model of diffusion[END_REF][START_REF] Nguyen | Parameter estimation using macroscopic diffusion MRI signal models[END_REF][START_REF] Regan | Mean residence time of molecules diffusing in a cell bounded by a semi-permeable membrane: Monte Carlo simulations and an expression relating membrane transition probability to permeability[END_REF]. The MRT is also commonly used in pharmacokinetics and geophysics [START_REF] Mcguire | A review and evaluation of catchment transit time modeling[END_REF][START_REF] Pontrelli | Mass diffusion through two-layer porous media: an application to the drug-eluting stent[END_REF][START_REF] Veng-Pedersen | Mean time parameters in pharmacokinetics. Definition Computation and Clinical Implications (Part I)[END_REF].

We propose a way to estimate the MRT in media containing thin layers with Monte Carlo simulations using continuous processes instead of random walks as in [START_REF] Regan | Mean residence time of molecules diffusing in a cell bounded by a semi-permeable membrane: Monte Carlo simulations and an expression relating membrane transition probability to permeability[END_REF]. We simulate paths of a stochastic process the density of which solves the PDE giving the concentration of the species of interest. To avoid the large computational cost induced by the layers, which forces to use a small time step around it, we replace the process by a simpler one called the snapping out Brownian motion [START_REF] Lejay | The snapping out Brownian motion[END_REF]. This is equivalent to replacing the thin layer by a semi-permeable membrane. We then show how to estimate the MRT from the Monte Carlo simulation. By the ways, we improve the estimation procedure of a first eigenvalue problem proposed in [START_REF] Lejay | Computing the principal eigenvalue of the Laplace operator by a stochastic method[END_REF]. Finally, we provide a one-dimensional benchmark test which shows the effectiveness of our approach unless the permeability of the membrane is too low.

A Monte Carlo estimation of the mean residence time

The concentration C(t, x) of the diffusive species (water, ...) follows an equation of the type

∂ t C(t, x) = ∇ • (D(x)∇C(t, x)), t > 0, x ∈ Ω (1)
in a box Ω of dimension d = 1, 2, 3 with a scalar diffusivity and periodic boundary conditions. Equivalently, we could solve (1) on a periodic medium composed of property translated copies of Ω. We consider here a simplified situation where Ω contains an interior part Ω intra separated from an exterior part Ω extra by a layer Ω layer of constant width ε > 0. The surface separating Ω layer and Ω intra (resp. Ω extra ) is denoted by Γ -(resp. Γ + ). We assume that Γ ± are of class C 2 (See Figure 1). For D 0 , µ > 0, the diffusivity is

D(x) = D 0 for x ∈ Ω intra ∪ Ω extra , while D(x) = D layer = µε for x ∈ Ω layer .
Under the above conditions, (1) is well posed and has a unique (up to an additive constant), periodic solution for any initial condition which is bounded and measurable.

The transmission condition

Since D is constant except on the interfaces Γ -and Γ + , we recast (1) into a transmission problem.

We denote by n -(x) (resp. n + (x)) the normal derivative at a point x ∈ Γ -(resp. x ∈ Γ + ) which is directed toward Ω intra (resp. Ω extra ). The concentration C is of class C 1 on the left and right of Γ ± . Besides, it satisfies

∂ t C(t, x) = D(x)△C(t, x), t > 0, x ̸ ∈ Γ ± , (2)
with

[n ± (x) • D(x)∇C(t, x)] Γ ± = 0 and [C(t, x)] Γ ± = 0, (3)
where at a point x ∈ Γ ± ,

[ f (x)] Γ ± = lim δ →0 f (x + δ n ± (x)) -f (x -δ n ± (x)).
The condition (3) is called the transmission condition (see [START_REF] Ladyženskaja | Linear and quasilinear equations of parabolic type[END_REF]§III.13,p. 224]). It specifies the continuity of the flow n ± • D(x)∇C(t, x) across the interface Γ ± as well as the continuity of the concentration. The variational formulation of ( 2)-( 3) is given in Section A.2.

The associated stochastic process

For the diffusivity D, Ω, Γ ± , . . . above extended to R d by periodicity, let us consider the stochastic process X associated to the divergence form operator L = ∇ • (D∇•). This process is a strong Markov process with continuous paths as proved in [START_REF] Stroock | Diffusion semigroups corresponding to uniformly elliptic divergence form operators[END_REF]. It is actually constructed in an abstract way from its transition semigroup whose densities are bounded from above and from below by a Gaussian kernel. When D is C 1 , the process X is solution to the Stochastic Differential Equation X t = x + t 0 σ (X s ) dB s + t 0 ∇D(X s ) ds, where σ σ T (x) = 2D(x). In our case, ∇D is not defined on Γ ± . Therefore, no such representation holds for X.

The probability p cross that the particle moves at distance at least R during a time step ∆t decreases exponentially fast with R 2 /∆t (this is a consequence of the Gaussian estimates on the density transition function of X [START_REF] Stroock | Diffusion semigroups corresponding to uniformly elliptic divergence form operators[END_REF]). For numerical simulations, we substitute to {X t+s } s≥0 a process {Z t+s } s≥0 with the same distribution until it leaves the ball of radius R centered on X t . In practice, R = 4 √ 2D 0 ∆t is sufficient, as the probability p cross is of order 6 × 10 -5 as discussed in [START_REF] Lejay | Simulating diffusion processes in discontinuous media: a numerical scheme with constant time steps[END_REF]. For away from the interfaces, the particle could be moved using Gaussian steps, which consists only in neglecting the probability to cross them.

This explains why, up to negligible errors, we focus on the dynamics of the process only at the vicinity of the layers.

Let φ be a C 2 diffeomorphism of R d . Then Y t = φ (X t ) is also a strong Markov process with continuous paths. Its infinitesimal generator A has domain Dom(A ) = φ (Dom(L )). From standard computations, for x ̸ ∈ Γ ± and f ∈ Dom(L ),

A f (x) = L ( f ∘ φ )(φ -1 (x)) = 1 2 D(φ -1 (x)) ∑ i, j,k=1,...,d ∂ 2 f ∂ x j ∂ x k (x) ∂ φ k ∂ x i (x) ∂ φ j ∂ x i (x) + D(φ -1 (x)) ∑ i, j=1,...,d ∂ f ∂ x j (x) ∂ 2 φ j ∂ x 2 i (x).
When φ is a translation or a rotation,

A f (x) = L f (φ (x)).
Recall that the layer Ω layer is surrounded by two C 2 -interfaces Γ -and Γ + at distance ε. We use an argument similar to the one of [START_REF] Bossy | Probabilistic interpretation and random walk on spheres algorithms for the Poisson-Boltzmann equation in molecular dynamics[END_REF] to reduce the problem to a "nearly 1-dimensional" one. Proposition 1. At the vicinity of any point x in Γ -, there exists a C 2 diffeomorphism φ such that φ (Γ -) and φ (Γ + ) are the hyperplanes orthogonal to (1, 0, . . . , 0) passing through φ ( x) = 0 and (ε, 0, . . . , 0).

Proof. Up to compositions with a rotation and a translation, we may assume that our reference point is

x = 0 and Γ -= U ∩ {(ω(x ′ ), x ′ ) | x ′ ∈ R d-1
} for a smooth map ω : R d-1 → R and a neighborhood U of 0. Since Γ -is orientable, we choose this map so that points in Ω intra have locally coordinates (z, x ′ ) which satisfy z ≤ ω(x ′ ), x ′ ∈ R d-1 . Let d : R d → R + be the distance between a point and Γ -. At the vicinity of Γ -, the distance is a well defined function with the same regularity as Γ - [START_REF] Gilbarg | Elliptic partial differential equations of second order[END_REF]Sect. 14.6,p. 354]. In addition, there exists a unique point y in

Γ -such that d(x) = |x -y|. Let n(x) = 1 -∇ω(x ′ )
be the unit vector normal to Γ -oriented toward the exterior of the cell at the

point x = (ω(x ′ ), x ′ ). Write φ (d, x ′ ) = x + dn(x) with x = (ω(x ′ ), x ′ ) for d ∈ R and x ′ ∈ R d-1 . Thus, φ is locally a diffeomorphism and φ -1 maps Γ -∩V to {(0, x ′ ) | x ′ ∈ R d-1 } ∩U as well as Γ + ∩V to {(ε, x ′ ) | x ′ ∈ R d-1
} ∩U for two neighborhoods V and U of 0. This function φ is the diffeomorphism we are looking at.

Using this, we focus on the component which is normal to the interface Γ ± . We are then is the same situation as the one described in [START_REF] Lejay | Simulation of a stochastic process in a discontinuous layered medium[END_REF].

Corollary 1. Using the diffeomorphism φ given above at the vicinity of a point x of the layer, we set X

φ = φ (X). The first component (X φ ) 1 of X φ is a diffusion process. Write d(x) = D(φ -1 (x, 0, . . . , 0)) for x ∈ R. The infinitesimal generator of (X φ ) 1 is P = ∇(d∇•).

The thin layer approximation

The particles have difficulties in traversing the thin layer which is then seen as a semi-permeable barrier. Thanks to a homogenization result (see e.g. [START_REF] Sánchez-Palencia | Nonhomogeneous media and vibration theory[END_REF]Chap. 13] or [START_REF] Perrussel | Asymptotic expansion of steady-state potential in a high contrast medium with a thin resistive layer[END_REF]), our diffusion equation could be replaced by a simpler one where the layer enclosed between Γ + and Γ -is reduced to a single interface Γ to which is associated a permeability µ (in m s -1 ). The condition at this single interface Γ is

µ (C(t, x+) -C(t, x-)) = D 0 n + (x) • ∇C(t, x) and [n + (x) • ∇C(t, x)] Γ = 0, x ∈ Γ, (4)
where C is the solution to the PDE

∂ t C(t, x) = D 0 △C(t, x) for x ̸ ∈ Γ and t > 0. (5)
This approximation then simplifies theoretical or numerical computations in modeling [START_REF] Coatléven | A macroscopic model including membrane exchange for diffusion MRI[END_REF][START_REF] Fieremans | Monte Carlo study of a twocompartment exchange model of diffusion[END_REF][START_REF] Perrussel | Asymptotic expansion of steady-state potential in a high contrast medium with a thin resistive layer[END_REF][START_REF] Ureña | Tissue-blood exchange of extravascular longitudinal magnetization with account of intracompartmental diffusion[END_REF]. The variational formulation of ( 4)-( 5), which is an instance of the Kedem-Katachalsky equation [19, Sect. 7.2], is given in Section A. [START_REF] Coatléven | A macroscopic model including membrane exchange for diffusion MRI[END_REF].

In dimension 1, a similar analysis can be carried out at the level of the stochastic process. Therefore, X may be replaced by a simpler process Y , called a SNapping Out Diffusion (SNOD), living in a medium where the layers are replaced by semi-permeable barriers of vanishing width which we call membranes [START_REF] Lejay | The snapping out Brownian motion[END_REF]. The coefficient µ in ( 4) is called the membrane permeability.

The state space of such a process is then a disjoint union of closed intervals, where each finite endpoint corresponds to one side of a membrane. The density transition function p(t, x, y) of X is solution to (4)-( 5) in the variables (t, x).

When a SNOD reaches a membrane, a local time clock drives a switching in which the particle is moved randomly on one or the other side of the membrane.

Using the Lamperti's type transform presented below in Proposition 2 to reduce its diffusivity to 1/2 (the one of the Brownian motion), a SNOD living in R with one membrane at 0 becomes a SNapping Out Brownian motion (SNOB) [START_REF] Bobrowski | Families of Operators Describing Diffusion Through Permeable Membranes[END_REF][START_REF] Lejay | The snapping out Brownian motion[END_REF], a process which depends on a single parameter κ > 0 related to D 0 and µ.

A SNOB behaves like an elastic Brownian motion on (-∞, 0-] ∪ [0+, ∞). When its local time at 0± goes above an independent exponential random variable of mean 1/κ, the process is displaced with probability 1/2 to either at 0+ or 0-. It starts afresh with the same dynamics. The behavior at 0 models a semi-permeable barrier the particle goes through after a random time in a Markovian manner.

The infinitesimal generator of the SNOB is S = 1 2 △. The domain of S contains as a core the functions satisfying the transmission conditions ∇u(0+) = ∇u(0-) and κ(u(0+) -u(0-)) = 2∇u(0) and which are bounded and of class C 2 elsewhere with a bounded first and second order derivatives.

Proposition 2 ([9, Sect. 6]). Let Y be a SNOD in a one-dimensional media of constant diffusivity D 0 with one membrane at 0 of permeability µ. With ψ(x) = y/ √ 2D 0 , the process Y ψ = ψ(Y ) is a SNOB of parameter κ = 2 3/2 µ/ √ D 0 when restricted to the vicinity of 0.

When f is a bounded, continuous function and Z a SNOB of parameter κ, it was shown in [START_REF] Lejay | The snapping out Brownian motion[END_REF] that

E x [ f (Z t )] = E x 1 + e -κL 0 t (B) 2 f (|B t |) + E x 1 -e -κL 0 t (B) 2 f (-|B t |) for x ≥ 0, ( 6
)
where B is a Brownian motion and L 0 t (B) its symmetric local time at 0. Using the representation (6) of the semigroup and some explicit formula for the simulation of the local time [START_REF] Lépingle | Euler scheme for reflected stochastic differential equations[END_REF], a simulation algorithm could be given for the SNOB which is presented in [START_REF] Lejay | Estimation of the mean residence time in cells surrounded by semi-permeable membranes by a Monte Carlo method[END_REF][START_REF] Lejay | The snapping out Brownian motion[END_REF].

The tangential components behave like a Brownian motion and are then easily simulated.

A Monte Carlo estimator of the mean residence time

The mean residence time (MRT) is an important macroscopic quantity that could be used to construct two-scale models such as in [START_REF] Coatléven | A macroscopic model including membrane exchange for diffusion MRI[END_REF].

In the box Ω, the operator Ľ = ∇(D∇•) with domain Dom( Ľ

) = { f ∈ Ľ2 (Ω) | Ľ f ∈ Ľ2 (Ω)
} is self-adjoint with a compact resolvent (here, Ľ2 (Ω) is the space of square integrable Ω-periodic functions, see Section A.2). In particular, it admits a spectral decomposition 0 ≤ λ 0 < λ 1 < • • • with eigenfunctions η 0 , η 1 , . . . chosen orthonormal for the Ľ2 (Ω) scalar product (regarding the sign of the eigenvalues, we use the convention L η k = -λ k η k ). The first eigenvalue is λ 0 = 0 and its associated eigenfunction is η 0 (x) = 1/ V per , where V per is the volume of the box Ω.

Proposition 3 ([16, Eq. ( 6)] or [9, § 2.1]). The mean residence time τ in the cell satisfies the approximation 1

τ ≈ V intra V per λ 1 , (7)
where V intra is the volume of interior part of the cell.

Determining the first positive eigenvalue λ 1 is then crucial for estimating the MRT. We explain now its relationship with the process X defined as the projection of the original process X on the periodic box Ω.

The density transition function p with respect to the Lebesgue measure of X admits the spectral expansion

p(t, x, y) = +∞ ∑ k=0 e -λ k t η k (x)η k (y).
Since η 0 (x) is constant and for a fixed point x ∈ Ω,

P(t) = Ω intra p(t, x, y) dy = V intra V per + Φ(x)e -λ 1 t + o(exp(-λ 1 t)) (8)
with Φ(x) = ∑ ks.t.λ k =λ 1 η k (x) Ω intra η k (y) dy, since the multiplicity of λ 1 is not necessarily 1.

As t increases, P(t) converges to the ratio P ∞ = V intra /V per . Thus, λ 1 is actually the exponential rate of convergence toward equilibrium.

On the other hand, P(t

) = P x [ Xt ∈ Ω intra ]
. By simulating a large number N of paths of the stochastic process X, we estimate P(t) at some given times {t k } k=1,...,M using the empirical estimators

P N (t k ) = 1 N N ∑ i=1 1 X(i) t k ∈Ω intra , (9)
where X(1) , . . . , X(N) are N independent samples of the paths X.

Remark 1. For practical computation, we use the SNOD approximation Y of X so that any thin interface is reduced to a single membrane. A similar analysis, carried out on its infinitesimal generator and eigensystem, shows that P x [ Yt ∈ Ω intra ] ≈ P(t) for t ≥ 0.

From the knowledge of {P N (t k )} k=1,...,M , we estimate λ 1 through a linear regression on {log(P N (t k )-P ∞ )} k=1,...,M . However, the first and last times t 1 and t M shall be carefully chosen. For this, we use the procedure presented in Section 3.

3 A procedure to estimate some exponential rate from Monte Carlo measurement

Let P(t) be the probability of an event, which decreases to P ∞ with the time t at some exponential rate λ . This is typical for first eigenvalues problems, such as [START_REF] Ladyženskaja | Linear and quasilinear equations of parabolic type[END_REF] and the estimation of survival probability of a Brownian motion killed when it exits from a domain [START_REF] Lejay | Computing the principal eigenvalue of the Laplace operator by a stochastic method[END_REF].

More precisely, we are willing to estimate λ under the assumption that

P(t) = P ∞ + c exp(-λt) + r(t) with r(t) = o t→∞ (exp(-λt)), (10)
where P ∞ ≥ 0 is known, yet c > 0, λ > 0 and r are unknown. In particular, we ignore how fast r(t) exp(-λt) decreases to 0. At different times 0 < t 1 < • • • < t M , this quantity P(t i ) is estimated by P N (t i ) from a Monte Carlo estimation with N independent samples. A Gaussian confidence interval allows one to assert that with a probability α,

|P(t) -P N (t)| ≤ c α P(t)(1 -P(t)) √ N
for N large enough and t > 0, with, e.g., c α = 2.57 for α = 99 %. From [START_REF] Lejay | The snapping out Brownian motion[END_REF], log(P N (t) -P ∞ ) = log cλt + ξ (t), with ξ (t) = log 1 + e λt r(t) + e λt (P N (t) -P(t)) ≈ e λt r(t) + e λt (P N (t) -P(t)).

The difficulty is then to find a time interval [t -,t + ] so that ξ (t) is small in front of log c and λt, so that a linear regression could be used.

When t -and t + are fixed, c est and λ est are the values which minimize

SQ(t -,t + ; c, λ ) = ∑ i=1,...,M t -≤t i ≤t + (log c -λt i -log(P N (t i ) -P ∞ )) 2 = ∑ i=1,...,M t -≤t i ≤t + ξ (t i ) 2 .
The residual sum of squares is

RSS(t -,t + ) = SQ(t -,t + ; c est , λ est ) = min c,λ SQ(t -,t + ; c, λ ).
Several criteria have been proposed in [START_REF] Lejay | Computing the principal eigenvalue of the Laplace operator by a stochastic method[END_REF] to find the best time interval [t -,t + ] for estimating λ through a least squares procedure. Yet they partially rely on a graphical procedure, which we want to avoid.

The Monte Carlo error can be neglected unless c α P(t)(1 -P(t))/ √ N, the magnitude of P N (t) -P(t), has the same order as P(t) -P ∞ ≈ c exp(-λt). Hence, it is useless to consider times t at which

c α P N (t)(1 -P N (t)) √ N ≥ β (P N (t) -P ∞ ) with say β = 1 10 .
This criterion fixes a limit for the maximum sampling time t M to consider. After this time, the estimations are no more relevant as being "flooded" by the Monte Carlo error. Getting rid of r(t) is more cumbersome, as this function is generally unknown. In the report [START_REF] Lejay | Estimation of the mean residence time in cells surrounded by semi-permeable membranes by a Monte Carlo method[END_REF], we have proposed a new methodology in which t -and t + are determined adaptively using the Akaike Information Criteria (AIC) [18,Example 13.6,p. 179]. This approach is suitable as the number of values of P N (t i ) to use varies with t -and t + .

Since there are two parameters (c and λ ) to estimate, we use the customary simplified formula for the AIC AIC(t -,t + ) = n log(RSS(t -,t + )) -n log n + 4, with n = #{t -≤ t i ≤ t + ; i = 1, . . . , M}.

The procedure is: We draw randomly some pairs (t -,t + ) with 0 ≤ t -< t + ≤ t M . We select the one for which AIC(t -,t + ) is minimized. 0 2,000 4,000 6,000 8,000 0.9 show that the values which are close to the minimum of AIC provide pretty close estimated values of λ and use a lot of sample times. We have used N = 100 000 paths and t i+1 -t i = 20 µs with the parameters of Table 1. This procedure, whose result is illustrated in Figure 2, no longer relies on a graphical interpretation as in [START_REF] Lejay | Computing the principal eigenvalue of the Laplace operator by a stochastic method[END_REF].

A computational benchmark test in dimension one

We test our simulation scheme by estimating a MRT, which we compare to the numerical, deterministic estimation of an eigenvalue problem through [START_REF] Ladyženskaja | Classical solvability of diffraction problems for equations of elliptic and parabolic types[END_REF]. For this, we place ourselves in a one-dimensional

Name

Quantity

From To

D 0 diffusivity in Ω intra ∪ Ω extra 2 × 10 -3 µm 2 µs -1 3 × 10 -
Ω extra = Ω ∖ (Ω intra ∪ Ω layer ). The diffusivity is D(x) = D 0 , for x ∈ Ω intra ∪ Ω extra .
In each layers of width ε, the diffusivity is D(x) = µε.

Using the thin layer approximation, Ω layer is just reduced to just two points at the positions of the membranes replacing the layers.

Still using this approximation, we consider the process Y living in the box with two membranes at y 1 and y 2 , and diffusivity D 0 . This process is obtained by projecting a process Y living in the periodic medium. The process Z = Y / √ 2D 0 has diffusivity 1/2 outside and inside the cell. Around each membrane, Y acts 1 like a SNOB with parameter κ = 2 3/2 µ/ √ D 0 (Proposition 2). We simulate Y with the algorithm given in [START_REF] Lejay | The snapping out Brownian motion[END_REF] to deal locally with each membrane. The impact of this choice is studied in Section B.

According to [START_REF] Ladyženskaja | Classical solvability of diffraction problems for equations of elliptic and parabolic types[END_REF], the MRT is then approximated by τ ≈ λ -1 V per /V intra . We use a realistic range of parameters given in Table 1. We perform N independent simulations. With (9), we compute empirical estimations P N (t k ) of P(t k ) given by ( 8) at given times 0 < t 1 < • • • < t M ≤ T . We then use the procedure of Section 3 to estimate λ .

Being the first non-zero eigenvalue of Ľ , λ could also be computed deterministically as an eigenvalue of a Sturm-Liouville problem. Numerically, it is easily checked that the eigenvalues with the thin layer approximation are close to the ones of the original problem [START_REF] Lejay | Estimation of the mean residence time in cells surrounded by semi-permeable membranes by a Monte Carlo method[END_REF]. The process is simulated only using the thin layer approximation, which leads to faster computations. Without the thin layer approximation, the simulations could be performed using the algorithms given in [START_REF] Lejay | Simulating diffusion processes in discontinuous media: a numerical scheme with constant time steps[END_REF]. Yet the thinness of the layer forces to use a very small time step, burdening the computational cost.

Numerical results are given in Figure 3. Unless µ is too small, the Monte Carlo results are in good agreements with the deterministic approach. A small value of µ means that the particles have difficulties to pass through the semi-permeable barrier. This problem is related to a rare events simulation problem. This shall be subject to future investigations.

A The thin layer problem: variational formulation A.1 Notations

We denote by Ľ2 (Ω) be the space of Ω-periodic functions which are square-integrable, that is

‖ f ‖ 2 Ľ2 (Ω) = ⟨ f , f ⟩ < +∞ with ⟨ f , g⟩ = Ω f (x)g(x) dx.
1 That is, by neglecting the excursions which are so large that Y crosses the other interface during the small time step. For an integrable Ω-periodic function f , we write

f Ω = 1 |Ω| Ω f (x) dx,
where |Ω| is the volume of Ω.

The Sobolev space Ȟ1 (Ω) is the completion of Ω-periodic functions of class C 1 (Ω) with respect to the norm ‖ f ‖ Ȟ1 (Ω) defined by

‖ f ‖ 2 Ȟ1 (Ω) = Ω | f (x)| 2 dx + Ω |∇ f (x)| 2 dx.
We then define

V = f ∈ Ȟ1 (Ω) f Ω = 0 .
As Ω is compact and owing to the periodic boundary conditions, Ȟ1 (Ω) is compactly embedded in Ľ2 (Ω), and so is V. This result is an instance of the Rellich-Kondrachev theorem.

A.2 Variational formulation for the initial diffusion problem

We consider the bilinear form

E( f , g) = Ω ⟨D(x)∇ f (x), ∇g(x)⟩ dx, ∀( f , g) ∈ Ȟ1 (Ω)
where D is the diffusivity defined in Section 2 in the medium with layer of width ε.

From classical results, the Poincaré inequality holds for each ε > 0: there exists a constant C ε such that

Ω ( f (x) -f Ω ) 2 dx ≤ C ε Ω |∇ f (x)| 2 dx, ∀ f ∈ Ȟ1 (Ω).
The proof of the next proposition proceeds from standard results since E is bilinear, symmetric and coercive thanks to the Poincaré inequality (See e.g. [START_REF] Zeidler | Nonlinear functional analysis and its applications. II/A. Linear monotone operators[END_REF]Chap. 22

])

Proposition 4. There exists a self-adjoint operator Ľ with domain Dom( Ľ

) ⊂ V such that E( f , g) = -⟨ Ľ f , g⟩, ∀( f , g) ∈ Dom( Ľ ) × V.
Besides, for any λ ≥ 0, λ -Ľ is invertible and (λ -Ľ ) -1 : Ľ2 (Ω) → Dom( Ľ ) is a bounded operator which is compact from Ľ2 (Ω) to Ľ2 (Ω). In particular, Ľ has an eigensystem {λ k , η k } k≥0 such that Ľ η k = -λ k η k , {λ k } k≥0 has no accumulation point and {η k } k≥0 is an orthonormal basis of Ľ2 (Ω).

Results on the local regularity of the solution through the interfaces Γ + and Γ -may be found e.g., in [START_REF] Ladyženskaja | Classical solvability of diffraction problems for equations of elliptic and parabolic types[END_REF].

A.3 Variational formulation for the thin layer approximation

With the thin layer approximation, we define by Γ = Γ + = Γ -, the layer of zero width which splits Ω to Ω extra and Ω intra .

We define

W = f ∈ Ȟ1 (Ω ∖ Γ) f is Ω-periodic , f Ω = 0 .
Applying the Rellich-Kondrachev inequality in both Ω extra and Ω intra proves that W is compactly embedded into Ľ2 (Ω).

Let us recall the trace theorem: the application γ

+ : f ∈ C (Ω extra ) → C (Γ) may be extended continuously from H 1 (Ω extra ) to L 1/2 (Γ), the functions g : Γ → R such that Γ g(x) 2 dσ is finite. The trace γ -on Γ -is defined similarly. Hence, for f ∈ W, [ f ] Γ = γ + f -γ -f is well defined as a function in L 1/2 (Γ). Therefore Γ [ f ] 2 Γ dσ < +∞,
where σ is the measure on the surface Γ.

Let us define the bilinear form

S( f , g) = Ω ⟨D 0 ∇ f (x), ∇g(x)⟩ dx + µ Γ [ f ] Γ • [g] Γ dσ , ∀ f , g ∈ W.
We now prove an inequality of Poincaré type.

Lemma 1. There exists a constant C such that

‖ f ‖ 2 Ľ2 (Ω) ≤ CS( f , f ), ∀ f ∈ W. ( 11
)
Proof. Assume that [START_REF] Lejay | Simulation of a stochastic process in a discontinuous layered medium[END_REF] is not true, meaning that there exists a sequence

{ f n } such that f n ∈ Ľ2 (Ω) with ‖ f n ‖ Ľ2 (Ω) = 1 and 1 ≥ C n S( f n , f n ) with C n ---→ n→∞ +∞. Applying the Poincaré inequality to both ( f n ) |Ω extra and ( f n ) |Ω intra implies that f n converges strongly to f that takes a constant value f + (resp. f -) on Ω extra (resp. Ω intra ). Yet Γ [ f n ] 2
Γ dσ converges to 0. By continuity of the traces γ + and γ -, f + = f -. Since ( f n ) Ω = 0, f Ω = 0 and thus f = 0.

From Lemma 1 and standard results (See e.g. [START_REF] Zeidler | Nonlinear functional analysis and its applications. II/A. Linear monotone operators[END_REF]Chap. 22]), we get the following proposition which is similar to Proposition 4. Proposition 5. There exists a self-adjoint operator S with domain Dom(S ) ⊂ W such that S( f , g) = -⟨S f , g⟩, ∀( f , g) ∈ Dom(S ) × W.

Besides, for any λ ≥ 0, λ -S is invertible and (λ -S ) -1 : Ľ2 (Ω) → Dom(S ) is a bounded operator which is compact from Ľ2 (Ω) to Ľ2 (Ω). In particular, S has an eigensystem {λ k , η k } k≥0 such that S η k = -λ k η k , {λ k } k≥0 has no accumulation point and {η k } k≥0 is an orthonormal basis of Ľ2 (Ω).

We consider now Ω-periodic functions f , u and v which are respectively of class C , C 1 and C 2 on Ω extra and on Ω intra (assuming now that the width of the layer is ε = 0) such that D 0 △u = f on Ω extra ∪ Ω intra .

Assume that for any such v as above, S(u, v) = -⟨ f , v⟩, that is

Ω f (x)v(x) dx = D 0 Ω ∇u(x) • ∇v(x) dx + µ Γ [u] Γ × [v] Γ dσ .
Using an integration by parts on Ω extra and Ω intra , since n + = -n -,

D 0 Γ + ∂ u ∂ n + (x)v(x) dσ (x) -D 0 Γ - ∂ u ∂ n + (x)v(x) dσ (x) = µ Γ [u] Γ × [v] Γ dσ ,
where Γ + (resp. Γ -) denotes the sides of Γ oriented toward Ω extra (resp. Ω intra ). Taking v continuous around Γ implies that ∂ u

∂ n Γ = 0. Taking v = 1 or Γ + and v = 0 on Γ -implies that D 0 ∂ u ∂ n + Γ + (x) = D 0 ∂ u ∂ n + Γ - (x) = µ[v] Γ (x), x ∈ Γ. (12)
This proves that if u is a classical solution to the variational problem S(u, v) = -⟨ f , v⟩ for any v ∈ W, then u satisfies (12) on Γ and thus S u = f . Conversely, if u satisfies ( 12) on Γ as well as D 0 △u = f , then u is solution to the variational problem S(u, v) = -⟨ f , v⟩ for any v ∈ W.

B Exact and approximate simulation

Let T be a time horizon and ∆t a time step so that ∆t = T /n for some integer n. Our estimation of the MRT relies on the simulation of the successive positions {X k∆t } k=0,...,n of the stochastic process in the medium with the membranes. We consider the particle living in the one-dimensional periodic media of Section 4: in the periodic box Ω of length 2(α + α ′ ), there are two membranes at distance 2α enclosing an inner cell Ω intra . The diffusivity is D 0 both in the inner and the outer cell.

Let B be a Brownian motion. We set β k = (B t ) k∆t≤t≤(k+1)∆t . We also define for each k ≥ 0, ξ k = {ξ k,i } i≥0 , an array of independent uniform random variables that are used to decide whether or not the particle crosses the membrane. Thanks to the Markov property, the successive positions X k∆t are computed iteratively by X (k+1)∆t = p(X k∆t + Ψ(X k∆t , β k , ξ k )), where p denotes the projection on the box Ω and Ψ is an abstract map.

Around each membrane at ±α, we define an interface layer by [±αδ , ±α + δ ] with δ = 4 √ 2D 0 ∆t. A numerical approximation {Z k } k=0,...,n of {X k } k=0,...,n is computed iteratively Z (k+1)∆t = p(Z k∆t + Ξ(Z k∆t , β k , ξ k )) using an approximation Ξ of Θ. Two rules are used to move the particle (See Figure 4). When Z k∆t does not belongs to an interface layer, we use a Gaussian step Z (k+1)∆t = Z k∆t + √ 2D 0 (B (k+1)∆t -B k∆t ). When Z k∆t belongs to an interface layer, we use a SNOB step around the interface (in which the other membrane is not "seen"). An algorithm for this step is easily deduced from (6) (see [START_REF] Lejay | Estimation of the mean residence time in cells surrounded by semi-permeable membranes by a Monte Carlo method[END_REF][START_REF] Lejay | The snapping out Brownian motion[END_REF] for the precise algorithm). A SNOB step is similar to a Gaussian step except for a possible symmetric replacement of the particle around the membrane according to the value of the local time and the one of ξ k,1 .

For t ≥ 0 , we write

ζ k (t) = X k∆t + √ 2D 0 (B t+k∆t -B k∆t ). Thus, Ξ(X k∆t , β k , ξ k ) = Θ(X k∆t , β k , ξ k ) on the
α -α -α -α = α + α α -4 √ 2D 0 ∆t α + 4 √ 2D 0 ∆t Ω intra Ω extra
Gaussian step SNOB step

Gaussian step

Figure 4: The periodic media and its interface layer (in gray).

We neglect the multiple membranes hits (I) since it implies that the range of ( √ 2D 0 B t ) t∈[0,∆t] is larger than the distance between the membranes. We then focus on (II).

On the event ∩ n-1 i=0 E i , Z n∆t = X n∆t and the simulation is exact. We now quantify the probability of such an event. If A denotes the complement of an event A,

P ∩ n-1 i=0 E i = P ∪ n-1 i=0 E i ≤ n-1 ∑ i=0 P E i . (13)
Using the Markov property, we set ζ (t) = x + √ 2D 0 B t for some point x and we consider the events: The event A(x) may then be decomposed as the (non disjoint) union of four events A 1 (x) = {ζ 0 = x ∈ (-α, α), ζ ∆t ̸ ∈ (-α, α) and α has been hit} A 2 (x) = {ζ 0 = x ∈ (-α, α), ζ ∆t ̸ ∈ (-α, α) and -α has been hit} A 3 (x) = {ζ 0 = x ∈ (-α, α), X ∆t ∈ (-α, α) and α has been hit} and A 4 (x) = {X 0 = x ∈ (-α, α), X ∆t ∈ (-α, α) and -α has been hit}.

A similar decomposition holds for A ′ (x). Using a symmetry argument, A 3 (x) = A 4 (-x) and A 1 (x) = A 2 (-x).

The probability that a Brownian bridge hits a barrier is known, so that by conditioning over the value of ζ ∆t = y, it holds that with ∆ = 2D 0 ∆t,

P(A 4 (x)) ≈ α-δ -α+δ exp -2(x -α)(y -α) ∆ exp -(y -x) 2 2∆ dy √ 2π∆ .
After some algebra, and neglecting multiple membranes hits,

P(A 4 (x)) = P(A 3 (-x)) ≈ α-δ -α+δ exp -(2α -x -y) 2 2∆ dy √ 2π∆ = Φ -α + x √ ∆ -δ -Φ -3α + x √ ∆ -δ .
Again, by neglecting multiple membranes hits,

P(A 2 (x)) = P(A 1 (-x)) ≈ Φ -α -x √ ∆ -Φ -α -2α ′ -x √ ∆ .
The probabilities P(A ′ i (x)) are estimated the same way by exchanging α ′ and α. Using [START_REF] Lejay | Simulating diffusion processes in discontinuous media: a numerical scheme with constant time steps[END_REF], by neglecting multiple membranes hits during one time step, p = P ∩ n-1 i=0 E i n-1 ∑ k=0 P(A(X k∆t )) + P A ′ (X k∆t ) .

To simplify the computations, we assume that the initial distribution of the process follows the invariant measure, which is the Lebesgue's one. Hence, each X k∆t is uniformly distributed over Ω. Thus, for i = 1, . . . , 4,

P(A i (X k∆t )) = 1 2α + 2α ′ α-δ -α+δ A i (x) dx and P A ′ i (X k∆t ) = 1 2α + 2α ′ α ′ -δ -α ′ +δ A ′ i (x) dx.
For a numerical result, we use the following parameters taken from Table 1: α = V intra /2 = 10 µm, α ′ = V extra /2 = 1 µm and D 0 = 3 × 10 -3 µm 2 µs -1 for a simulation up to T = 5 × 10 4 µs. This choice that maximizes the number of times a Gaussian step is performed. The time step is ∆t = 1 µs so that 5 × 10 4 steps are performed. Thus, p 5 × 10 4 × 2.2 × 10 -6 = 0.11. This implies that for more than 89% of the particles, the exact and the approximate schemes lead to the same sequence of positions. The actual value is greater than this because we have only estimated the number of time the particle hit one of the membranes when it starts away from the interface layer. However, this does necessarily not lead to a bad replacement. A finer analysis could then be performed using the distribution of the local time of a Brownian bridge.

Figure 1 :

 1 Figure 1: A box Ω of the periodic medium.

  -6 µm, D 0 = 2 × 10 -3 µm 2 µs -1 2 = 10 -2 µm, L intra = 10 µm, L extra = 2 µm

Figure 2 :

 2 Figure 2: Estimation of λ using a least square procedure. The curve in red represents the estimated function c exp(-λt) from the time interval [t -,t + ] delimited by the triangle. The two figures at theshow that the values which are close to the minimum of AIC provide pretty close estimated values of λ and use a lot of sample times. We have used N = 100 000 paths and t i+1 -t i = 20 µs with the parameters of Table1.

10 -6 10 - 5 10 -4 10 -DFigure 3 :

 105103 Figure 3: Estimation of λ by a deterministic method and a Monte Carlo one with N = 100 000 particles.

  following event E k = {ζ k does not hit any membrane during [0, ∆t]} ∪ {ζ k hit only the membrane at ±α during [0, ∆t] and |X k∆t ∓ α| ≤ δ }.In other words, Ξ(X k∆t , β k , ξ ) and Θ(X k∆t , β k , ξ k ) may give different values if (I) the two membranes are hit during [0, ∆t] (multiple membranes hits) or (II) X k∆t is not in an interface layer yet ζ k reaches one membrane during [0, ∆t].

  0

A

  (x) = {ζ hits -α or α during [0, ∆t] with ζ 0 = x ∈ Ω intra } A ′ (x) = {ζ hits -α or α during [0, ∆t] with ζ 0 = x ∈ Ω extra }.

Table 1 :

 1 3 µm 2 µs -1 A realistic range of parameters. medium, which allows explicit computations of the first positive eigenvalue λ = λ 1 . We consider a cell Ω intra surrounded by two thin layers at positions y 1 and y 2 in a box Ω. The set of points in the layers is Ω layer . We set

	µε ε µ = V intra V extra V per T	D layer ε	diffusivity in Ω layer layer width membrane's permeability intracellular size extracellular size medium's period max. observations' time	10 -8 µm 2 µs -1 10 -6 µm µs -1 5 µm ≈ V intra /10 10 -6 µm 2 µs -1 0.01 µm 10 -4 µm µs -1 20 µm V intra +V extra 20 000 µs 50 000 µs

Acknowledgement. The author wishes to thank J.-R. Li (Inria) who have proposed this research subject and provided the data of Table 1. This work has been supported by the ANR projects SIMUDMRI (ANR-10-COSI-0008) and H2MNO4 (ANR-12-MONU-0012-01).

The author also wishes to thank the referees for their careful reading that led to numerous improvements.


F López-Mir

V Naranjo

S Morales

J Angulo

PROBABILITY DENSITY FUNCTION OF OBJECT CONTOURS USING REGIONAL REGULARIZED STOCHASTIC WATERSHED

Keywords: mathematical morphology, segmentation, stochastic watershed, probability density function

In this paper, a probability density function of object contours based on the stochastic watershed transform is carried out. The watershed transform produces an over-segmentation of the image due to noise, illumination problems, low contrast, etc., because each regional minimum of the image gives place to a region in the output image. To solve this problem, the efforts are focused on the definition of markers to impose new minima in the image, and enhancing the gradient image. The stochastic watershed performs a probability density function (pdf) of the object contours based on a MonteCarlo simulation of random markers. A variation on the method for defining this pdf based on regional regularization of the image is carried out. The objective is to obtain a pdf of the object contours with less noise and better contrast than that produced by the stochastic watershed to use it as a new gradient image for segmentation purposes.

INTRODUCTION

The watershed transform is an image segmentation technique based on mathematical morphology with application in multiple situations [START_REF] Meyer | The watershed concept and its use in segmentation : a brief history[END_REF]. This technique interprets the grayscale image as a topographic surface where the pixel intensity is equal to the altitude, Fig. 1. Each regional minimum of this surface represents a catchment basin, and the idea of the watershed transformation is to flood the surface from these basins. To prevent the merging of two different basins the watershed lines or watersheds are built. These watersheds lines will conform the result of the watershed transform. The standard framework of watershed transform is based on the flooding algorithm, but there are other alternative frameworks: based on a continuous formulation using topographic distance [START_REF] Najman | Watershed of a continuous function[END_REF]; the topological watershed based on discrete geometry tools [START_REF] Bertrand | On topological watersheds[END_REF]; the graph-based watershed using minimum spanning-tree algorithms [START_REF] Meyer | Watersheds, waterfalls, on edge or node weighted graphs[END_REF]; the power watershed algorithm [START_REF] Couprie | Power watershed: A unifying graph-based optimization framework[END_REF]; the viscous watershed [START_REF] Vachier | The viscous watershed transform[END_REF]; etc. A segmentation based on the watershed transform is successfully achieved if the minima of the image are located inside of the objects of interest and the edges of these objects are the maxima of the image. For this reason, the gradient image is usually used as input of this transform. The oversegmentation is produced by the great number of minima that appear in the gradient image. The reduction of these minima is the objective of several solutions to improve the watershed transform. One approach is the marker-watershed transform where the minima are imposed by a set of markers introduced artificially in the image gradient [START_REF] Meyer | Morphological segmentation[END_REF]. The automatic definition of these markers is not an easy task in some cases but manual markers are user-dependent and the robustness of the watershed transform is limited. The hierarchical watershed paradigms and the waterfall algorithm are based on the strategy of region fusion and the relative weights of the watershed lines that separate the initial basins [START_REF] Beucher | Watershed, hierarchical segmentation and waterfall algorithm[END_REF]. Other solutions try to achieve an equivalent image of the object edges that reduces the minima inside of the objects and enhances the boundaries. This is the case of the stochastic watershed [START_REF] Angulo | Stochastic watershed segmentation[END_REF]. Original stochastic watershed is based on Monte Carlo simulations. It was shown in [START_REF] Meyer | A stochastic evaluation of the contour strength[END_REF] that the corresponding probability density function of contours obtained by stochastic watershed can be calculated using graph algorithms; this latter approach is out of the scope of the paper. More recently, other variants of stochastic watershed has been also proposed [START_REF] Bernander | Improving the stochastic watershed[END_REF].

The objective of this paper is to obtain a probability density function of the object edges in a similar way of the original stochastic watershed but taking into account some ideas of the hierarchical and waterfall paradigms in order to achieve an equivalent gradient image less noisy where the regional minima are considerable reduced in comparison to the stochastic watershed or the morphological gradient image. This image will be the input of the watershed algorithm or other variation as the marker-controlled watershed.

The paper is divided as follows. In Section 2, some tools used in this work, the methodology and the mathematical formulation of our approach is explained. Results of several experiments that demonstrate the advantages of this method are presented in Section 3. Finally, in Section 4 the main conclusions are discussed.

METHOD

Tools

Let I be an image defined as the mapping

I(x) : E → T , (1)
where x ∈ E is the pixel position in the support space of pixels E, e.g., for 2D images E ⊂ N 2 . In the case of discrete images, T = {t min , t min + 1, ..., t max } represents the pixel intensity.

The main operators define in mathematical morphology are the dilation (δ) and the erosion (ε) and the combination of them produces the opening (γ) and the closing (ϕ), the two basic filters in mathematical morphology [START_REF] Soille | Morphological Image Analysis[END_REF]. The morphological gradient is defined as ρ(I) = δ(I) -ε(I). Typically, it is the input image of the watershed algorithm and its variations. With the marker-controlled watershed the minima are imposed by a set of seeds that mark the objects of interest and the background, so the output regions can be controlled according the number of makers defined. The marker-controlled watershed transform is defined as Eq. 2:

W S(ρ(I)) m(x) = W S(R ε (ρ(I)∧m(x)) (m(x))), (2)
where m(x) is a set of markers:

m(x) = 0, if x belongs to a marker, 255, otherwise, (3)
and R ε (ρ∧m(x)) (m(x)) is a reconstruction by erosion that imposes the minima of the gradient image, ρ(I), to the pixels set to zero in m(x). Thus, the output regions are limited to the number of connected components of marker image m(x).

In marker-controlled watershed, marker definition is essential to achieve an accurate segmentation. The stochastic watershed uses a different strategy due to the random definition of the markers. This variation is based on an interactive process whose objective is obtain a probability density function of the object contours.The stochastic watershed transform is based on applying M times the marker-controlled watershed with N markers (m i (x)) that change randomly in each iteration. Using this approach a probability density function can be computed with the Parzen method [START_REF] Duda | Pattern Classification and Scene Analysis[END_REF] as Eq. 4.

pdf(x) = 1 M M i=0 W S i (x) * G(x; s), (4)
where W S i (x) = W S(ρ(I)) mi(x) and G(x; s) represents typically a Gaussian kernel of variance s 2 (in our case, s = 3) and mean m ( m = 0), that is defined as:

G(x; s) = 1 2πs 2 e -( x 2 2s 2 ) . (5)
Fig. 2 shows some images that produce the algorithm. Fig. 2a and 2b are a synthetic image and its respective gradient. This image has internal low gradients that divided the cube in four different regions. This situation tries to simulate when the illumination or contrast is not uniform in an image. Fig. 2c and2d show two outputs of the marker-controlled watershed with two different set of markers created randomly. In these output images the watershed lines has the same weight and low gradients are enhanced. Fig. 2e is the average of the 30 output images (the watershed lines) of the markercontrolled watersheds carried out following Eq. 4 to obtain the pdf of object contours. In this image, the low gradients that separate the square in four regions are high contrasted while artificial edges that are created in each iteration are reduced by the stochastic procedure. Fig. 2f is the average of Fig. 2b and2e proposed in [START_REF] Angulo | Stochastic watershed segmentation[END_REF] ant it is observable that local gradients intensity are not taken into account and the final segmentation may be worse that in the case of the gradient image is directly used.

Regional graylevel-based model

Stochastic watershed enhances low gradients in the image that can be produced by problems as not uniform illumination or noise. This is justified since the pdf is directly extracted from the gradient image, which contains information about all the local changes in the intensity of the pixels. The stochastic procedure fulfils to reduce irrelevant contours that the watershed itself introduces but it enhances significant contours and also low gradients of the image, Fig. 2. For this reason, in this paper a novel strategy in the way of obtaining the pdf of object contours is presented. In this strategy, irrelevant edges

with 1≤n≤N C n = E and C p ∩ C q = ∅, ∀p = q. (6)
The stochastic procedure gives in each iteration an estimation of the gradient image which only measures the local energy of the object contours. In our approach and according to the representation of Eq. 6, each output image of the marker-controlled watershed transform can be represented as an image partition Π, i.e.,

I(x) → ρ(I(x)) → W S(x) → Π(W S(x)). ( 7
)
Each n-connected class C n of the partition Π is then evaluated with the mean intensity of pixels belonging to this class from the initial image I in order to construct a "mosaic image", denoted by p(x), and formally defined as of the image can be interpreted as a regional edgeness energy term, Eq. 9 (Fig. 3c).

p(x) = {µ Cn (I) : x ∈ C n }, (8)
E(x) = ρ(p(x)) (9)
This regional edgeness energy term is introduced in each of the M realizations of the stochastic watershed. That is, by integrating the gradient of the mean-based simplified model of each realization instead of directly the watershed lines, Eq. 10. The objective is to obtain a new probabilistic density function of object contours taking into account local properties of the image, Fig. 3d.

pdf(x) = 1 M M i=0 E i (x) * G(x; s). (10)

RESULTS

Several experiments demonstrate the benefits of the proposed method in comparison with other approaches. The effect of the two parameters (M , N ) of the stochastic watershed have been analysed to obtain the number of minima that different methods produced in the pdf of object contours. Fig. 4 shows the image and the gradient used in the next experiments. In this gradient image, the number of regional minima was calculated (the number of output regions of the watershed transform). After, the N parameter was fixed (N = 50) and the M parameter was changed between five and one hundred with a step of five. In each value of M , the stochastic watershed was applied five times independently in the gradient image and five pdf of object contours were obtained. Then, The number of minima of each pdf was calculated and the average in each value of M was obtained. This procedure tries to mitigate the random procedure of the whole algorithm. The same steps were carried out with our approach. In other experiment, the same protocol was applied but the M parameter was fixed (M = 30) and the N parameter was changed applying the stochastic watershed and our approach. Fig. 4c summarizes these results, representing the number of minima that our approach and the stochastic watershed produces in relation to the variation of the parameter M (fixing N ) and N (fixing M ) independently.

CONCLUSIONS

In this paper some of the problems of the watershed transform have been analysed. The problem associated with the number of regional minimum (since each one represents a region in the output image) is the most important to obtain an accurate segmentation with this procedure. This great number of minima is associated with a process of noise, low contrast or non-uniform illumination often introduced at the stage of image capture and it is hardly predictable. To solve this problem, different approaches try to reduce the number of regional minima merging regions or imposing artificial minima in the image. Other approaches seek to improve the image gradient to do it less noisy and/or better enhanced to get a lower number of minima. The stochastic watershed transform creates a probability density function of the object contours based on an iterative process and in the creation of random markers. This process allows to filter the artificial edges that can generate the watershed process itself and generally reduces the number of regions in comparison to the gradient image (Fig. 4c) but low gradients sometimes caused by illumination, contrast or noise are enhanced (Fig. 2). To mitigate this problem, some variations in the process of obtaining the pdf of the object contours are proposed in this paper. With this approximation, the number of regional minimum decreases respect to the stochastic watershed (Fig. 4c), since the noise of the pdf is lower than in the original stochastic watershed.

The effects of the parameters M, N in our approach are analysed in Fig. 4c. The M parameter can increase, but there is a value where the pdf is stable and if the M value increases results are similar and computational cost increases. The N parameter has a direct relationship with the number of minima detected in each iteration and the pdf, so the higher value of N, the higher number of regions/minima there is in the pdf.

In [START_REF] Bernander | Improving the stochastic watershed[END_REF] an improvement of the stochastic watershed based on stratified random markers and noise addition is presented. Fig. 5 shows the pdf of our approach and the result presented in [START_REF] Bernander | Improving the stochastic watershed[END_REF]. Our approach is less noisy and some objects are better detected than [START_REF] Bernander | Improving the stochastic watershed[END_REF]. In all cases under study, the way in which the random markers are generated can be decisive for the final result. In our work, the worst case was taken into account where such markers were generated uniformly and randomly throughout the image in comparison to [START_REF] Bernander | Improving the stochastic watershed[END_REF]. In addition, if a priori information was considered the improvement of our approximation can be more significantly than when markers are randomly generated. This would be case of medical images, where much information can be known a priori (organ position, intensity of organs under study, etc.).

Fig. 1 .

 1 Fig. 1. Image represented as topographic surface: (a) original image; (b) gradient image; (c) representation as topographic surface of gradient image.

Fig. 2 .

 2 Fig. 2. Stochastic watershed procedure (M=30, N=50): (a) original image; (b) gradient image; (c) j-th output image of the marker-controlled WS; (d) i-th output image of the marker-controlled WS; (e) pdf of stochastic watershed; (f) (pdf +gradient)/2.

Fig. 3 .

 3 Fig. 3. Regional graylevel-based model applied in image of Fig.1 (M=30, N=50): (a) pdf proposed in [9]; (b) i-th mosaic image, p i (x); (c) E i (x); (d) pdf(x).

Fig. 4 .

 4 Fig. 4. (a) original image; (b) gradient image; (c) Number of minima that different methods produce according M, N parameters.

Fig. 5 .

 5 Fig. 5. (a) image used in [11]; (b) pdf extracted from [11]; (c) pdf of our approach.

* Thanks to CDTI for funding this work under the project BRAIM (IDI-20130020).


Tom Cheesman

Avraham Roos

Version Variation Visualization (VVV): Case Studies on the Hebrew Haggadah in English

Keywords: visualization, algorithm, Hebrew, English, translation, haggadah, corpus

The following case study illustrates how 'Version Variation Visualisation' (VVV) software tools were successfully used to highlight previously unnoticed linguistic features in a corpus of 50 differing English-language (re)translations of the Hebrew Passover Haggadah. VVV's visualizations facilitate overviews, close reading, and navigation among versions. Its 'Eddy and Viv' algorithms enabled the authors to identify outliers (unusual versions) within the corpus, and pointed them to specific segments in the source text that generated the most variation among the translations. Such patterns could then be explained using close reading techniques with the help of cultural-historical background information. VVV is shown to be a useful tool in analysing multi-translation corpora.

INTRODUCTION

The 'Version Variation Visualization' project has developed online tools to support comparative, algorithm-assisted investigations of a corpus of multiple versions of a text, e.g. variants, translations, adaptations [START_REF] Cheesman | Reading Originals by the Light of Translations[END_REF][START_REF] Cheesman | Six Maps of Translations of Shakespeare[END_REF][START_REF] Cheesman | Translation Sorting with Eddy and Viv[END_REF][START_REF] Cheesman | Translation Sorting with Eddy and Viv[END_REF]Cheesman et al., -13, 2016;;[START_REF] Thiel | Visualizing Translation Variation: Designing Tools for Literary Scholars in Translation Studies and Linguistics[END_REF]links: www.tinyurl.com/vvvex). A segmenting and aligning tool allows users to 1) define arbitrary segment types, 2) define arbitrary text chunks as segments, and 3) align segments between a 'base text' (a stabilised iteration of the translated source text), and translated versions of the text. (The 'source text' of multiply translated texts is typically unstable, with historical variations in differing manuscripts and printed editions; our system requires the initial establishment of a 'base text' to serve as a stable basis for comparisons among translated versions.) The alignment tool can automatically align recurrent defined segment types in sequence. Several visual interfaces in the prototype installation enable exploratory access to parallel versions, to comparative visual representations of versions' alignment with the base text, and to the base text visually annotated by an algorithmic analysis of variation among versions of segments. Data can be filtered, viewed and exported in diverse ways. Many more modes of access and analysis can be envisaged. The tool is language neutral. Experiments so far mostly use modern texts: German Shakespeare translations. Roos is working on a collection of approx. 100 distinct English-language translations of a Hebrew text with ancient Hebrew and Aramaic passages: the Haggadah [START_REF] Roos | Making a Clean Breast of English Passover Haggadah Translations: Data Visualization of Bowdlerization in Haggadah Translations of Ezekiel[END_REF] I THE HAGGADAH On the evening before Passover (Pesach), Jews gather at home to celebrate a festive ceremony and meal with family and friends, to commemorate the biblical Exodus of the Jewish people out of Egypt. They eat the traditional matza, drink the prescribed four glasses of wine, and read from the Haggadah. This is a Hebrew text with instructions for a 15-phase ceremony: what to say or sing, which acts to perform, in what order, when to eat or drink what, etc. For example, in phase 1 the blessing over the wine is made, phase 11 is the meal and in phase 12 a specially broken piece of a Matzah, called Afikoman, is consumed. All participants hold a printed copy of the Haggadah. Typically, many different versions (plural: Haggadot) are present in the room.

The Hebrew Haggadah text is a compilation of Bible quotes, excerpts from traditional rabbinical teachings (Mishnah, Midrash), Exodus narrative, explanations of the festival's history, and Passover 'laws'. The text probably dates back to 200-300 CE. The oldest complete manuscript dates to the 10th century CE. Thousands of variant Hebrew-language versions are extant, in manuscript and print. There are translations in over 40 different languages [START_REF] Yudlov | 'Otz̲ ar ha-haggadot: Bibliografiah šel hagadot pesaḥ me-rešit ha-defus ha-ʻivri ʻad šenat tš[END_REF]. The first English-language version appeared in London in 1770. Countless more have since appeared. [START_REF] Yudlov | 'Otz̲ ar ha-haggadot: Bibliografiah šel hagadot pesaḥ me-rešit ha-defus ha-ʻivri ʻad šenat tš[END_REF] catalogues 823 English-language editions to 1960. The rate of production of new ones has since been accelerating exponentially. Most of these are retranslations, variously dependent on precursors.

Roos is compiling a digital corpus of English-language Haggadah translations, and using digital tools to compare them and visualize the differences. He aims to explain the differences in terms of their cultural historical contexts, and so shed light on translators' minds and motives.

II VERSION VARIATION VISUALIZATION (VVV) -Eddy and Viv Algorithms

VVV compares multiple retranslation documents at segment level, and visualizes the similarities and differences, in order to facilitate overviews, close reading, and navigation among versions. An algorithm called 'Eddy' ('∑D') quantifies variation among versions of a base text segment, in order to distinguish more and less predictable or distinctive versions. An algorithm called Viv ('variation in variation') aggregates Eddy metrics, and projects the result onto the base text segment, in order to distinguish more and less variously translated segments. The algorithms can be applied to the aligned corpus or any selected sub-corpus.

Eddy

The Eddy value assigned to a particular version or section indicates its "strangeness" as compared to other versions. Essentially, Eddy assigns lower metrics to wordings which are closer to the notional average, and higher metrics to more distant ones. So, Eddy ranks versions on a cline from low to high distinctiveness, or originality, or unpredictability. It sorts common-or-garden translations from interestingly different ones.

Eddy can be implemented in various ways. Our standard approach is: Each word in the corpus word list [where corpus means the corpus of aligned segment versions] is considered as representing an axis in N-dimensional space, where N is the length of the corpus word list. For each version, a point is plotted within this space whose co-ordinates are given by the word frequencies in the version word list for that version. (Words not used in that version have a frequency of zero.) The position of a notional 'average' translation is established by finding the centroid of that set of points. An initial 'Eddy' variation value for each version is calculated by measuring the Euclidean distance between the point for that version and the centroid. [Flanagan in [START_REF] Cheesman | Translation Array Prototype 1: Project Overview[END_REF] No stop words are excluded; no stemming, lemmatisation or parsing is performed. Users can also select a more primitive arithmetical formula, and one using Dice's co-efficient.

In the VVV 'Eddy and Viv' view, when a base text segment is selected, the segment-versions are displayed in a scrollable list in Eddy order, with associated metrics, and with a visual representation of relative Eddy value. The list can be re-ordered to display by date, translator name, or segment length in characters. Eddy values can also be displayed, explored, and exported in the form of charts and tables.

Examples of Eddy use will be provided in section 4.1.

Viv

Viv aggregates the Eddy values for a segment. In our standard approach, Viv is the average of Eddy values of version-segments. Users can also select Viv as the standard deviation of Eddy values. Viv indicates where translators differed most or least, in relation to the base text. (This function is comparable with the amount of layering associated with a word or string of words in the TRAViz visualization: [START_REF] Jänicke | TRAViz: A Visualization for Variant Graphs[END_REF])

In the VVV 'Eddy and Viv' view, Viv is represented on the base text by a tonal underlay, varying with the relative value of each segment. Metrics can be viewed by brushing a segment. Floor and ceiling values can be altered to facilitate surveying the base text.

Segments can be filtered in various ways (text search, Eddy/Viv ranges, segment lengths, etc), in the base text and in the version corpus or subcorpora, and texts exported in CSV tables with associated Eddy and Viv metrics.

Examples of Viv use will be provided in section 4.3.

As one reviewer commented, Eddy/Viv is not the only possible approach to comparing differing translations/versions. Measuring the overlap of words (or lemmas) among segments would achieve the same effect. Such a method would also need to calculate a centroid and distances from it.

VVV is specifically created to compare numerous retranslations of the same source text, making it ideal for research into Haggadah version variation. It can help a researcher identify variations, and present them to an audience.

III RELATED WORK

There has been some digital work on larger retranslations corpora, involving works of wide intrinsic interest, but none designed to facilitate access to multiple translations, and the translated work, together with algorithmic analyses. [Janicke et al., 2015] take an in some ways similar approach, but their 'TRAViz' interface offers a very different mode of text visualization, is monolingual (shows no translated text), and works best with more limited variation and shorter texts. Similarly, Juxta-Commons, CollateX and Stemmaweb are monolingual, do not rate the "strangeness" of variants in comparison to all others (Eddy value) and do not create a heat map in the source text (Viv Value) revealing which ST instances generated most target text variants, two of the most powerful VVV features. It is especially these comparisons between the ST and the TT that this research focuses on. Whereas the above software highlights differences in the versions, VVV highlights how these differences are connected to the ST, thus helping us the explain the reasons for the variants.

[Lapshinova-Koltunski, 2013] describes a parallel multi-translation corpus designed to support computational linguistic analyses of differences between professional translations, student translations, MT outputs and edited MT outputs. [Shei and Pain, 2002] proposed a similar parallel corpus, with an interface designed for translator training. These projects only offer access to filtered segments of the text corpus, and do not envisage exploring variation among retranslations. [Altintas, Can, and Patton, 2007] used two time-separated (c.1950, c.2000) collections of published translations of the same seven English, French or Russian literary classics into Turkish, in order to quantify aspects of language change. This raises the question whether such translations 'represent' their language. Corpus-based Translation Studies [Baker, 1993;Kruger et al., 2011] has established that translated language differs from untranslated language. We also know from decades of work in Descriptive Translation Studies [Morini, 2014;Toury, 2012] that retranslations vary for complex genre-, market-, subculture-specific and institutional factors, and individual psychosocial factors, involving the translators and others with a hand in the work (commissioners, editors), and their uses of resources including source versions and prior (re)translations.

IV USING VVV WITH HAGGADAH SAMPLES

One section of the Haggadah concerns four sons who represent four different attitudes to Judaism. They each ask a question which characterizes their attitude towards the Passover festivities, and the text then suggests how to respond to these questions. This section of the Hebrew source text has 126 words and is divided into six parts: (1) introduction; (2) characterization of each of the four sons; (3)-( 6) one paragraph for each son, with his question and the response. The source text contains twelve manually defined segments: units of meaning.

60 different translations of this section were uploaded to VVV, segmented and manually aligned with the Hebrew base text. Each translation contains between five and twelve of the source text segments because translations sometimes disregard segments of the ST or merge two ST segments into one TT translation.

Exploring with Eddy

In part (2), characterizing the four sons, most translators use straightforward terms: 'wise', 'wicked', 'simple', and 'one who does not know how to ask'. Some are more creative. Eddy highlights certain translations as 'strange'. VVV automatically rates each version, thus ranking all 60 versions from most common to most exotic. In The corpus includes C18 and C19 versions, but none appear in table 1. Almost all high Eddy versions date from the 1940s and after. The general retranslation trend is towards greater variation, probably at least in part because of copyright issues and a need to differentiate in order to stand out in the ever growing crowd.

As a distant reading tool, VVV's Eddy values reveals to us that in comparison to other versions of the same period, the 1906 translation (WILROS) is an early outlier, therefore worth further investigation. Close reading reveals that the language use in this particular translation is indeed quite extraordinary, with phrases such as "Israel's exode from Egypt", "and took cognizance of them", "of which Jerusalem is emblematic", "cut the sea in twain", etc. Historically and culturally, the translator William Rosenau was a radical leader of Reform Judaism, and he would eventually edit the revised edition of the Reform Union Haggadah, with a thoroughly rewritten source text. The version in our corpus still adheres to the traditional ST, but Rosenau's radicalism clearly already shines through in his translation and is picked up by VVV.

It is also intriguing that no version is consistently in the highest 5 for all four sons (see Table 2). A translation's relative Eddy varies, as we will see in the next section.

Eddy Variation Chart

The poet Abraham Regelson published several Haggadot. Roos's collection includes one from 1944 (REGFORST1) and another from 1952 (REGFORST2). VVV's Variation Chart ('Eddy Overview') helps us compare these two translations (see Figure 1). This chart plots each version's Eddy values on the y-axis, for segments in sequence on the x-axis. The user can select which versions' graphs to display or hide, select an area to zoom in on, and hover over a node to display base text and version (as is shown in Figure 1). In Figure 1 we see Regelson using higher Eddy-value language (more distinctive language in relation to the corpus) in 1952 than in 1944. One exception is highlighted, in part (2) of the passage (discussed in Section 4.3). In Table 3 What says the Simple Son? -"What is all this about?" Therefore, say to him: "With might of hand, Praised be God, praised be He. Blessed be He who gave the Torah to His people, Israel. Blessed be He. On the subject of the Passover service the Torah speaks of -FOUR SONS

One is intelligent, one is ill-mannered, one is indifferent, and one is not even able to ask a question.

1. The INTELLIGENT son asks: "What is the meaning of all the Passover customs and ceremonies, the rules and rites which God has commanded?" You will explain to him all the traditions of Pesach down to the last detail of the Afikoman.

The

Viv: Variation in Variation

In VVV's 'Eddy and Viv' interface, source text segments are highlighted according to their Viv value: the higher the value, the darker the underlay tone. We can visually identify which source text segments produced the most variant translations, or 'read the original by the light of the translations' [START_REF] Cheesman | Reading Originals by the Light of Translations[END_REF]. The segment with highest Viv value is the one beginning: 'We may not eat an afikoman…' (underlined above). This is a quote from the Mishna (oral laws compiled about 200 CE by Rabbi Judah HaNasi, the basis for the later Talmud). Already in the Talmud (c. 500 CE) the correct meaning of the term afikoman had become obscure and was disputed. In Talmudic traditions, afikoman (Hebrew )אפיקומן is said to derive from Greek epikomen or epikomion (επί Κομός), 'that which comes after', variously interpreted as (A) 'dessert', or (B) 'afterdinner entertainment/revelry', and additionally as (C) a metaphor. The translation from the Sefaria site cited above clearly favors the option that it is a form of dessert.

But our corpus contains five different interpretations in the context of the answer to the wise son: (A1) a proscribed dessert; (A2) the prescribed dessert eaten at the 13 th phase of the Passover ceremony (the piece of a matza called Afikoman); (B1) proscribed excessive subsequent entertainment (distinguishing Passover from pagan celebrations); (C1) prescribed teaching of all of the (Passover) law: because the afikoman is the last law in the section on Pesach, so "We don't leave anything until after the afikoman" could mean, "we don't stop studying until we have learned everything"; (C2) prescribed sacrifice of a Passover lamb.

There is also a sixth option for translation: leaving afikoman to stand in the target text, uninterpreted.

This range of options explains the segment's high Viv value. Some of the variant English versions, low in the Eddy value list, are shown in the VVV 'Eddy and Viv' view in Figure 3. Viv and Eddy values are calculated according to manually set segments (meaningful units). These can be single words, phrases, sentences, paragraphs or even the whole text. By creating a whole text segment, one could easily aggregate data for a quantitative comparison of whole texts. This way, although no one version consistently deviates from the norm for all 4 sons (as shown above), it is possible to ascertain which translation as a whole is furthest from the base texts. It should be noted, however, that by doing this, the Viv value becomes irrelevant. The version that scored the highest overall Eddy value was the 1974 Polychrome Haggadah by Jacob Freedman, whose translation is extremely verbose and elaborate. 1) expands to an extraordinary degree. This is the verbose and elaborate translation by Jacob Freedman mentioned at the end of section 4.3. The 1955 and 1967 versions are also expansive. They did not appear in Table 1. Figure 4 now explains why: both omit the segments which ST's part 2 with the 4 segments of the names of the four sons, shown in Table 1.

Conclusion

Using VVV can yield valuable insights when comparing multiple variants, and is also useful for presenting findings visually. Manually comparing different versions becomes difficult with larger corpora. When Viv is highlighted in base text segments, even researchers with no knowledge of a language (in this case Hebrew) can identify the parts that warrant closer inspection.

VVV offers a useful range of visualization modes, but many more can be developed. Future research planned on the Haggadah includes comparing the language use of translators when translating Hebrew and Aramaic text passages, comparing the translations of biblical Hebrew versus more modern Hebrew text passages, and translated transitive versus intransitive verb forms.

Thanks

We would like to thank VVV developer Kevin Flanagan for continuing technical support, and Eran Hadas for his endeavors to streamline uploading and segmenting texts in VVV.

Figure 1 .

 1 Figure 1. Comparing REGFORST1 and REGFORST2 in the Eddy Variant Chart

Figure 2 .

 2 Figure 2. Partial Screenshot of Eddy and Viv view of the Four Sons Passage

Figure 3 .

 3 Figure 3. Eddy and Viv view of the Four Sons Passage

4. 4

 4 Parallel View Visualization: Alignment Maps Parallel view visualizations include a distant overview of segment alignments between base text and versions: an 'alignment map'. Successive segments of the base text are represented as a vertical 'barcode': the thickness of a bar represents segment length in words. Segments of a version are represented in the same format. Alignments are represented by lines connecting base text and version. This enables rapid identification of translators' editing decisions: omission, addition, reduction, expansion, and transposition.

Figure 4 (Figure 4 .

 44 Figure 4 (created from screenshots) shows ten examples of the 'Four Sons' passage. The unchanging base text is on the left, the version on the right of each 'map'. The afikoman segment is highlighted.

Table 1 .

 1 Table 1, the left column gives the original Hebrew and the commonest translation (lowest Eddy value). The second column gives the five translations with highest Eddy values (rounded Eddy figures given in column 3): outliers in comparison to all other versions in the corpus, worth further exploration. Names of the Four Sons: Translations with lowest and Highest Eddy.

	Base text and
	lowest translation Eddy	5 highest Eddy translations	Eddy value	Version reference	Version date
					One is intelligent	1.28 REGFORS	1952
						T2
		חכם	אחד	One is understanding	1.30 POLYHH	1974
	One wise			The intelligent child	1.38 WILROS	1906
					A clever son	1.45 TCH	1954
					The first is sensible	1.52 MSAM	1942
					The second mean	1.38 MSAM	1942
					The rebellious child	1.38 ANIM	2005
	One wicked רשע	ואחד	One is ill-mannered	T2 1.50 REGFORS	1952
					and one who is stubborn	1.63 GUT	1956
					one is recalcitrant and scornfully insolent 1.81 POLYHH	1974
					One Artless	1.17 NAH	2012
					One is indifferent	1.32 REGFORS	1952
		תם	ואחד		T2
	One simple			A simpleton	1.32 LEHM	1972
					The naïve son	1.33 HOS	2009
					A dull son	1.45 TCH	1954
					and a fourth, a child that does not yet	1.90 GLATZ	1989
					know how to ask
					and the fourth incapable of even asking a	1.91 POLYHH	1974
	יודע	שאינו	ואחד	question
		לשאול	and one who is too young to ask any	1.92 TCH	1954
	and one who	questions about the things he sees
	does not know how to ask	to make inquiry of his own accord and the child who does not know enough	1.94 WILROS	1906
					and the child still too young to even	2.13 NSEX	1983
					inquire of the Pesakh

Table 2 .

 2 Highest Eddy Scorers from Table1.

	2 son1 REGFORST	POLYHH TCH	WILROS MSAM
	2 son2 REGFORST	POLYHH	MSAM ANIM	GUT
	2 son3 REGFORST	TCH		NAH LEHM HOS
	son4	POLYHH TCH	WILROS	GLATZ	NSEX

Table 3 .

 3 With a strong hand God took us out of Egypt the Lord hath taken us out of Egypt, from the house of slavery." But the One Who Wits Not To Ask-sit is for thee to open talk with him, as it is said: "And thou shalt tell thy son on that day, saying: 'This is on account of what the Lord did for me when I went forth from Egypt.''' where we were slaves." 4. The INCOMPETENT one-get him started by quoting the words from the Bible: "In that day you shall tell your son saying: (Point to the ceremonial dishes.) "All of this is because of what the Lord did for me when I came out of Egypt." Comparing Regelson 1944 and Regelson 1952.It makes sense to assume that, having in 1944 already translated the Haggadah in a quite straightforward manner, Regelson decided for the 1952 translation to try out different translation techniques, more off the beaten path. That would explain why his 1952 translation scores higher Eddy values. Examples of this can be seen by the alliterated names of the four sons (Intelligent, ill-mannered, indifferent, and incompetent), something not found in any of the other versions. On the other hand, having used a etymology-based final comment for the wise son in 1944, Regelson opts for a very specific Jewish-sources based closing for his 1952 intelligent son. We will comment on this further in the next section.

	ILL-MANNERED one asks: "What's the sense of
	all this business of yours?" Yours, he says; and none of
	his. By refusing to identify himself with his people he
	denies a basic principle of our religion. You may fling this
	in his teeth: "I do this because of what the Lord did for me
	when He rescued me from Egypt." Me not him. Let him
	know that had he been there, by denying his brothers he
	could not have been saved.

3. The INDIFFERENT one merely asks: "What is this?" Tell him: "

Journal of Data Mining and Digital Humanities http://jdmdh.episciences.org

http://jdmdh.episciences.org

Funding

The VVV project was supported in 2011-13 by Swansea University (Research Incentive Fund and Bridging the Gaps). The main phase of software development was funded by a 6-month Research Development Grant under the Digital Transformations theme of the Arts and Humanities Research Council (UK), reference AH/J012483/1, in 2012.


Virginie Majchrzak

published or not. The documents may come

DC and AC Magnetic Flux Interaction in a Coupling

Transformer with a Virtual Air Gap Principle of a Dynamic Voltage Restorer


Virginie Majchrzak

research documents, whether they are published or not. The documents may come

Coupling Transformer Operation of a Dynamic Voltage

Restorer Under Electrical Grid Conditions