Direct computation of myocardial deformation and strain from
tagged cine MRI
Citation for published version (APA):
Florack, L. M. J., & Assen, van, H. C. (2008). Direct computation of myocardial deformation and strain from
tagged cine MRI. (CASA-report; Vol. 0822). Technische Universiteit Eindhoven.
Document status and date:
Published: 01/01/2008
Document Version:
Publisher’s PDF, also known as Version of Record (includes final page, issue and volume numbers)
Please check the document version of this publication:
• A submitted manuscript is the version of the article upon submission and before peer-review. There can be
important differences between the submitted version and the official published version of record. People
interested in the research are advised to contact the author for the final version of the publication, or visit the
DOI to the publisher's website.
• The final author version and the galley proof are versions of the publication after peer review.
• The final published version features the final layout of the paper including the volume, issue and page
numbers.
Link to publication
General rights
Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright owners and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights. • Users may download and print one copy of any publication from the public portal for the purpose of private study or research. • You may not further distribute the material or use it for any profit-making activity or commercial gain
• You may freely distribute the URL identifying the publication in the public portal.
If the publication is distributed under the terms of Article 25fa of the Dutch Copyright Act, indicated by the “Taverne” license above, please follow below link for the End User Agreement:
www.tue.nl/taverne
Take down policy
If you believe that this document breaches copyright please contact us at:
openaccess@tue.nl
providing details and we will investigate your claim.
and Strain from Tagged Cine MRI
Luc Florack1,2 and Hans van Assen2
1 Department of Mathematics and Computer Science, Eindhoven University of
Technology, The Netherlands, L.M.J.Florack@tue.nl
2 Department of Biomedical Engineering, Eindhoven University of Technology,
The Netherlands, H.C.v.Assen@tue.nl
Abstract. Myocardial deformation and strain can be investigated us-ing suitably encoded cine MRI that admits disambiguation of material motion. Practical limitations restrict the analysis to in-plane motion in cross-sections of the heart (2D+time), but the proposed method readily generalizes to 3D+time. Time evolution of the deformation tensor is gov-erned by a first order ordinary differential equation, which is completely determined by the velocity gradient tensor. We solve this matrix-ODE analytically, and present results obtained from healthy volunteers. The proposed method is robust and requires only off-the-shelf algorithms.
1
Introduction
Cine MRI [1], combined with (C)SPAMM3 encoding technology [2–4], admits disambiguation of local tissue motion, thus enabling the extraction of myocardial deformation and strain [5]. This can be achieved either without explicit a priori regularization, or through exploitation of sparse constraints combined with in-terpolation and/or regularization [3, 6–31]. Possible encodings are DENSE (Dis-placement ENcoding with Stimulated Echoes), and HARP (HARmonic Phase). Tagging based methods using HARP technology form our point of departure. Given a dense motion field within the myocardium, our aim is to devise an op-erational procedure for direct extraction of myocardial deformation and strain. By “direct” we mean that we seek to obviate sophisticated preprocessing steps, such as segmentation of, or interpolation between tag lines, finite element meth-ods explicitly coupled to the tagging pattern, etc. Although such sophisticated “indirect” procedures exist and have been proven powerful, they require spe-cific algorithmics that is neither trivially implemented nor readily available. In addition we aim to minimize the number of extrinsic control parameters. We
2
believe that the parsimony of our method facilitates applicability and optimiza-tion, since only off-the-shelf algorithms (linear filtering and inversion of linear systems) are needed in our computation of myocardial deformation and strain.
2
Theory
2.1 Deformation
The velocity gradient tensor, with components Lα
βrelative to a coordinate frame,
relates the rate of change of a momentary infinitesimal line element d ˙xα to the
line element dxβitself. It can be shown that (employing summation convention)
d ˙xα= Lα
βdxβ with Lαβ =
∂vα
∂xβ (α, β = 1, . . . , n) . (1)
The numbers Lα
β constitute a matrix L with row and column index4 α, β. To
get an estimate of this tensor field we have applied the algorithm of Van Assen et al. and Florack et al. [9, 11], solving a linear system of algebraic equations. Tissue deformation can be described by a map f : M → N : x(X, t0) 7→ x(X, t), in which x(X, t) denotes the spatial position of a material point at time t, with reference position X = x(X, t0) at time t0(“Lagrange picture”). Domain M and codomain N are copies of the deformable tissue medium (a subset of IRn) at
times t0, respectively t ≥ t0. We fix t0so as to correspond to end-diastole. The associated differential map,
f∗≡ df : T MX −→ T Nx: v 7→ f∗(v) = viFiαfα, (2)
relative to a basis {fα} of T Nx, provides a local linearization of tissue
deforma-tion, the deformation tensor. Here T MXis the tangent space of M at the fiducial
point X = x(X, t0), and T Nxthat of N at the f -mapped point x = x(X, t). The
matrix F (t, t0) of this linear map relative to the local coordinate charts {xα, N }
and {Xi, M } is given by the Jacobian
Fα
i =
∂xα
∂Xi. (3)
By virtue of the chain rule, the relation between deformation and velocity gra-dient tensors, Eqs. (1) and (3), is given by the first order ODE [32]
˙
F = L F , (4)
4 Our analysis is confined to a single short-axis plane, whence n = 2, i.e. we only
subject to an initial condition, viz. F (t = t0, t0) = I. Equivalently, F (t, t0) = I +
Z t
t0
L(s) F (s, t0) ds . (5)
This induces an expansion known as the matricant, cf. Gantmacher [32]: F (t, t0) = I + Z t t0 L(τ ) dτ + Z t t0 L(τ ) Z τ t0 L(σ) dσ dτ + . . . (6) The matricant has the following property:
F (t, t0) = F (t, t1) F (t1, t0) (t0< t1< t) . (7) It follows that, if we split the interval [t0, t] into n parts by using intermediate points t1, . . . , tn−1separated by ∆tk = tk− tk−1(k = 1, . . . , n, with tn= t),
F (t, t0) = F (t, tn−1) . . . F (t1, t0) with t0< t1< . . . < tn−1< tn= t. (8)
For an infinitesimally narrow time interval [tk−1, tk] we have by approximation
F (tk, tk−1) = I + L(t∗k) ∆tk+ higher order terms in ∆tk, (9)
in which t∗
k is any point satisfying tk−1 ≤ t∗k ≤ tk. Eqs. (8–9) give rise to a
representation in terms of a so-called multiplicative integral [32]: F (t, t0) =∼ Z t t0 (I + L(τ ) dτ )def= lim ∆tk→0(I +L(t ∗ n) ∆tn). . .(I +L(t∗1) ∆t1) . (10) One recognizes the multiplicative counterpart of the Riemann sum approxima-tion for ordinary (“additive”) integrals. One can show that this is identical to
F (t, t0)[32]= ∼ Z t
t0
exp (L(τ ) dτ )def= lim
∆tk→0exp (L(t
∗
n) ∆tn). . .exp (L(t∗1) ∆t1) . (11) Several properties of the deformation tensor are manifest in this representation. For instance, since for square matrices A, B, one has (i) det AB = det A det B, (ii) det(I + ²A) = 1 + ²tr A + O(²2), and (iii) det exp A = exp tr A, it shows that
det F (t, t0) =∼ Z t t0 (1 + tr L(τ ) dτ ) =∼ Z t t0 exp (tr L(τ ) dτ ) . (12) In particular, a divergence free velocity field (div v = tr L = 0) preserves volumes: det F (t, t0) = 1. Furthermore, exp A exp B = exp(A + B) if [A, B] = 0, whence for a stationary velocity field (L(t) = L0pointwise constant) one finds F (t, t0) = exp ((t − t0)L0). Finally, the multiplicative integral suggests a straightforward numerical approximation, viz. by using either Eq. (10) or (11) without limiting procedure. (In this case the two representations are of course no longer identical.)
4
2.2 Strain
On the basis of the differential map f∗, Eq. (2), and its transpose f∗T, one defines
an intrinsic mapping on T MX, known as the Lagrangian strain tensor [33],
E ≡1 2(f
T
∗ ◦ f∗− idT MX) : T MX→ T MX: v 7→ E(v) = viEijej, (13)
relative to a basis {ej} of T MX, with mixed tensor components
Eij= 1 2 ³ g`jFα `hαβFiβ− δij ´ . (14)
Here gij are the components of the dual (Euclidean) metric tensor in domain
co-ordinates {Xi, M }, and h
αβ those pertaining to codomain coordinates {xα, N }.
Note that if f∗is an isometry, then E vanishes identically by construction. Thus
E captures nonrigid deformations only.
Below we consider a single Cartesian coordinate system for both domain and codomain. (In such a system the representations of the various metric tensors, gij and hαβ and their duals, gij and hαβ, all simplify to identity matrices.)
3
Experiment
Figure 1 illustrates various scalar fields extracted from the strain tensor field, Eqs. (13–14), by contraction with a pair of local unit vectors, viz. radial (r=radial) and azimuthal (c=circumferential) basis vectors of the polar coordinate system centered at the midpoint of the region of interest, and those defining the strain tensor’s eigensystem. If u, v ∈ T MX are two such unit vectors, then the local
scalar quantity derived from the strain tensor is given by the inner product (u, E(v)) = (E(u), v) = gjkukEijvi. Tensor components are evaluated in
Carte-sian coordinates. Figure 2 illustrates temporal evolutions of these quantities spatially averaged over the respective regions of interest, together with their standard deviations. Table 1 shows statistics for three healthy volunteers.
4
Conclusion
We have proposed a simple and robust linear model for extracting myocardial deformation and strain. Material motion and gradient velocity are determined by a multi-scale linear system of algebraic equations. We have analytically solved the linear matrix-ODE governing myocardial deformation. By discretizing the closed-form solution we have subsequently solved for the induced Lagrangian strain tensor field, yielding results that are typical for healthy volunteers, cf.
Garot et al. [12], thus demonstrating the feasibility of our method. An advantage of our method is that only off-the-shelf algorithms are needed. This serves clarity, facilitates optimization, and enables implementation on dedicated hardware. Of course, our quantitative results demand evaluation by a more comprehensive benchmarking on in vivo as well as carefully controlled synthetic or phantom data. Moreover, by virtue of the transparent mathematical framework, it is pos-sible to assess tolerances on the basis of theoretical error propagation models. Future investigation will be needed to pursue these recommendations, and com-pare theoretically predicted and experimentally observed performance measures. Acknowledgments. The Netherlands Organisation for Scientific Research (NWO) is gratefully acknowledged for financial support. Jos Westenberg has provided the MR data used in our experiments.
base mid apex
I II III I II III I II III
−0.16±0.13 −0.17±0.11 −0.19±0.11 −0.19±0.09 −0.17±0.08 −0.15±0.09 −0.16±0.16 −0.20±0.09 −0.19±0.13
0.20±0.25 0.16±0.21 0.14±0.25 0.04±0.20 0.14±0.15 0.03±0.15 0.07±0.23 0.20±0.23 −0.02±0.18 0.02±0.15 −0.04±0.12 −0.07±0.13 0.01±0.07 −0.02±0.09 −0.06±0.11 0.07±0.14 −0.03±0.11 −0.05±0.12 0.25±0.22 0.19±0.20 0.18±0.23 0.07±0.17 0.16±0.15 0.09±0.16 0.15±0.21 0.22±0.22 0.05±0.19
−0.22±0.11 −0.21±0.09 −0.24±0.09 −0.23±0.09 −0.19±0.08 −0.21±0.07 −0.26±0.12 −0.23±0.09 −0.27±0.09
Table 1. Average and standard deviation of strains from three healthy volunteers (labels I,II,III) over base, mid, and apex short-axis cross-sectional ROIs at that time frame at which Emaxhas attained its maximum in mid-slice (tI= 16, tII= 15, tIII= 14).
Rows correspond to, from top to bottom, Ecc, Err, Ecr, Emax, and Emin.
References
1. Zerhouni, E.A., Parish, D.M., Rogers, W.J., Yang, A., Shapiro, E.P.: Human heart: Tagging with MR imaging—a method for noninvasive assessment of myocardial motion. Radiology 169(1) (1988) 59–63
2. Axel, L., Dougherty, L.: MR imaging of motion with spatial modulation of mag-netization. Radiology 171(3) (1989) 841–845
3. Axel, L., Montillo, A., Kim, D.: Tagged magnetic resonance imaging of the heart: a survey. Medical Image Analysis 9(4) (2005) 376–393
4. Fischer, S.E., McKinnon, G.C., Maier, S.E., Boesiger, P.: Improved myocardial tagging contrast. Magnetic Resonance in Medicine 30(2) (1993) 191–200
5. Axel, L., Dougherty, L.: Heart wall motion—improved method of spatial modula-tion of magnetizamodula-tion for MR imaging. Radiology 172(2) (1989) 349–350
6. Aletras, A.H., Balaban, R.S., Wen, H.: High-resolution strain analysis of the human heart with fast-DENSE. Journal of Magnetic Resonance 140(1) (1999) 41–57
6
7. Aletras, A.H., Ding, S., Balaban, R.S., Wen, H.: DENSE: Displacement encoding with stimulated echoes in cardiac functional MRI. Journal of Magnetic Resonance 137(1) (1999) 247–252
8. Aletras, A.H., Wen, H.: Mixed echo train acquisition displacement encoding with stimulated echoes: An optimized DENSE method for in vivo functional imaging of the human heart. Magnetic Resonance in Medicine 46(3) (2001) 523–534 9. Assen, H.v., Florack, L., Suinesiaputra, A., Westenberg, J., ter Haar Romeny, B.:
Purely evidence based multiscale cardiac tracking using optic flow. In Miller, K., Paulsen, K.D., Young, A.A., Nielsen, P.M.F., eds.: Proceedings of the MICCAI Workshop on Computational Biomechanics for Medicine II (Brisbane, Australia, October 29, 2007). (2007) 84–93
10. Axel, L., Chen, T., Manglik, T.: Dense myocardium deformation estimation for 2D tagged MRI. In Frangi, A.F., Radeva, P.I., Santos, A., Hernandez, M., eds.: Proceedings of the Third International Conference on Functional Imaging and Modeling of the Heart—FIMH 2005 (Barcelona, Spain, June 2–4 2005). Volume 3504 of Lecture Notes in Computer Science., Berlin, Springer-Verlag (2005) 446– 456
11. Florack, L., van Assen, H., Suinesiaputra, A.: Dense multiscale motion extraction from cardiac cine MR tagging using HARP technology. In Niessen, W., Westin, C.F., Nielsen, M., eds.: Proceedings of the 8th IEEE Computer Society Workshop on Mathematical Methods in Biomedical Image Analysis, held in conjuction with the IEEE International Conference on Computer Vision (Rio de Janeiro, Brazil, October 14–20, 2007), Omnipress (2007) Digital proceedings.
12. Garot, J., Bluemke, D.A., Osman, N.F., Rochitte, C.E., McVeigh, E.R., Zerhouni, E.A., Prince, J.L., Lima, J.A.C.: Fast determination of regional myocardial strain fields from tagged cardiac images using harmonic phase MRI. Circulation 101(9) (March 2000) 981–988
13. Gilson, W.D., Yang, Z., French, B.A., Epstein, F.H.: Complementary displacement-encoded MRI for contrast-enhanced infarct detection and quantification of myocar-dial function in mice. Magnetic Resonance in Medicine 51(4) (2004) 744–752 14. Gupta, S.N., Prince, J.L.: On variable brightness optical flow for tagged MRI. In
Bizais, Y., ed.: Proceedings of the Fourteenth International Conference on Informa-tion Processing in Medical Imaging–IPMI 1995 (Ile de Berder, France), Dordrecht, Kluwer Academic Publishers (1995) 323–334
15. Haber, I., Kikinis, R., Westin, C.F.: Phase-driven finite element model for spatio-temporal tracking in cardiac tagged MRI. In: Proceedings of the Fourth International Conference on Medical Image Computing and Computer-Assisted Intervention—MICCAI 2001 (Utrecht, The Netherlands, October 14–17 2001). Vol-ume 2208 of Lecture Notes in Computer Science., Berlin, Springer-Verlag (2001) 1332–1335
16. Kraitchman, D.L., Sampath, S., Castillo, E., Derbyshire, J.A., Boston, R.C., Bluemke, D.A., Gerber, B.L., Prince, J.L., Osman, N.F.: Quantitative ischemia detection during cardiac magnetic resonance stress testing by use of FastHARP. Circulation 107 (March 2003) 2025–2030
17. Kuijer, J.P., Hofman, M.B., Zwanenburg, J.J., Marcus, J.T., van Rossum, A.C., Heethaar, R.M.: DENSE and HARP: Two views on the same technique of phase-based strain imaging. Journal of Magnetic Resonance Imaging 24(6) (2006) 1432– 1436
18. Osman, N.F., Kerwin, W.S., McVeigh, E.R., Prince, J.L.: Cardiac motion track-ing ustrack-ing CINE harmonic phase (HARP) magnetic resonance imagtrack-ing. Magnetic Resonance in Medicine 42(6) (1999) 1048–1060
19. Osman, N.F., McVeigh, E.R., Prince, J.L.: Imaging heart motion harmonic phase MRI. IEEE Transactions on Medical Imaging 19(3) (2000) 186–202
20. Ozturk, C., Derbyshire, J.A., McVeigh, E.R.: Estimating motion from MRI data. Proceedings of the IEEE 91(10) (2003) 1627–1648
21. Pan, L., Prince, J.L., Lima, J.A.C., Osman, N.F.: Fast tracking of cardiac motion using 3D-HARP. IEEE Transactions on Biomedical Engineering 52(8) (2005) 1425–1435
22. Prince, J.L., McVeigh, E.R.: Motion estimation from tagged MR image sequences. IEEE Transactions on Medical Imaging 11(2) (1992) 238–249
23. Ryf, S., Spiegel, M.A., Gerber, M., P., B.: Myocardial tagging with 3D–CSPAMM. Journal of Magnetic Resonance Imaging 16(3) (2002) 320–325
24. Sampath, S., Derbyshire, J.A., Atalar, E., Osman, N.F., Prince, J.L.: Real-time imaging of two-dimensional cardiac strain using a harmonic phase magnetic res-onance imaging (HARP-MRI) pulse sequence. Magnetic Resres-onance in Medicine 50(1) (2003) 154–163
25. Sampath, S., Derbyshire, J.A., Osman, N.F., Atalar, E., Prince, J.L.: Real-time imaging of cardiac strain using an ultra-fast HARP sequence. In: Proceedings of the 9th International Society for Magnetic Resonance in Medicine, ISMRM (Glasgow, Scotland, April 21–27 2001). (2001) 111
26. Spottiswoode, B.S., Zhong, X., Hess, A.T., Kramer, C.M., Meintjes, E.M., Mayosi, B.M., Epstein, F.H.: Tracking myocardial motion from cine DENSE images using spatiotemporal phase unwrapping and temporal fitting. IEEE Transactions on Medical Imaging 26(1) (2007) 15–30
27. Steenstrup-Pedersen, K., Nielsen, M.: Computing optic flow by scale-space inte-gration of normal flow. In Kerckhove, M., ed.: Scale-Space and Morphology in Computer Vision: Proceedings of the Third International Conference, Scale-Space 2001, Vancouver, Canada. Volume 2106 of Lecture Notes in Computer Science., Berlin, Springer-Verlag (July 2001) 14–25
28. Suinesiaputra, A., Florack, L.M.J., ter Haar Romeny, B.M.: Multiscale optic flow analysis of MR tagging heart image sequences. European Journal of Medical Physics 19 (2003) 48
29. Suinesiaputra, A., Florack, L.M.J., Westenberg, J.J.M., ter Haar Romeny, B.M., Reiber, J.H.C., Lelieveldt, B.P.F.: Optic flow computation from cardiac MR tag-ging using a multiscale differential method—a comparative study with velocity-encoded mri. In: Proceedings of the Sixth International Conference on Medical Image Computing and Computer-Assisted Intervention—MICCAI 2003 (Montr´eal, Canada, November 15–18 2003). Volume 2878–2879 of Lecture Notes in Computer Science., Berlin, Springer-Verlag (2003) 483–490
30. Young, A.A.: Model tags: Direct three-dimensional tracking of heart wall motion from tagged magnetic resonance images. Medical Image Analysis 3(4) (1999) 361– 372
31. Young, A.A., Axel, L.: Three-dimensional motion and deformation of the heart wall: Estimation with spatial modulation of magnetization–a model-based ap-proach. Radiology 185(1) (1992) 241–247
32. Gantmacher, F.R.: The Theory of Matrices. American Mathematical Society (2001)
33. Marsden, J.E., Hughes, T.J.R.: Mathematical Foundations of Elasticity. Dover Publications, Inc., Mineola (1994)
8 -0.3 0.3 0 -0.3 0.3 0 -0.3 0.3 0 -0.3 0.3 0 -0.3 0.3 0
Fig. 1. Color-coded strain fields in short axis mid-slice cross-sections, time frames 1, 6, 11, 16, 21 (a full heart cycle subtends 42 frames), regularized through Gaussian blurring with spatial scale σ = 1.0, for one healthy volunteer. First row: Circumferential strain Ecc. Second row: Radial strain Err. Third row: Shear strain Ecr. Fourth row:
Minimal strain eigenvalue Emin. Fifth row: Maximal strain eigenvalue Emax.
5 10 15 20 -0.3 -0.2 -0.1 0.0 0.1 0.2 0.3 Time Frame Strain MID slice Ecc Err Ecr 5 10 15 20 -0.3 -0.2 -0.1 0.0 0.1 0.2 0.3 Time Frame Strain MID slice Emax Emin
Fig. 2. Temporal evolution of scalar strain quantities over the first half of the heart cycle (i.e. mostly systolic) in mid-slice cross-section. Error bars indicate standard de-viations over the spatial ROI in each time frame, thus capture all sources of variation due to noise, numerics, and true spatial variability. Legends explain the various graphs. Notice the strong correlation between the extrinsic (polar system related) and instrinsic (eigensystem related) strains. (In the eigensystem, shear strain vanishes identically.)
