# Analysis of Large-Amplitude Pulses in Short Time Intervals: Application to Neuron Interactions

Authors: Gianni Mattioli; Massimo Scalia; Carlo Cattani
Journal: Mathematical Problems in Engineering
 (2010)

---

## Abstract

This paper deals with the analysis of  a nonlinear dynamical system which characterizes  the axons  interaction and is  based on a generalization of FitzHugh-Nagumo system. The parametric domain of stability is investigated for both the linear and third-order approximation. A further generalization is studied in presence of  high-amplitude (time-dependent) pulse. The corresponding numerical solution for some given values of parameters are analyzed through the wavelet coefficients, showing both the sensitivity to local jumps and some unexpected inertia of neuron's as response to the high-amplitude spike.


---

## Body

## 1. Introduction

Theclassical model of one excitatory neuron [1–11] is considered under a spike train (time-dependent) excitation. The dynamical system as well as its evolution are investigated under the dependence on some featuring parameters and in presence of a high-amplitude time-dependent pulse. It was already shown, in a previous paper [12], dealing with Hodgkin-Huxley model, that neuron's multiple firing does not reflects immediately on a direct response and there is some kind of time delay (inertia) before the pulse becomes effective. It should be noticed that there are different models of neurons (see, e.g., [1] and references therein) such as the Integrate-and-fire, FitzHugh-Nagumo [2, 9], Morris-Lecar, and the more general Hodgkin-Huxley model [6–8]. In the early sixties, FitzHugh [2] proposed, as a simplified model, a generalization of the Van der Pol equation thus showing the existence of a limit cycle and some periodicity. In all models, if we consider a simple system consisting of a neuron and a synapse, it is known that the activity of stimulated axons can be detected by an abrupt change in the electrical potential. These pulse in a short time are called spikes or axons firing. However, as shown later these abrupt changes do not appear immediately, thus showing, even in this simplified model of neurons interactions, the existence of some kind of neurons inertia, which can be physically explained by their cooperative responce to firing.In the following we will propose a suitable generalization of the FitzHugh-Nagumo model (FHN) and we will study both on the linear and third-order approximation. Equilibrium points, parametric domain of stability will be outlined in absence of time-dependence. The presence of high-amplitude time-dependent source of pulses can be studied on the nonlinear time series of the numerical approximation. The pulse action of the firing process, based on FitzHugh suggestion of delta Dirac function, is proposed in the form of a high-amplitude localized function with compact support in a short interval, which replicates a few times before it disappears. According to FitzHugh any change in the electrical potential are localized in a short time interval and within this interval has a significant amplitude.


## 2. Simplified FitzHugh-Nagumo Model

The most successful and general model in neuroscience describes the neuron activity in terms of two conductances and electric potentials. The Hodgkin-Huxley model [6–8] considers the neuron activity as an electrical circuit. Cells membrane store charges, electrochemical forces arise because of the imbalanced ion concentration inside and outside the cell. It can be written, in adimensional form, as  (2.1)dudt=a1n4(a2-u)+a3m3h(a4-u)+a5(a6-u)-Iext(t),dndt=αn(u)(1-n)-βn(u)n,dmdt=αm(u)(1-m)-βm(u)m,dhdt=αh(u)(1-h)-βh(u)h,
being
(2.2)αn(u)=ϕ·a7(u-a8)1-e-(u-a8)/a9,βn(u)=ϕ·a10e-(u-a11)/a12,αm(u)=ϕ·a13(u-a14)1-e-(u-a14)/a15,βm(u)=ϕ·a16e-(u-a17)/a18,αh(u)=ϕ·a19e-(u-a20)/a21,βh(u)=ϕ·a221+e-(u-a23)/a24.An alternative simplified 2-dimensional model (FHN) that has been proposed by FitzHugh-Nagumo [2, 9] can be substantially represented by the system(2.3)dxdt=a(x-x33)-by+γz,dydt=αx-βy,
where z=z(t) is a pulse function defined in a very short-range interval (pulse function), having as limiting case the delta Dirac. In the following we will consider a spike train. Parameter γ defines the amplitude of the pulse function z(t).This dynamical system depends [2, 9] on a,b,α,β,andγ in the sense that the evolution would be completely different, starting from some critical values. In the following we will show that the solution tends asymptotically to a constant value.A suitable generalization of system  (2.3) is(2.4)dxdt=asinx-bh(y)+γz(t),dydt=αf(x)-βg(y),
with nonnegative parameters  (2.5)a≥0,b≥0,β≥0,
and odd functions
(2.6)f(x)=-f(-x),h(y)=h(-y),g(y)=g(-y),f(x0)=f′′(x0)=f′′′(x0)=0,h(x0)=h′′(x0)=h′′′(x0)=0,g(x0)=g′′(x0)=g′′′(x0)=0.
Up to the first-order it is(2.7)dxdt=ax-by+γz(t),dydt=αx-βy,
and to the third-order it is(2.8)dxdt=a(x-x36)-by+γz(t),dydt=αx-βy,
so that the FHN system  (2.3) might be considered as the third-order approximation of system  (2.4).


### 2.1. Autonomous Systemγ=0

#### 2.1.1. Linear Case

Whenγ=0, and assuming by analogy with system (2.3)  that(2.9)a=1,b=1,
from the linear system  (2.7), we have(2.10a)dxdt=x-y,(2.10b)dydt=αx-βy. There is only one equilibrium point at the origin:
(2.11)(x̅,y̅)=(0,0).
The eigenvalues are
(2.12)λ1,2=12(1-β±Δ),Δ≡1-4α+2β+β2.
The parametric domain of feasible values for α and β is defined according to  (2.5) and by the parabolic curve as in Figure 1: (2.13)α=14(β+1)2.
Thus we have  (2.14)λ1=λ2∈ℝ,α=β24+β2+1,
where β>1, stable equilibrium in (0,0), and β<1, unstable equilibrium in (0,0);        (2.15)λ1≠λ2∈ℝ,α>β24+β2+1,
where λ1<0,λ2<0, stable equilibrium in (0,0),λ1<0,λ2>0, unstable equilibrium in (0,0) (saddle), and λ1>0,λ2>0, unstable equilibrium in (0,0);       (2.16)λ1=λ2∈ℂ,α<β24+β2+1,
where β>1,stable equilibrium in (0,0), and β<1, unstable equilibrium in (0,0).Figure 1
Parametric domain for the linear FHN model (8).As a special case, by assuming thatt0=0,x(0)=x0,andy(0)=y0 with α=0, β≠0, from (2.10b) we have (2.17)y(t)=y0e-βt,
and from (2.10a) we have (2.18)x(t)={x0et+y01+β[e-βt-et],β≠-1,et(x0-y0t),β=-1.
Analogously, by assuming that t0=0,x(0)=x0,andy(0)=y0 with β=0, α≠0, from (2.10b) we have (α≠1/4) (2.19)x(t)=et/21-4α[x01-4αcosh(1-4α2t)+(x0-2y0)sinh(1-4α2t)],y(t)=et/21-4α[y01-4αcosh(1-4α2t)+(2αx0-y0)sinh(1-4α2t)],x(t)=et/22[x0(2+t)-2y0t],y(t)=et/24[x0t-2y0(t-2)],
when α=1/4.


#### 2.1.2. Third-Order Approximation

Assuming that(2.20)a=1,b=1,
from the third-order system  (2.8), we have(2.21)dxdt=(x-x36)-y,dydt=αx-βy.
The null clines intersect at the points (2.22)(x̅1,y̅1)=(0,0),(x̅2,y̅2)=(-6(β-α)β,-αβ6(β-α)β),(x̅3,y̅3)=(6(β-α)β,αβ6(β-α)β),
so that, being β>0, there are 3 disjoint equilibrium points (Figure 2) when (2.23)β>α,
and only one when (2.24)α≤-β.
The Jacobian of  (2.21) is (2.25)J(x,y)=(1-x22-1α-β)
so that (2.26)J(x̅2,y̅2)=J(x̅3,y̅3)=(1β(3α-2β)-1α-β).
The eigenvalues are(2.27)λ1,2=12β(3α-2β-β2±Δ1),Δ1≡(3α-2β)2+β2(β2-4β+2α).
The parametric domain of feasible values for α and β is defined by the curve (2.28)(3α-2β)2+β2(β2-4β+2α)=0
as in Figure 3.Figure 2
Null clines for (2.21).Figure 3
Parametric domain of (2.21).


#### 2.1.3. The General Case  (2.4)

Whenγ=0, assuming that (2.29)a=1,b=1,f(x)=x,g(y)=y
from the system  (2.4), we have (2.30)dxdt=sinx-y,dydt=αx-βy.
The null clines have at least two intersections (Figure 4) when (2.31)α≤βorβ≥-4.6α
andonly one when
(2.32)α≤β.
otherwise, there are more than one intersection.Figure 4
Null clines for (2.8).The Jacobian of  (2.30) and eigenvalues coincide with the linear case.


## 2.1. Autonomous Systemγ=0

### 2.1.1. Linear Case

Whenγ=0, and assuming by analogy with system (2.3)  that(2.9)a=1,b=1,
from the linear system  (2.7), we have(2.10a)dxdt=x-y,(2.10b)dydt=αx-βy. There is only one equilibrium point at the origin:
(2.11)(x̅,y̅)=(0,0).
The eigenvalues are
(2.12)λ1,2=12(1-β±Δ),Δ≡1-4α+2β+β2.
The parametric domain of feasible values for α and β is defined according to  (2.5) and by the parabolic curve as in Figure 1: (2.13)α=14(β+1)2.
Thus we have  (2.14)λ1=λ2∈ℝ,α=β24+β2+1,
where β>1, stable equilibrium in (0,0), and β<1, unstable equilibrium in (0,0);        (2.15)λ1≠λ2∈ℝ,α>β24+β2+1,
where λ1<0,λ2<0, stable equilibrium in (0,0),λ1<0,λ2>0, unstable equilibrium in (0,0) (saddle), and λ1>0,λ2>0, unstable equilibrium in (0,0);       (2.16)λ1=λ2∈ℂ,α<β24+β2+1,
where β>1,stable equilibrium in (0,0), and β<1, unstable equilibrium in (0,0).Figure 1
Parametric domain for the linear FHN model (8).As a special case, by assuming thatt0=0,x(0)=x0,andy(0)=y0 with α=0, β≠0, from (2.10b) we have (2.17)y(t)=y0e-βt,
and from (2.10a) we have (2.18)x(t)={x0et+y01+β[e-βt-et],β≠-1,et(x0-y0t),β=-1.
Analogously, by assuming that t0=0,x(0)=x0,andy(0)=y0 with β=0, α≠0, from (2.10b) we have (α≠1/4) (2.19)x(t)=et/21-4α[x01-4αcosh(1-4α2t)+(x0-2y0)sinh(1-4α2t)],y(t)=et/21-4α[y01-4αcosh(1-4α2t)+(2αx0-y0)sinh(1-4α2t)],x(t)=et/22[x0(2+t)-2y0t],y(t)=et/24[x0t-2y0(t-2)],
when α=1/4.


### 2.1.2. Third-Order Approximation

Assuming that(2.20)a=1,b=1,
from the third-order system  (2.8), we have(2.21)dxdt=(x-x36)-y,dydt=αx-βy.
The null clines intersect at the points (2.22)(x̅1,y̅1)=(0,0),(x̅2,y̅2)=(-6(β-α)β,-αβ6(β-α)β),(x̅3,y̅3)=(6(β-α)β,αβ6(β-α)β),
so that, being β>0, there are 3 disjoint equilibrium points (Figure 2) when (2.23)β>α,
and only one when (2.24)α≤-β.
The Jacobian of  (2.21) is (2.25)J(x,y)=(1-x22-1α-β)
so that (2.26)J(x̅2,y̅2)=J(x̅3,y̅3)=(1β(3α-2β)-1α-β).
The eigenvalues are(2.27)λ1,2=12β(3α-2β-β2±Δ1),Δ1≡(3α-2β)2+β2(β2-4β+2α).
The parametric domain of feasible values for α and β is defined by the curve (2.28)(3α-2β)2+β2(β2-4β+2α)=0
as in Figure 3.Figure 2
Null clines for (2.21).Figure 3
Parametric domain of (2.21).


### 2.1.3. The General Case  (2.4)

Whenγ=0, assuming that (2.29)a=1,b=1,f(x)=x,g(y)=y
from the system  (2.4), we have (2.30)dxdt=sinx-y,dydt=αx-βy.
The null clines have at least two intersections (Figure 4) when (2.31)α≤βorβ≥-4.6α
andonly one when
(2.32)α≤β.
otherwise, there are more than one intersection.Figure 4
Null clines for (2.8).The Jacobian of  (2.30) and eigenvalues coincide with the linear case.


## 2.1.1. Linear Case

Whenγ=0, and assuming by analogy with system (2.3)  that(2.9)a=1,b=1,
from the linear system  (2.7), we have(2.10a)dxdt=x-y,(2.10b)dydt=αx-βy. There is only one equilibrium point at the origin:
(2.11)(x̅,y̅)=(0,0).
The eigenvalues are
(2.12)λ1,2=12(1-β±Δ),Δ≡1-4α+2β+β2.
The parametric domain of feasible values for α and β is defined according to  (2.5) and by the parabolic curve as in Figure 1: (2.13)α=14(β+1)2.
Thus we have  (2.14)λ1=λ2∈ℝ,α=β24+β2+1,
where β>1, stable equilibrium in (0,0), and β<1, unstable equilibrium in (0,0);        (2.15)λ1≠λ2∈ℝ,α>β24+β2+1,
where λ1<0,λ2<0, stable equilibrium in (0,0),λ1<0,λ2>0, unstable equilibrium in (0,0) (saddle), and λ1>0,λ2>0, unstable equilibrium in (0,0);       (2.16)λ1=λ2∈ℂ,α<β24+β2+1,
where β>1,stable equilibrium in (0,0), and β<1, unstable equilibrium in (0,0).Figure 1
Parametric domain for the linear FHN model (8).As a special case, by assuming thatt0=0,x(0)=x0,andy(0)=y0 with α=0, β≠0, from (2.10b) we have (2.17)y(t)=y0e-βt,
and from (2.10a) we have (2.18)x(t)={x0et+y01+β[e-βt-et],β≠-1,et(x0-y0t),β=-1.
Analogously, by assuming that t0=0,x(0)=x0,andy(0)=y0 with β=0, α≠0, from (2.10b) we have (α≠1/4) (2.19)x(t)=et/21-4α[x01-4αcosh(1-4α2t)+(x0-2y0)sinh(1-4α2t)],y(t)=et/21-4α[y01-4αcosh(1-4α2t)+(2αx0-y0)sinh(1-4α2t)],x(t)=et/22[x0(2+t)-2y0t],y(t)=et/24[x0t-2y0(t-2)],
when α=1/4.


## 2.1.2. Third-Order Approximation

Assuming that(2.20)a=1,b=1,
from the third-order system  (2.8), we have(2.21)dxdt=(x-x36)-y,dydt=αx-βy.
The null clines intersect at the points (2.22)(x̅1,y̅1)=(0,0),(x̅2,y̅2)=(-6(β-α)β,-αβ6(β-α)β),(x̅3,y̅3)=(6(β-α)β,αβ6(β-α)β),
so that, being β>0, there are 3 disjoint equilibrium points (Figure 2) when (2.23)β>α,
and only one when (2.24)α≤-β.
The Jacobian of  (2.21) is (2.25)J(x,y)=(1-x22-1α-β)
so that (2.26)J(x̅2,y̅2)=J(x̅3,y̅3)=(1β(3α-2β)-1α-β).
The eigenvalues are(2.27)λ1,2=12β(3α-2β-β2±Δ1),Δ1≡(3α-2β)2+β2(β2-4β+2α).
The parametric domain of feasible values for α and β is defined by the curve (2.28)(3α-2β)2+β2(β2-4β+2α)=0
as in Figure 3.Figure 2
Null clines for (2.21).Figure 3
Parametric domain of (2.21).


## 2.1.3. The General Case  (2.4)

Whenγ=0, assuming that (2.29)a=1,b=1,f(x)=x,g(y)=y
from the system  (2.4), we have (2.30)dxdt=sinx-y,dydt=αx-βy.
The null clines have at least two intersections (Figure 4) when (2.31)α≤βorβ≥-4.6α
andonly one when
(2.32)α≤β.
otherwise, there are more than one intersection.Figure 4
Null clines for (2.8).The Jacobian of  (2.30) and eigenvalues coincide with the linear case.


## 3. FHN Model with a High-Amplitude Spike

Let us consider the general nonautonomous system  (2.4) when γ≠0:(3.1)dxdt=sinx-y+γz(t),dydt=αx-βy.
The time-dependent function z(t) is a high-amplitude function (Figure 5), based on Haar wavelets  (3.2)z(t)≡∑k=-3310ψ84(t-k)+5ψ94(t-k)+ψ140(t-k)
which replicates itself in a finite interval (Figure 6).Figure 5
High-amplitude function10ψ84(t)+5ψ94(t)+ψ104(t).Figure 6
z(t).The basic wavelet functionψkn(t) is defined as
(3.3)ψkn(t)≡2n/2{1,k2n≤t<k+0.52n,-1,k+0.52n≤t<k+12n,0,elsewhere.


### 3.1. Numerical Simulation

Let us assume for the values of parameters the following:γ=102, β=2 and two different values for α, that is, α=-0.01, α=2 in the interval t∈[0,8]. As initial conditions it is assumed that x(0)=0, y(0)=0. By using the Runge-Kutta 4th-order method, with the accuracy 10-6, we obtain in the interval (0<t≤8), the solution in correspondence with different initial conditions.From a direct inspection of the solution (see Figures7 and 8) it can be seen that under a spike firing there is some delay effect so that the perturbation show its influence only after some time delay. This perturbation of the system makes the orbits nearby the origin unstable. The uniqueness in phase space is going to be lost at least in the initial time interval. The orbit in the phase plane is self-crossing. The uniqueness of motion is lost. When α=-0.01, the origin behaves as an attractor, while in the case α=2,  x, y diverge asymptotically.Numerical solution of (3.1) with γ=102, β=2, and α=-0.01 (a, c, e), α=2 (b, d, f) in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)(c)(d)(e)(f)Phase orbits of (3.1) with γ=102, β=2,and α=-0.01 (a), α=2 (b) in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)


## 3.1. Numerical Simulation

Let us assume for the values of parameters the following:γ=102, β=2 and two different values for α, that is, α=-0.01, α=2 in the interval t∈[0,8]. As initial conditions it is assumed that x(0)=0, y(0)=0. By using the Runge-Kutta 4th-order method, with the accuracy 10-6, we obtain in the interval (0<t≤8), the solution in correspondence with different initial conditions.From a direct inspection of the solution (see Figures7 and 8) it can be seen that under a spike firing there is some delay effect so that the perturbation show its influence only after some time delay. This perturbation of the system makes the orbits nearby the origin unstable. The uniqueness in phase space is going to be lost at least in the initial time interval. The orbit in the phase plane is self-crossing. The uniqueness of motion is lost. When α=-0.01, the origin behaves as an attractor, while in the case α=2,  x, y diverge asymptotically.Numerical solution of (3.1) with γ=102, β=2, and α=-0.01 (a, c, e), α=2 (b, d, f) in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)(c)(d)(e)(f)Phase orbits of (3.1) with γ=102, β=2,and α=-0.01 (a), α=2 (b) in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)


## 4. Wavelet Analysis

TheHaar scaling functionφ(t) is the characteristic function on [0,1]. By translation and dilation we get the family of functions defined (in [0,1])  as(4.1)φkn(t)≡2n/2φ(2nt-k),(0≤n,0≤k≤2n-1),φ(2nt-k)={1,t∈Ωkn,0,t∉Ωkn,Ωkn≡[k2n,k+12n).
The Haar wavelet family {ψkn(t)} is the orthonormal basis [13]:  (4.2)ψkn(t)≡2n/2ψ(2nt-k),∥ψkn(t)∥L2=1,ψ(2nt-k)≡{-1,t∈[k2n,k+1/22n),1,t∈[k+1/22n,k+12n),0,elsewhere.(0≤n,0≤k≤2n-1),
Without loss of generality, we can restrict ourselves to 0≤n, 0≤k≤2n-1⇒Ωkn⊆[0,1]. Let Y≡{Yi},  (i=0,…,2M-1,2M=N<∞,M∈ℕ), be a finite energy time-series; ti=i/(2M-1), is the regular equispaced grid of dyadic points.Thediscrete Haar wavelet transform is the operator 𝒲N which maps the vector Y into the vector of the wavelet coefficients{α,βkn}:  (4.3)𝒲NY={α,β00,…,β2M-1-1M-1},Y={Y0,Y1,…,YN-1},
However in order to reduce the computational complexity it has been proposed [14] the short wavelet transform as follows. Let the set Y={Yi} of N data be segmented into σ segments (in general) of different length. Each segment Ys,s=0,…,σ-1 is made of ps=2ms, (∑sps=N), data: (4.4)Y={Yi}i=0,…,N-1=⨁s=0σ-1{Ys},Ys≡{Ysps,Ysps+1,…,Ysps+ps-1},
being, in general, ps≠pr. The short discrete Haar wavelet transform of Y is (see [14]) 𝒲ps,σY, (4.5)𝒲ps,σ≡⨁s=0σ-1𝒲sp,Y=⨁s=0σ-1Ys,𝒲ps,σY=(⨁s=0σ-1𝒲ps)Y=(⨁s=0σ-1𝒲psYs),𝒲2msYs={α00(s),β00(s),β01(s),β11(s),…,β2ms-1-1ms-1(s)},
with 2ms=ps,∑s=0σ-1ps=N.There follows that, the matrix of the wavelet transform is expressed as a direct sum of lower-order matrices so that the short transform is a sparse matrix [14]. When the short wavelet transform maps short interval values into a few set of wavelet coefficients, it can be considered as a first-order approximation. However, since the wavelet transform maps the original signal into uncorrelated sequences, the short wavelet transform describes, for each sequence of detail coefficients, its local behavior. When ps=p=N,σ=1, the above coincides with the ordinary wavelet transform. We assume, in the following, that ps=p=N/σ,s=0,…,σ-1,(σ>1).


### 4.1. Critical Analysis through the Wavelet Coefficients

Let us take the short Haar wavelet transform for the two time series, obtained by using the Runge-Kutta4th order method, with the accuracy 10-6, of system  (3.1) with the following values of parameters-initial conditions: (4.6)γ=102,β=2,x0=0,y0=0,t∈[0,8].
The two time series correspond to the two values of α:α=-0.01 and α=2. It can be seen from Figures 9 and 10 that the jump is more visible in some coefficients. For instance both in Figure 9 and in Figure 10 the highest value of the amplitude of wavelet coefficients for x(t) is in β01 while for y(t) is in β00. Probably this difference is due to the fact that x(t) is less regular than y(t).Wavelet coefficients of (3.1) with γ=102, β=2, and α=-0.01 in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)(c)(d)(e)(f)Wavelet coefficients of (3.1) with γ=102, β=2, and α=2 in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)(c)(d)(e)(f)


## 4.1. Critical Analysis through the Wavelet Coefficients

Let us take the short Haar wavelet transform for the two time series, obtained by using the Runge-Kutta4th order method, with the accuracy 10-6, of system  (3.1) with the following values of parameters-initial conditions: (4.6)γ=102,β=2,x0=0,y0=0,t∈[0,8].
The two time series correspond to the two values of α:α=-0.01 and α=2. It can be seen from Figures 9 and 10 that the jump is more visible in some coefficients. For instance both in Figure 9 and in Figure 10 the highest value of the amplitude of wavelet coefficients for x(t) is in β01 while for y(t) is in β00. Probably this difference is due to the fact that x(t) is less regular than y(t).Wavelet coefficients of (3.1) with γ=102, β=2, and α=-0.01 in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)(c)(d)(e)(f)Wavelet coefficients of (3.1) with γ=102, β=2, and α=2 in the interval t∈[0,8]; the initial conditions are x0=0, y0=0.
(a)(b)(c)(d)(e)(f)


## 5. Conclusion

In this paper the FitzHugh-Nagumo model has been considered with a high-amplitude pulse. The analysis was done by using wavelet coefficients and it has been shown that the dynamical system shows some kind of inertia against the rapid jumps. In fact, the jumps are detected with some delay with respect to the time they appear.



---

*Source: 895785-2010-06-15.xml*