Numerical Computation for
MRE (Magnetic Resonance Elasticity) by
Applying Numerical Differentiation Method
Kazuaki Nakane
Osaka University Graduate School of Medicine, Division of Health Sciences, Japan
Abstract. Palpation is a standard medical method to detect of abnormalities of
the human body, because a pathological state of soft tissues is often correlated
with changes in stiffness. However, a pathological lesion may be undetectable
by palpation if it is located deep in the body or if it is too small. Recently, new
technics are developing to identify stiffness (Young’s modulus). The observed
data contain noise, it prevent to substitute the data for the mathematical model
equation of the human body. Because the model equation includes the second
differential terms. “Numerical differentiation” is a numerical method to deter-
mine the derivatives of an unknown function from the given noisy values of the
unknown function at the scattered points. In this talk, we will apply this method
to the noisy data, and introduce numerical results.
1 Introduction
Palpation is a standard medical method to detect of abnormalities of the human body,
because a pathological state of soft tissues is often correlated with changes in stiffness.
However, a pathological lesion may be undetectable by palpation if it is located
deep in the body or if it is too small. Recently, new technics are developing to identify
stiffness(Young’s modulus). The principle of them is that:“By giving vibration from
outside of the body, we observethe propagationof the wave. By processing the observed
data, the stiffness can be identified.”
To observe the propagation of the wave, there are two ways. One is by using ul-
tarsonic device, the other is MRI device. For either technique, it is usually derived by
substituting the observed data for a mathematical model expressing the human body.
The human body can be approximated by the elastic model, we have to calculate the
second derivative of the data.
Because the observed data contain the noise, it is impossible to derivate them, di-
rectly. So far, by modifying data, for example taking average of neighbourhood, we
get the Young’s modulus. As for measuring a value of elasticity or the identification of
the border of the pathological part, we can not get correct information by this method.
We have to develop the method to calculate the derivation of the noisy data without
modifying of data.
“Numerical differentiation” is a numerical method to determine the derivatives of an
unknown function from the given noisy values of the unknown function at the scattered
points. at the scattered points.
Nakane K. (2008).
Numerical Computation for MRE (Magnetic Resonance Elasticity) by Applying Numerical Differentiation Method.
In Proceedings of the 2nd International Workshop on e-Health Services and Technologies, pages 73-77
DOI: 10.5220/0001896600730077
Copyright
c
 SciTePress

The higher order for one-dimensional and the first order for two dimensional nu-
merical differntiaton along the line of this method were given in [2] and [3], respec-
tively. For the two dimensional case, the new ingredient was that the variational prob-
lem for the regularised minimisation problem is solved by using Green’s function for
the Laplacian with Dirichlet boundary condition and a scheme for computing the first
order derivative was given in [3]. The numerical example showed that this method was
efficient. But in many applications, it is necessary to compute higher order derivatives.
In this talk, we will treat a problem concerned with MRE(Magnetic Resonance Elas-
ticity). It give an example to apply a numerical differntiation for the second derivatives
from noisy scattered data.
Mathematical setting and a mathematical result are introduced in the following.
2 Problem
Suppose we know Ω ⊂ R
2
is a bounded domain with piecewise C
2
boundary and
ρ = ρ (x) ∈ H
4
(Ω) is a function defined in Ω. Let N be a natural number and {x
i
}
N
i=1
be a group of points in Ω. We assume that Ω is divided into N parts {Ω
i
}
N
i=1
, and there
is only one points of {x
i
}
N
i=1
in each part. For simplicity we also assume that the areas
|Ω
i
| of all Ω
i
(1 ≤ i ≤ N ) are the same. We denote by d
i
the diameter of Ω
i
and let
d = max{d
i
}.
We will discuss the following problem:
Suppose that we know the approximate value ˜ρ
i
at point x
i
i.e.
| ˜ρ
i
− ρ(x
i
)| ≤ δ, i = 1, 2, · · · , N,
where δ > 0 is a given constant called the error level.
We want to find a function f
∗
(x) which approximates function ρ(x) such that
lim
d→0,δ→0
kf
∗
− ρk
H
2
(Ω)
= 0.
We treat this problem as the following optimisation problem by using Tikhonov
regularisation method.
Problem 2.1. Define a cost function Φ(f )
Φ(f) =
1
N
N
X
j=1
(f(x
i
) − ˜ρ
j
)
2
+ αk∆
2
fk
2
L
2
(Ω)
, f ∈ H
where H = {f : f ∈ H
4
(Ω), f |
∂Ω
= ∆f|
∂Ω
= 0}, and α > 0 is a regularisation
parameter. Then, the problem is to find f
∗
∈ H such that Φ(f
∗
) ≤ Φ(f ) for every
f ∈ H.
The existence and uniqueness of the minimiser of Problem is shown in [1].
Theorem 2.1. Problem 2.1 is equivalent to finding a unique solution f
∗
∈ H for the
following variational problem:
Z
Ω
∆
2
f
∗
∆
2
hdx = −
1
αN
N
X
j=1
(f(x
i
) − ˜ρ
j
)h(x
j
) (1)
74

for all h ∈ H. This equation is the Euler equation of Problem 2.1. Moreover, the min-
imiser of Problem 2.1 is unique.
To solve numerical differentiation problem, it is necessary to provide a scheme for
construction f
∗
. For that, by a formal argument using Green’s function of bi-harmonic
operator, we derive a method how to construct f
∗
. It will be shown as a theorem that
the constructed f
∗
by this method is the solution of (1).
We get the equation (2) by (1) (cf. [1])
αf
∗
(x) +
1
N
N
X
j=1
(f(x
i
) − ˜ρ
j
)
Z
Ω
G(x
j
, y)G(y, x)dy
= 0, (2)
where G(x, y) is a bi-harmonic Green’s function. By defining
a
j
(x) =
Z
Ω
G(x
j
, y)G(y, x)dy
and
c
j
= −
1
αN
(f
∗
(x
j
) − ˜ρ
j
),
(2) becomes
f
∗
(x) =
N
X
j=1
c
j
a
j
(x).
Now the problem of constructing f
∗
reduces to computing the coefficient c
j
from ˜ρ
j
.
From the definition of a
j
and c
j
with x = x
j
(j = 1, 2, · · · N), we obtain
c
j
= −
1
αN
(f(x
j
) − ˜ρ
j
)
= −
1
αN
(
N
X
k=1
a
k
(x
j
)c
k
− ˜ρ
j
). (3)
Let A be a (N, N ) matrix which is defined by
A =

αNδ
ij
+ a
i
(x
j
)

,
where δ
ij
is Kroneker’s delta. And let c and b be vectors
c = (c
j
), b = (b
j
).
Then (1) becomes the linear equations
Ac = b.
Solving this equations, we will obtain coefficients c
j
. Then we get f
∗
. To make sure f
∗
is the solution of our problem, we give the following theorem.
75

Theorem 2.2. Suppose function
f
∗
(x) =
N
X
k=1
a
k
(x)c
k
where {c
j
}
N
j=1
is the solution of linear system (3), then f
∗
is the solution of Problem
2.1.
3 Numerical Results
Here, we will see our method, the numerical differentiation, is effective. We insert a
noise to the numerical data that are given by numerical computation of direct problem.
Figure 1 is comparison with the elastic coefficient(stiffness) that we set in the direct
problem and the numerical result of our method. Reconstruction is performed with good
precision.
0
2e-07
4e-07
6e-07
8e-07
1e-06
1.2e-06
1.4e-06
1.6e-06
0 5e-09 1e-08 1.5e-08 2e-08
’muline’
72*x
Fig.1. Without noise.
Figure 2, we insert the noise(1%). For comparison, we show the result which is given
by ordinary differential method Figure 3.
0
2e-07
4e-07
6e-07
8e-07
1e-06
1.2e-06
1.4e-06
1.6e-06
0 5e-09 1e-08 1.5e-08 2e-08
’muline’
72*x
Fig.2. Numerical differential method with 1% noise.
Figure4, 10 noise is inserted to the data. We can reconstruct the elastic modulus. We
can see our method is robust to the noise.
76

0
2e-07
4e-07
6e-07
8e-07
1e-06
1.2e-06
1.4e-06
1.6e-06
0 5e-09 1e-08 1.5e-08 2e-08
’mulineexa’
72*x
Fig.3. Ordinary differential method with 1% noise.
0
2e-07
4e-07
6e-07
8e-07
1e-06
1.2e-06
1.4e-06
1.6e-06
0 5e-09 1e-08 1.5e-08 2e-08
’muline’
72*x
Fig.4. Numerical differential method with 10% noise.
4 Conclusions
The real observed data contain the noise. It is impossible to derive the stiffness(Young’s
modulus), directly. Usually, by modifying the data, it is given. As for measuring a cor-
rect value of elasticity or the identification of the border of the pathological part, we
have to develop the method to calculate the differential coefficient of the noisy data
without modifying of data.
In this note, we proposed the numerical differential method and introduced numeri-
cal results. We can see our method is robust to the noise. However, we do not know the
noise level of the real observed data. It is necessary to devise how to choose parameters
with the real data.
References
1. G. Nakamura, S. Wang and Y. Wang, Numerical differentiation for the second order deriva-
tives of functions with two variables, preprint.
2. Y. B. Wang, Y. C. Hon and J. Cheng, Reconstruction of High Order Derivatives from Input
Data, Submitted to Journal of Inverse and Ill-posed Problems(2004).
3. Y. B. Wang and T. Wei, Numerical Differentiation for two-dimensional scattered data, Sub-
mitted to Journal of Mathematical Analysis and Applications(2004).
77