Dynamic model of the left ventricle for use in simulation of myocardial perfusion SPECT and gated SPECT

Simulation is a useful tool in cardiac SPECT to assess quantification algorithms. However, simple equation-based models are limited in their ability to simulate realistic heart motion and perfusion. We present a numerical dynamic model of the left ventricle, which allows us to simulate normal and anomalous cardiac cycles, as well as perfusion defects. Bicubic splines were fitted to a number of control points to represent endocardial and epicardial surfaces of the left ventricle. A transformation from each point on the surface to a template of activity was made to represent the myocardial perfusion. Geometry-based and patient-based simulations were performed to illustrate this model. Geometry-based simulations modeled ~1! a normal patient, ~2! a well-perfused patient with abnormal regional function,~3! an ischaemic patient with abnormal regional function, and ~4! a patient study including tracer kinetics. Patient-based simulation consisted of a left ventricle including a realistic shape and motion obtained from a magnetic resonance study. We conclude that this model has the potential to study the influence of several physical parameters and the left ventricle contraction in myocardial perfusion SPECT and gated-SPECT studies. © 2003 American Association of Physicists in Medicine.@DOI: 10.1118/1.1589497 #


I. INTRODUCTION
Simulation is a useful tool for elucidating the effect of radiation-matter interactions on images and assessing the accuracy of quantification algorithms.Although anthropomorphic phantoms constituted an early approach to the simulation of cardiac SPECT studies, they were static models that did not permit a simulation of the cardiac cycle.Some attempts to simulate movements in anthropomorphic phantoms have led to simple models of the cardiac cycle. 1,2However, these models failed to accurately simulate realistic movements such as rotation, change of position due to respirationinduced movement, or anomalous movements caused by ischaemia.Moreover, these models are adapted to a particular anatomy and physiology and are not sufficiently versatile to simulate changes in shape and perfusion.
Numerical simulation is another possibility.In essence, numerical simulation consists of mathematically representing two physical distributions: the geometrical characteristics of the radiation source ͑activity map͒ and the spatial distribution of attenuating material ͑attenuation map͒.Information included in activity and attenuation maps is used by deterministic or Monte Carlo projectors during the image formation process.Existing models fall into two groups: geometrybased and voxel-based.Geometry-based phantoms [3][4][5] are derived from simple equations.The main advantage of this type is that different anatomical variations and images can be generated at multiple resolutions.Voxel-based phantoms 6,7 are generally extrapolated from patient data, and are fixed to a particular anatomy and resolution.It is possible to make geometry-based phantoms reasonably realistic but they will not be as realistic as voxel-based phantoms.Recently some attempts at building hybrid phantoms that preserve the best characteristics of voxel-based and geometry-based models have aroused considerable interest. 8,9Both kinds of phan-toms have been used to analyze the effects of anatomical variations in cardiac SPECT images 10 and gated SPECT 11,12 images.
Despite the large number of recent improvements, some features of cardiac SPECT and gated-SPECT studies remain nonrealistically simulated.Even though clinical results have clearly demonstrated that ischaemia causes immediate and severe dysfunctional contraction, 13,14 simultaneous simulation of ischaemia and anomalous movements remain underexplored.This work is focused on reducing the limitations in two areas.The first area is the creation of a model of the left ventricle ͑LV͒, which allows the simulation of left ventricular motion in diseased states, taking care to simulate not only the movement but also the perfusion.The second area is that, when a stable tracer is being simulated, the total and the regional activity within the LV should remain constant throughout the cardiac cycle.

A. Characteristics of the model
In order to model a dynamic LV, a number of characteristics should be taken into account.In this section we seek to provide a summary of the main characteristics of the numerical phantom.
͑1͒ Anatomy.We assume that in myocardial perfusion SPECT and gated-SPECT studies, only the left ventricle is clearly observable, and only this will be relevant in the simulation.6][17] Such a simplification will be acceptable in most perfusion and gated-SPECT studies and allows us to concentrate on improving the left ventricle modeling.This approximation will not be acceptable in cases where high activity in the right ventricle can influence LV activity.For example, after injection of the tracer, a right ventricle blood pool can affect images of the LV.͑2͒ Left ventricle contraction.LV movement consists of the contraction of the cardiac walls.During this contraction the apex remains almost static, whereas the base moves toward the apex, forcing blood into the circulatory system.As a result of contraction, the myocardial wall increases its thickness, with a significant decrease in the perimeter of the endocardial wall and a smaller decrease in the epicardial wall perimeter.In addition, other movements such as rotation or movements due to respiration could be included.An indepth study of the left ventricle movement can be found in Refs.18 -21.͑3͒ Malfunctional contraction.From a kinetic point of view, a diseased state of the left ventricle results in anomalous movements. 13,14Some interesting cases to simulate are akinesis, hypokinesis, and diskinesis and their associated changes in thickening.͑4͒ Rotational movement.It has been established that during systole the base of the LV moves clockwise, and the apical region moves anticlockwise.3][24][25] The midventricular region remains almost stationary during systole.͑5͒ Respiratory movement.Breathing affects the motion of the heart.The diaphragm moves up and down between 1-10 cm during breathing 9,26 and the heart moves in a similar fashion.Respiratory motion can alter left ventricular function and cause artifacts in myocardial SPECT, leading to a misinterpretation of images, especially in the area of the inferior wall of the LV. 27͑6͒ Activity.A fundamental condition to be modeled-when a stable tracer is being simulated-is that the total and regional activities of the heart during the cardiac cycle should remain constant.This means that, although the LV moves and thickens, the amount of activity in a region of the LV does not change during the cardiac cycle.

B. Input parameters
To generate the phantom a number of parameters must be supplied to the algorithms.

͑1͒ Parameters representing endo-and epicardial surfaces.
In this model, the endocardial and epicardial surfaces of the LV are described using a number of control points ͑CP͒.Each CP is specified using four coordinates: three for space and one for time (r,,,t).CP must be carefully placed since regions of the surface without CP will not be satisfactorily modeled.The choice of CPs could be obtained from equations representing the surfaces, from images with good anatomical resolution such as magnetic resonance ͑MR͒ or be specified by the operator.Figure 1 shows two examples of different ways of defining the coordinates of CPs ͑both are used in the following section ''Illustrative models''͒.The upper half of the image shows two short-axis slices of a MR image of a heart at end-systole and end-diastole.CPs used to represent endo-and epicardial surfaces are presented.Different short-axis slices from base to apex should be used to determine the three-dimensional coordinates of the CPs.The lower half of the image shows a geometrical model to represent endo and epicardial surfaces.CPs used to represent endo and epicardial surfaces are distributed uniformly over the surfaces.͑2͒ Parameters representing perfusion.Perfusion is represented in a two-dimensional template.As the ventricular wall moves and thickens, the CPs move in space, but the activity remains fixed to the myocardial wall.It is necessary for each CP, regardless of its position during the cardiac cycle, to transform into the same point on the template.A mathematical transformation between CP and points on the template is used to determine the amount of activity in the myocardium.This transformation is found by using as input parameters pairs of coordinates on the surface (,,t) and their corresponding coordinates (x,y) on the template.The mathematics used to calculate this transformation is detailed in the Appendix.Note that this transformation does not need to be bijective and some points can be transformed to more than one point on the template, as occurs at the apex.
This will not produce any problems if care is taken not to introduce incongruencies in the template.To make the model easier to compute, we consider the transformed points on the template ordered in the vertex of a rectangular array.This restriction helps to simplify the implemented algorithms but does not introduce limitations to the shape and contraction of the modeled LV.
Figure 2 shows a graphical representation of the transformation using a geometrical model.A LV with a number of CPs on the surface is shown on the left.A template of activity and the transformed CP are presented on the right.The template of activity can be fixed by the operator or obtained from patient data.In this paper, the template was represented by a 32ϫ32 matrix and stored in a text file.The matrix was generated automatically to represent uniform perfusion.Ischaemia was simulated by manual modification of some elements of the matrix.

C. Activity and attenuation maps
The numerical phantom of the LV consists of threedimensional matrices representing activity and attenuation distributions.Activity and attenuation maps are needed for each moment in the cardiac cycle to simulate a gated-SPECT study.Three steps are followed to fix the numerical activity to each voxel: ͑1͒ The equations describing the endocardial and epicardial surfaces are determined for each gate.These equations were obtained from the CP and the transformation by using bicubic splines: 28 r en t ϭr en t ͑ ,͒ϭr en t where r en t and r ep t are the equations of the endocardial and epicardial surfaces at time t.͑2͒ For each voxel of coordinates (r,,,t), we verify whether or not it falls within the LV myocardium "r en t (,)ϽrϽr ep t (,)….The coordinates of all these voxels in the myocardium wall are transformed to the template, as described in the Appendix.The activity for these voxels is obtained from the template.͑3͒ Finally, the activity map is normalized to the regional activity of the LV.To this end, the template is divided into sections.Each section corresponds to a part of the LV and can be normalized individually.
Thorax information should be included in order to obtain the total activity map.The mathematical cardiac torso phantom ͑MCAT͒ 5 was employed in this work.Thus, the total activity was obtained by inserting our LV inside the MCAT phantom in the correct position and orientation.The relative activity per voxel assigned to the lung and the rest of the body ͑including ventricular chamber͒ was 15% of the normal myocardial wall.The liver and kidney were not included.
The attenuation map used is the one provided by MCAT.Attenuation coefficients for soft tissue, lung, and spine were 0.150, 0.044, and 0.189 cm Ϫ1 .
To generate a nongated study, all the activity maps can be averaged together during simulation.To generate a gated study, individual time samples are used separately or small groups of adjacent time samples are averaged to simulate movement in every time interval.

D. Projections and reconstruction
Simulation of projections was performed with Monte Carlo methods by using the SimSET code. 29,30Photons of 140 keV were used to simulate a Technetium ( 99m Tc) acquisition.Thallium ( 201 Tl) decays by electron capture, emitting x-ray and gamma photons, with 95% of the emission between 69 and 82 keV.For the sake of simplicity, photons of 71 keV were simulated in a thallium acquisition.In both cases, a low-energy high-resolution collimator and an energy window of 20% were employed.The projections were simulated on a 64ϫ64 matrix with a pixel size of 0.634 cm to simulate standard SPECT and gated-SPECT formats.The postprocessing of all the simulations was performed by using a filtered-backprojection algorithm with a Butterworth filter ͑order: 3.5, cut-off frequency: 0.44 cm Ϫ1 ͒.

E. Illustrative models
To illustrate this model, five studies were simulated.These were created with 48 CPs distributed on the surface.A 128ϫ128ϫ128 matrix with a voxel size of 0.31 cm was employed for activity and attenuation maps.The main characteristics of each simulation were the following.
͑a͒ Case A. Well-perfused patient with normal regional function.In this case, a geometric simulation of the LV of a patient with normal perfusion and uniform wall thickening and contraction was performed.The geometry of this LV was fixed by the operator and consisted of two ellipsoids truncated by the valvular plane.The data used to create the geometric model ͑size, radial, and longitudinal contraction, thickness,...͒ was fixed to be similar to those of a real left ventricle.The contraction, thickness, and ejection fraction can be calculated analytically from the equation of the truncated ellipsoids used to generate the CPs. Figure 3  To simulate normal tissue, the parameters chosen were c max ϭ5, ␣ϭ2, ␤ϭ8, c 0 ϭ1000, ϭ5.To simulate a reversible defect, the parameters selected were c max ϭ3.6, ␣ϭ2, ␤ϭ8, c 0 ϭ900, ϭ5.To simulate an irreversible defect, the parameters chosen were c max ϭ0, ␣ϭ2, ␤ϭ8, c 0 ϭ400, ϭ10.͑e͒ Case E. Patient-based simulation.This case simulates a normal LV with information obtained from a real patient.A gated MR image of a normal subject was acquired using a Tl sequence in eight intervals.The LV was reoriented to obtain the short axis.This simulation included contraction and rotational motion obtained from the MR image by placing the CPs on the surface of the left ventricle.Forty-eight CPs per interval were manually placed on the short axis.However, without specific landmarks on the heart surface, the real contraction and rotational motion were only approximated in the model.Only large geometrical characteristics were considered; papillary muscles were not simulated.
The template used to fix the perfusion was uniform.

III. RESULTS
The results of the illustrative models are shown in this section.Figure 6 gives the results of cases A, B, and C.This figure shows five slices of the simulated LV ͑three short axis, horizontal axes, and vertical axes͒.The first row ͑case A͒ shows the results of the simulation of a patient with a wellperfused myocardium.A decrease in perfusion can be seen from apex to base due to attenuation.The second row ͑case B͒ shows the results of the simulation of a patient with normal perfusion but reduced regional function in the septum.Despite simulating the activity of case A, a lower intensity is observed in the septum.This can be attributed to the fact that the partial volume effect has a higher impact on the septum than on the lateral wall.Thus, example B clearly shows that changes in movement or contraction patterns can yield mis-leading perfusion levels.Although this error would probably not be significant enough to result in a wrong diagnosis, it might influence the calculation of perfusion changes between two studies of the same patient.
The third row ͑case C͒ shows the simulation of a patient with ischaemia in the lateral wall and reduced regional function.This simulates a case that is often encountered in clinical practice.Further studies including varying degrees of hypokinesis and ischaemias in various regions should be undertaken in order to draw conclusions.
Figure 7 gives the results of simulation D. This figure shows a horizontal long axis, a vertical long axis and three short axes of the LV for each simulation.The upper half presents a LV with a reversible defect in the anterior wall and normal tissue in the rest of the myocardium.This left ventricle is represented in stress ͑S͒ and in redistribution ͑R͒ situations.The normal tissue shows a uniform uptake and also a uniform redistribution in all the walls of the LV.The reversible defect in the anterior wall displayed a reduced uptake in this wall in comparison with the normal region.In the redistribution, both regions are indistinguishable.The lower half of the figure shows a LV with an irreversible defect in the anterior wall and normal tissue in the rest of the myocardium.The irreversible defect displayed a lack of uptake and no redistribution.
Figure 8 gives the result of the simulation of case E obtaining the CPs from a gated MR image.The upper row shows the myocardial perfusion SPECT of the patient and the lower row displays the simulated LV.Large anatomical variations can be taken into account to generate a realistic model.This is the case of thinning of the myocardial wall and shortening of the septum.

IV. DISCUSSION
A dynamic LV model was designed and implemented and some examples are reported.Our initial results suggest that the phantom has the potential to simulate normal LV as well as LV affected by different cardiac disorders such as those affecting contraction, morphology, and tracer kinetics.
Geometry-based simulations modeled normal and pathological LVs with contraction and tracer kinetics disorders.In this paper, these disorders were simulated separately in different models.However, the methodology used to generate the LV does not exclude the possibility of including contraction, hypoperfusion, and tracer kinetics in the same model.
Respiratory-induced movement was not included.Nevertheless, it has been treated by other authors. 26Their methodology has the potential to be included in our model so as to simulate respiratory-induced movement.
A patient-based simulation was generated using anatomical information obtained from gated MR images.Figure 8 shows that this methodology has the potential to create a realistic model for the LV including the contraction and change in shape.Anatomical variations can be modeled using several sets of patient data.However, the lack of some information limits the realism of the simulation.Without specific landmarks on the heart surface, no realistic time correspondence can be established for the CPs from frame to frame.In order to realistically and accurately model the motion of a CP from frame to frame, tagged MRI data would have to be used.In addition, a uniform perfusion was simulated which is not the case in actual patient studies.Other sources of discrepancy between modeled and real images are that abdominal structures and papillary muscles were not considered.
Applications of this phantom could include the evaluation of ͑1͒ alterations of the apparent distribution of activity due to contraction of the heart; ͑2͒ the accuracy of methods for determining cardiac global and regional function and ͑3͒ four-dimensional reconstruction methods of gated cardiac SPECT.Moreover, PET and dynamic PET studies could be simulated using PET projectors.

ACKNOWLEDGMENTS
This work was supported in part by Ministerio de Ciencia y Tecnologia, MCYT ͑SAF99-0137 and SAF2002-04270-C02-01/02͒ and Fondo de Investigacio ´n Sanitaria, FIS ͑G03/ 185͒.The work of S. Bullich was supported by a grant from the Institut d'Investigacions Biome `diques August Pi i Sunyer ͑IDIBAPS͒.

APPENDIX: TRANSFORMATION BETWEEN 3-D SPACE AND 2-D TEMPLATE
When the ventricular wall moves and thickens, the activity remains fixed to the myocardial wall.It is necessary to carry out a transformation between points on the twodimensional template and CP.This transformation requires a CP to transform to the same point on the template regardless of its movement.This transformation can be represented mathematically as follows: This minimization has been implemented with the Simplex algorithm. 28sing this methodology, two transformations can be found: a transformation from the points of the endocardial

FIG. 1 .
FIG. 1.An example of the methods used to define the coordinates of CPs.The upper half of the image shows two short-axis slices of a MR image of a heart at the end systole and end diastole.The lower half of the image shows a geometrical model to represent endo-and epicardial surfaces.CPs are distributed uniformly over the surfaces.Different short-axis slices from the base to the apex should be used to determine the three-dimensional coordinates of the CPs.

FIG. 3 .
FIG.3.The horizontal long axis ͑HLA͒, vertical long axis ͑VLA͒, and mid-ventricular short axis ͑SA͒ diagrams of a geometric left ventricle.The continuous line shows the end systole and the dotted line shows the end diastole.The upper row of the figure ͑simulation A͒ represents normal contraction.The motion includes a small movement of the epicardial surface and a larger movement of the endocardial surface.The apex remains almost static and the base moves toward the apex.Rotational movement of the control points was included and is represented by arrows.The central row of the figure ͑simulation B͒ shows a patient with a reduced regional function in the septum.The lower row of the figure ͑simulation C͒ represents a patient with ischaemia and reduced regional function in the lateral wall ͑lateral akinesis͒.The light gray region represents the placement of the ischaemia.
Figure 5 illustrates 201 Tl kinetics for different degrees of defects ͑normal tissue, reversible defect and irreversible defect͒.The time where a SPECT acquisition was performed is marked in gray.Two simulations were performed.The first simulation included a geometrical model, as described in model A ͑but without contraction͒ with a reversible defect in the anterior wall.The second simulation modeled an irreversible defect in the anterior wall.In both simulations, the rest of the myocardium had normal tissue.Each SPECT acquisition was obtained by sampling the activity in five-time frames during the acquisition and by averaging them.

FIG. 4 .
FIG. 4. The figure representing the simulated LV volume during the cardiac cycle for a patient with a normal contraction ͑continuous line͒ and with abnormal contraction ͑dashed line͒.

FIG. 5 .
FIG.5.A figure representing the simulated thallium-201 kinetics in the myocardium.A continuous line represents tracer kinetics for normal tissue, a dashed line represents tracer kinetics for a tissue affected by a reversible defect, and a dotted line represents tracer kinetics for a tissue affected by a persistent defect.The time intervals where the SPECT acquisitions were performed are marked in gray.

FIG. 8 .
FIG. 8.A figure showing the results of simulation E. The upper row displays a myocardial perfusion SPECT of a patient and the lower row displays the simulated LV with the CPs obtained from a MR image of the same patient. a͒ 28 t :͑x i ,y i ͒→͑ en,i where (x i ,y i ) are points on the template, (r en , en , en ,t) are CP on the endocardial surface, and (r ep , ep , ep ,t) are points on the epicardial surface as represented in Fig.1.With x and y as the ordinates in the vertexes of a rectangular array, several bicubic splines28can be adjusted to determine the transformation for the endocardial surface " en,i t ϭ en,i t (x i ,y i ), en,i t ϭ en,i t (x i ,y i )… and epicardial surface" ep,i t ϭ ep,i t (x i ,y i ), ep,i t ϭ ep,i t (x i ,y i )….The inverse transformation to find (x,y) given ͑,͒ can be obtained by minimizing the following two-dimensional functions: g en,, t ͑ x,y ͒ϭ͓"Ϫ en,i t ͑ x,y ͒… 2 ϩ"Ϫ en,i t ͑ x,y ͒… 2 ͔ 1/2 , t ͑ x,y ͒ϭ͓"Ϫ ep,i t ͑ x,y ͒… 2 ϩ"Ϫ ep,i t ͑ x,y ͒… 2 ͔ 1/2 .