Influence of Poisson’s ratio on elastographic direct and inverse problems J Fehrenbach Laboratoire MIP, Universit´e Paul Sabatier, 31062 TOULOUSE CEDEX 09, France E-mail:
[email protected] Abstract. We consider the displacement of an elastic material under an external compression (axial or almost axial stress). We assume that only one component of the displacement is observed, in the direction of compression (axial displacement), or alternatively, that two components are observed in a plane. These hypotheses are in accordance with an imaging modality, namely ultrasonic elastography. In the case of a homogeneous medium we show that any value of Poisson’s ratio allows to predict the observed value of the axial displacement. When two components of the displacement are measured in a plane, the Poisson’s ratio of the plane strain model that predicts the observed displacement is not the same as the tri-dimensional material. These facts are illustrated by numerical experiments in the case of an inhomogeneous medium. We also present results on experimental phantom data, where the inverse problem of reconstructing the Young’s modulus is solved assuming different values for Poisson’s ratio.
Keywords: elastography, Poisson’s ratio, Young’s modulus, inverse problem Submitted to: Phys. Med. Biol.
1. Introduction Ultrasonic elastography is a medical imaging modality that provides the directional displacement of a tissue under motion (Ophir et al 2002). The motion can either be a quasi-static motion (Ophir et al 1991), or a time-dependent motion (Gao et al 1995, Tanter et al 2002). Sonograms are acquired at time intervals and the displacement of the tissue is estimated by cross-correlation algorithms (Ophir et al 1991). When a single transducer is used, the displacement in the direction of the transducer (axial displacement) is known. When an array of transducers is used, the displacement along the array is also available (lateral direction), see (Konofagou and Ophir 1998). The displacement in the orthogonal direction (elevational motion) cannot be estimated at present time. The clinical applications of elastography rely on the fact that elastic properties of tissues provide medical information: certain kinds of cancers tend to be stiffer than the surrounding tissue (Garra et al 1997, Krouskop et al 1998, Souchon et al 2003),
Influence of Poisson’s ratio
2
deep vein thrombosis tends to be stiffer when it ages (Emelianov et al 2002). Given the motion of the tissue, a rough estimate of the stiffness is the inverse of the strain (Ponnekanti et al 1995), and a more precise estimate requires the use of inverse elastic models (Kallel and Bertrand 1996). We address here quasi-static elastography, where the external compression is usually in the axial direction. Most biological tissues are incompressible, hence the three components of the displacement u are related by: ∇ · u = 0. On the other hand, the measured displacements are partial displacements (we know only one or two components of u). The measures are usually available in a plane, and a plane strain model (linear isotropic elasticity) is often used to solve the inverse elasticity problem. However, the plane strain assumption requires that the medium is invariant and infinite along the elevational direction, or equivalently that it is invariant and clamped along the elevational direction. In (Konofagou and Ophir 1998), the medium is confined in the elevational direction and the two components of the displacement are used to estimate Poisson’s ratio. Unfortunately the elevational confinement of the medium is not always easy to achieve. The questions we address here are: when the material is not confined in the elevational direction, how does the axial displacement depend on Poisson’s ratio ? What brings the information that the tissue is incompressible ? Does it help to use incompressible plane strain models to identify the properties of the tissue (shear modulus) in the case the plane strain assumptions are not met ? Several works address the inverse problem in static elastography. In (Kallel and Bertrand 1996) the authors use a 2D simulation with a Poisson ratio of 0.495 and they can identify inclusions in the domain. In (Doyley et al 2000), the simulations and inverse problem on experimental data use a plane strain model with a Poisson ratio of 0.495, the phantom is constrained not to move in the lateral faces (it is not the elevational but the lateral faces that are clamped, this is far from the plane strain assumptions), and accurate values of the contrast are recovered. In (Oberai et al 2004), 2D and 3D simulations are performed with ν ∼ 0.499 9995, and in the phantom experiment the lateral faces are constrained, in (Doyley et al 2005) the simulations are performed in 2D with ν = 0.495 and during the in-vitro experiment, the faces of the phantom seem to be free, in (Fehrenbach et al 2006) Poisson’s ratio is 0.45 in the plane strain simulation and in the in-vitro experiment, where the elevational faces of the phantom are not constrained. All these studies report interesting results on the point of view of identification of the Young’s modulus. It is surprising that a model (plane strain) is applied in a situation where the hypotheses of this model are not met, and however yields valuable results. Our contribution here is to explain why this abuse of the plane strain model is not serious in the case of axial measurements. We also show what care should be taken when processing displacement measurements in both the axial and lateral directions. We assume in this work that the tissue is not constrained in the elevational direction, or equivalently that the stress is along the axial direction. We tackle two situations: the first one is when only axial displacements are known. We show analytically (on a
3
Influence of Poisson’s ratio
simple model) and numerically that the axial component of the solution of the direct problem does not depend much on Poisson’s ratio. The inverse problem is also studied on experimental data: the shear modulus distribution estimated by the inverse problem does not depend much on Poisson’s ratio. The second situation is when both axial and lateral displacements are known. In this case, a plane strain model is often used to predict the displacement. However, the plane strain model assumes that there is no displacement in the elevational direction, this hypothesis is not fulfilled since the elevational faces are left free. We study the effect of elevational motion: the Poisson’s ratio of the plane strain model is not the same as the tri-dimensional tissue in order to predict the same displacement. We study analytically a simple case, and numerically a more realistic case. This paper is organized as follows: in section 2 we present analytical results for the forward problem. In section 3 a numerical experiment is performed with a non-homogeneous phantom and the displacements are compared with displacements predicted by plane strain models. In section 4 we present the results of the inverse problem applied to experimental data with several values for Poisson’s ratio. 2. Analytical results 2.1. Direct problem We model an experiment as follows: the domain Ω is a unit cube filled with an elastic material, the lower horizontal face is perfectly lubricated and the upper horizontal face is perfectly lubricated and moved downwards of a small amount η, see figure 1. The four other faces are free. In order to ensure uniqueness of the resulting displacement, we also prevent horizontal translations and rotations.
imaging plane x1
x2
x3
Figure 1. The model experiment.
The x1 axis is the axial direction, the x2 axis is the lateral direction and the x3 axis is the elevational direction. The axial faces of the cube are the faces orthogonal to the x1 axis (top and bottom face on figure 1), the lateral faces of the cube are the faces orthogonal to the x2 axis (right and left face on figure 1), and the elevational faces are the faces orthogonal to the x3 axis (front and rear face on figure 1). The points of
Influence of Poisson’s ratio
4
the cube are subject to a displacement u = (u1 , u2 , u3 ) and a stress distribution σ that satisfy the following boundary conditions:
u1 (0, x2 , x3 ) = 0, u1 (1, x2 , x3 ) = −η, σn = 0 in the lateral and elevational faces, σn is vertical in the axial faces.
(1)
Since no volume forces are considered, the equilibrium equations write: ∇ · σ = 0 in Ω.
(2)
Hooke’s law for a linear isotropic medium describes the relation between the displacement and the stress fields: E Eν ∇·uI + ε(u), (3) σ= (1 + ν)(1 − 2ν) 1+ν where I is the identity matrix, ε(u) is the linearized strain tensor: 1 ε(u) = (∇u + ∇uT ), (4) 2 the space-dependent coefficient E is the Young’s modulus and Poisson’s ratio ν is assumed to be constant in the domain. The case of an incompressible material is a limit case of our model, since it is obtained by letting ν → 0.5. On the numerical point of view, instabilities appear as ν → 0.5. However, it is shown in (Guillaume et al 2002) that the solution depends analytically on Poisson’s ratio, hence it is continuous as ν → 0.5. The incompressible solution is expressed in (Guillaume et al 2002) as the sum of a power series of compressible (hence numerically stable) solutions. 2.2. Homogeneous cube If the material is homogeneous (the Young’s modulus E is constant), the solution of (1)-(3) is ηx1 u= −νηx2 , −νηx3
this expression is still valid in the incompressible case. The first component of the displacement does not depend on Poisson’s ratio ν, in other words: all the values of ν ∈ [0, 0.5] give rise to the same axial displacement. 2.3. Inverse problem with a plane image
The imaging region is the plane of equation x3 = 0, see figure 1. We now assume that in this plane the two components of the displacement are known: the observed displacement is ! ηx1 obs u = . −νηx2
5
Influence of Poisson’s ratio
Let us consider a plane strain model in a square domain, for a homogeneous medium with Poisson’s ratio ν 0 . Under the prescribed boundary conditions, the displacement is u
plane
=
ηx1 ν0 − 1−ν 0 ηx2
!
.
The value of ν 0 that predicts the same displacement as uobs in the plane strain approximation is defined by ν ν0 = . 1+ν The difference between these values comes from the fact that in the model experiment, the elevational faces are free to move, and in the plane strain model they are clamped (or the material is infinite in the elevational direction, which amounts the same). For instance, in the case of an incompressible material, ν = 0.5 and ν 0 = 1/3. We can say that when an incompressible 3D-medium is observed in 2D and the elevational faces are left free, the observed displacement is not the displacement predicted by an incompressible plane strain model. It coincides with the displacement predicted by a plane strain model with a Poisson’s ratio of 1/3. The difference is observed on the lateral component of the displacement (along x2 ). 2.4. Case of an almost axial stress We consider here the case of an elastic material subject to an almost axial stress: we assume that the stress tensor is σ = αe1 ⊗ e1 + o(α), where (e1 , e2 , e3 ) is the orthonormal basis associated to the axes x1 , x2 , x3 and o(α) is negligible compared to α. This is likely to occur in situations analogous to the one described above, namely in case of axial compression, when the other faces of the material are left free, or when the other faces are subject to pressure forces that are negligible compared to the axial compression. The Young’s modulus distribution is E, and we study the effect of changing Poisson’s ratio on the solution. Let us consider two values of Poisson’s ratio ν0 and ν, we denote ε(ν0 ) and ε(ν) the strain tensors of the corresponding solutions. Hooke’s law 3 yields: Eε(ν0 ) = (1 + ν0 )σ − ν0 tr σI, and Eε(ν) = (1 + ν)σ − ν tr σI. The difference of these equations reads: E(ε(ν0 ) − ε(ν)) = (ν0 − ν)σ − (ν0 − ν) tr σI = −α(ν0 − ν)(e2 ⊗ e2 + e3 ⊗ e3 ) + o(α).
6
Influence of Poisson’s ratio
This shows that Poisson’s ratio has a negligible influence on the axial component of the strain, hence it has a negligible influence on the axial displacement. 3. Numerical results We have shown analytically that the axial displacement under an axial compression (other faces free) does not depend on Poisson’s ratio. This was shown on a homogeneous cubic phantom, and in the case of an axial stress. We have also shown that when observing a displacement in a plane, a 3D incompressible medium undergoes a 2D displacement that does not come from an incompressible medium, rather from a medium with Poisson’s ratio 1/3. The general case of a material with space-dependent Young’s modulus is hard to tackle analytically. In this section, we perform numerical experiments that tend to prove that the same phenomena still hold. 3.1. Dependence of the axial component on Poisson’s ratio We consider a cubic phantom with non-homogeneous elastic material. Poisson’s ratio is constant equal to ν, and Young’s modulus depends on the space variable as follows: the phantom comprises a hard cylindrical inclusion where Young’s modulus is 3 times higher that the background, and a hard ball where it is 10 times higher, see figure 2. The boundary conditions are still given by (1).
x1 x2
x3
Figure 2. The inhomogeneities for the numerical experiment: the cylinder is 3 times stiffer that the background, and the ball 10 times stiffer
A P1 finite element method was implemented on a 1.5 MHz computer with GetFem++ (www-gmm.insa-toulouse.fr/getfem/). The cube was discretized with 8339 points and 42709 convexes. The displacement u(ν) = (u1 (ν), u2 (ν), u3 (ν)) was computed for several values of the Poisson’s ratio: ν = 0.2, 0.35, 0.49, 0.4995. We show on figure 3 the axial, lateral and elevational displacements for ν = 0.35, and on figure 4 each component of the difference u(ν = 0.35) − u(ν = 0.4995). The axial component has not a significant variation, while the other components have a significant variation. We show in table 1, for each component (axial, lateral and elevational), the relative difference of the displacements between ν0 = 0.35 and each other value of ν.
7
Influence of Poisson’s ratio
Figure 3. The three components of the displacement for ν = 0.35: axial (left), lateral (center) and elevational (right).
Figure 4. The difference of the displacements between ν = 0.35 and ν = 0.4995: axial (left), lateral (center) and elevational (right) component Table 1. Relative difference ||ui (ν) − ui (ν0 )||/||ui (ν0 )|| between the displacement obtained for ν and for ν0 = 0.35, for i = 1 (axial component), i = 2 (lateral component) and i = 3 (elevational component).
ν = 0.2 ν = 0.49 ν = 0.4995
axial
lateral
elevational
0.005 0.01 0.03
0.412 0.413 0.425
0.433 0.432 0.478
We observe from table 1 and figure 4 that the axial displacement does not vary significantly as Poisson’s ratio varies in a wide range [0.2, 0.4995]. On the contrary, the lateral and elevational displacements depend heavily on Poisson’s ratio. This confirms the results of sections 2.2 and 2.4. 3.2. Influence of Poisson’s ratio on the lateral displacement In this section we illustrate the fact that the lateral displacement of a plane strain model depends on Poisson’s ratio, and we show which value of Poisson’s ratio predicts a displacement close to the displacement observed in the imaging plane. We use a 3-dimensional numerical phantom similar to the one presented in figure 2, except the fact that it has only one hard cylindrical inclusion. The material is assumed to be almost incompressible (ν = 0.4995), and it is submitted to the same boundary conditions (1). We are given two components of the displacement in the imaging plane
8
Influence of Poisson’s ratio Π of equation x3 = 0, see figure 5: uobs = (u1 , u2 ) in Π.
We want to find a plane material that gives rise to the same (plane) displacement field under the same boundary conditions, assuming a plane strain model, see figure 6. The shear modulus of the material is assumed to be known (it is the shear modulus of the 3D material in the imaging plane), and we try several values for Poisson’s ratio: ν = 0.2, 0.3, 0.32, 0.33, 0.35, 0.37, 0.4, 0.45 and 0.49. For each of these values, we compute uplane (ν) = (uplane (ν), uplane (ν)) the displacement predicted by a plane strain 1 2 model.
Figure 5. Axial (left) and lateral (right) displacements for the 3-dimensional medium in the imaging plane Π.
x1 ν x2
Figure 6. Left: the plane model: the shear modulus is known and Poisson’s ratio ν is a parameter ; center: axial displacement for the plane strain model for ν = 0.33 ; right: lateral displacement for the plane strain model for ν = 0.33.
In order to find the value of Poisson’s ratio that predicts the observed displacement u , we show on figure 7 the relative difference between the observed axial (resp. lateral) displacement, and the axial (resp. lateral) displacement predicted by the plane strain model: ||uobs − uplane (ν)||/||uobs i i i ||. The following observation can be made: the plane strain model predicts an accurate axial component for any value of ν (the relative difference is of order of a few percents), and it predicts an acceptable lateral component for values of ν between 0.3 and 0.35 (the relative difference is of order 10-15%). This confirms the result of section 2.3. obs
9
Influence of Poisson’s ratio 1 0.9
relative difference
0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0.2
0.25
0.3 0.35 0.4 0.45 Poisson ratio ν of the plane model
0.5
Figure 7. The relative difference of the axial (circles) and lateral (squares) components between the observations uobs in the imaging plane and the plane strain model uplane i i for different values of Poisson’s ratio ν.
4. Experimental results 4.1. Methods An experiment was performed by A. Kumar (Houston) on a tissue-mimicking phantom made of gelatin-agar-water mixtures. The phantom was made as a homogeneous block with a cylindrical inclusion running through the center. The background was prepared by mixing 5% by weight gelatin and 3% by weight agar in de-ionized water at 80o C (Kallel et al 2001). The cylindrical inclusion was made from reinforced sponge material (Srinivasan et al 2004). Sample pieces of the background and the inclusion were also made simultaneously when the phantoms were prepared. These samples were tested in a Nano-Indenter XP platform (MTS Systems Corporation, Oak Ridge, TN) to determine their respective Young’s moduli. The procedure to use the Nano-Indenter to determine the modulus is described in detail in (Srinivasan et al 2004). An HDI 1000 (ATL-Philips Inc, Bothell, WA) US scanner, with a 128-element array, 5 MHz center frequency and a 60% fractional bandwidth was used for data acquisition, the sampling rate was 20 MHz. The transducer was attached to a compressor plate. The axial displacement was estimated through cross-correlation of the pre- and post-compression RF A-lines. The inverse problem was implemented to estimate the relative Young’s modulus from the axial displacement estimation. The implementation of this inverse problem is described in (Fehrenbach et al 2006). A plane strain model is used, for several values of Poisson’s ratio : ν = 0.32, 0.37, 0.42. 4.2. Results The estimations of the Young’s modulus contrast are shown on figure 8. The contrast of the inclusion relatively to the background does not vary significantly for these 3 values of ν, it is between 2 and 2.5 (this is to be compared with the results obtained by the Nano-Indenter that gives a contrast of 2.1±0.4). This is in accordance with the results of section 2.4.
Influence of Poisson’s ratio
10
Figure 8. The Young’s modulus contrast estimated with a plane strain model for different values of ν: ν = 0.32 (left), ν = 0.37 (center) and ν = 0.42 (right)
5. Discussion We have studied the effect of Poisson’s ratio ν on the axial displacement in the case of an axial (or almost axial) stress. The analytical study and the numerical results show that Poisson’s ratio has a very small influence on the axial displacement (direct problem). It is therefore possible to choose any plane strain model (incompressible, almost incompressible or compressible) to study the inverse problem. This explains why inverse problems using only axial displacements give reasonnable answers whatever value is chosen for ν. However, for numerical stability reasons, it is not desirable to use almost incompressible models (locking phenomena may appear). When bidirectional displacements are measured (axial and lateral displacement), the following was observed: when the lateral and elevational faces are free, or when the stress is almost axial, the plane strain model that is close to what is observed in ν the imaging plane has a Poissson’s ratio of ν 0 = if the 3-dimensional material 1+ν has a Poisson’s ratio of ν. When the inverse problem is tackled with a plane strain approximation, it should not be adressed with a Poisson’s ratio of ν, rather with the apparent value ν 0 . In the case of an incompressible medium, the apparent Poisson’s ratio (for the equivalent plane strain approximation) is 1/3. Acknowledgments This work was supported in part by ACINIM-IMEG CNRS project (France). The author is grateful to Arun Kumar and Jonathan Ophir (University of Texas, Medical School at Houston) for the experimental data. References Doyley M, Meaney P and Bamber J 2000 Evaluation of an iterative reconstruction method for quantitative elastography Phys. Med. Biol. 45 1521–40 Doyley M, Srinivasan S, Pendergrass SA, Wu Z and Ophir J 2005 Comparative evaluation of strainbased and model-based modulus elastography Ultrasound Med. Biol. 31(6) 787–802 Emelianov S, Chan X, O’Donnell M, Knipp B, Myers D, Wakefield D and Rubin J 2002 Triplex ultrasound: elasticity imaging to age deep vein thrombosis Ultrasound Med. Biol. 28 757–67
Influence of Poisson’s ratio
11
Fehrenbach J, Masmoudi M, Souchon R and Trompette P 2006 Detection of small inclusions by elastography Inverse Probl. 22 1055–69 Gao L, Parker K and Alam S 1995 Sonoelasticity imaging: theory and experimental verification J. Acoust. Soc. Am. 97 3875–86 Garra B, Cespedes E, Ophir J, Spratt S, Zuurbier R, Magnant C and Pennanen M 1997 Elastography of breast lesions: initial clinical results Radiology 202 79–86 Guillaume P, Masmoudi M and Zeglaoui A 2002 From compressible to incompressible materials via an asymptotic expansion Numerische Mathematik 91 649–73 Kallel F and Bertrand M 1996 Tissue elasticity reconstruction using linear perturbation method IEEE T. Med. Imaging 15(3) 299–313 Kallel F, Prihoda CD and Ophir J 2001 Contrast-transfer efficiency for continuously varying tissue moduli: simulation and phantom validation Ultrasound Med. Biol. 27 1115–25 Konofagou E and Ophir J 1998 A new elastographic method for estimation and imaging of lateral displacements, lateral strains, corrected axial strains and Poisson’s ratios in tissues Ultrasound Med. Biol. 24(8) 1183–99 Krouskop T, Wheeler R, Kallel F, Garra B and Hall T 1998 Elastic moduli of breast and prostate tissues under compression Ultrasonic Imaging 20 260–74 Oberai A, Gokhale N, Doyley M and Bamber J 2004 Evaluation of the adjoint equation based algorithm for elasticity imaging Phys. Med. Biol. 49 2955–74 Ophir J, Cespedes I, Ponnekanti H, Yazdi Y and Li X 1991 Elastography: a quantitative method for imaging the elasticity of biological tissues Ultrasonic imaging 13 111–34 Ophir J, Alam S, Garra B, Kallel F, Konofagou E, Krouskop T, Merritt C, Righetti R, Souchon R, Srinivasan S and Varghese T 2002 Elastography: imaging the elastic properties of soft tissues with ultrasound J. Med. Ultrasonics 29 155–71 Ponnekanti H, Ophir J, Huang Y and Cespedes I 1995 Fundamental mechanical limitations on the visualization of elasticity contrast in elastography Ultrasound Med. Biol. 21 533–43 Srinivasan S, Krouskop T and Ophir J 2004 Comparing elastographic strain images with modulus images obtained using nanoindentation: preliminary results using phantoms and tissue samples Ultrasound Med. Biol. 30(3) 329–43 Souchon R, Rouviere O, Gelet A, Detti V, Srinivasan S, Ophir J and Chapelon JY 2003 Visualisation of HIFU lesions using elastography of the human prostate in vivo: preliminary results Ultrasound Med. Biol. 29 (7) 1007–15 Tanter M, Bercoff J, Sandrin L and Fink M 2002 Ultrafast compound imaging for 2D-motion vector estimation: application to transient elastography IEEE T. Ultrason. Ferroelectr. Freq. Control 49 1363–74.