Monte-Carlo simulation of a slot-scanning digital mammography system for tomosynthesis

1 Introduction

Full-field digital mammography (FFDM) has emerged as an improvement on the conventional screen-film mammography in the past two decades. FFDM performs better than screen-film mammography for women under 50 years of age, but its sensitivity decreases to less than 50% for denser breasts [23].

The two-dimensional (2D) images produced by FFDM may show overlapping breast tissue because three-dimensional (3D) anatomical data is projected on to two dimensions in the image. Small tumours can be obscured due to overlapping tissue, resulting in a high number of recalls, more screenings and possible biopsies. The false negative rate for mammography is estimated at one in five breast cancers that are present at time of screening [2]. It is also estimated that for women who undergo annual mammography screening for a 10-year period in the USA, approximately 50% of cases produce at least one false positive result [20] with the necessary additional tests such as diagnostic mammograms, ultrasound scans and biopsies creating anxiety in the patient and increasing the healthcare costs [2]. Digital breast tomosynthesis (DBT), a 3D alternative to FFDM, has been shown in several studies to improve detection rates compared to FFDM and to be superior in detecting multiple breast tumours [19]. Also, combining FFDM with DBT has been shown to improve detection rates in dense breasts compared to using FFDM alone, leading to reduced recall rates [11, 37]. DBT has the potential to improve screening sensitivity and specificity, and its clinical application is predicted to have the potential to become the gold standard for breast imaging [12, 13].

In tomosynthesis, multiple images of the compressed breast are acquired at different angles and the projections are reconstructed to display images of planar slices through the breast. The optimization of DBT parameters namely number of projections, total angular range and angular intervals remains an area of research interest [10] particularly to improve image quality while reducing patient dose. An increase in the number of projections and in angular range improves the image quality. Optimizing these parameters is useful for the design of DBT systems.

Several manufacturers apply different methods to acquire tomosynthesis data [19]. These differences make clinical comparisons difficult. Tomosynthesis machine manufacturers vary the arc of the movement, the number of exposures and projections and the angular intervals between projections. In addition to these tomosynthesis parameters, mammography parameters like the system resolution, source-detector geometry, binning of image data and X-ray source/filtered source can also be varied.

Various models have been developed to study FFDM and DBT systems, but they simulate only parts of the imaging chain like the X-ray source energy spectrum [28] and have limited flexibility to enable simulation of different phantoms or geometries [22]. Very few simulation tools are available for the simulation of slot-scanning imaging systems, although the Monte-Carlo method has been shown to be effective for modelling a slot-scanning system [24]. However, to the authors’ knowledge there are no simulation tools specifically designed for slot-scanning DBT. Helvie et al. [19] recognizes DBT as a new modality and concludes that favourable initial trial results need to be confirmed with more representative trials and the effect of DBT parameters on the detection of micro-calcifications needs to be explored further. A flexible simulation system that models the physics of X-ray photon transport accurately and includes the geometric requirements of DBT is lacking. This paper describes such a simulation system that can be used to model FFDM and DBT for slot-scanning imaging systems.

This paper presents a tissue-equivalent 3D breast phantom, its validation for slot-scanning digital mammography and its use in simulating digital breast tomosynthesis. The mathematical phantom simulated is based on the physical tissue-equivalent phantom designed for mammography (Model 011A) manufactured by CIRS. The model was validated by comparing images acquired using a slot-scanning digital mammography machine with simulated phantom images, using the edge-test method and modulation transfer function (MTF), noise power spectrum (NPS) and detective quantum efficiency (DQE). Different 2D projections of the 3D phantom were simulated and the phantom was reconstructed using filtered backprojection.

2 Background

Computer simulation techniques have been widely used to perform X-ray photon simulations [31]. Initially, algorithms based on the Boltzmann equation were used to model radiation transport problems to simulate the flight paths of X-ray photons, but had several limitations for simple geometries [3]. The use of more accurate but computationally intensive algorithms like the Monte-Carlo based simulation has been made possible by improvements in computer technology [3]. Several implementations of Monte-Carlo techniques are available; the PENetration and Energy LOss of Positrons and Electrons (PENELOPE) [31] is used in this paper.

Several 2D and 3D mathematical simulations of breast tissues and micro-calcifications have been developed. Bakic et al. [8] developed a three dimensional breast phantom based on a realistic distribution of tissue structures. A compression model was used to change the size and placement of structures in the phantom. A mammogram acquisition model was used to simulate synthetic mammograms using a monoenergetic parallel beam source applied to the compressed breast phantom. Mammograms simulated using this system were printed on film and shown to radiologists. Although the simulated images were not validated with real X-ray images it was concluded that the proposed model could be used to study the correlation between 3D breast composition and 2D mammograms.

Hussein et at. [22] defined a tissue-equivalent mathematical breast phantom, based on the commercially available CIRS phantom, and simulated the phantom using a prototype of a linear slot-scanning mammography system. Pixel intensities in different areas of the phantom were estimated using empirical equations for variables in the X-ray imaging chain. The simulated and real images of the phantom were compared using image contrast, contrast-to-noise ratio, signal-to-noise ratio and detective quantum efficiency. Monte-Carlo methods were not used by [22] and although images from a slot-scanning system were simulated, individual slots were not simulated. The effects of different methods of combining the slots to acquire the final image could therefore not be studied. Pixel intensities for the defined phantom were simulated, thus the simulation process was not independent of the phantom used. Furthermore, the phantom was two-dimensional so it was only useful for simulating mammography and did not accommodate different geometries, materials and phantom configurations and 3D modelling.

Hunt et al. [21] used high-resolution voxelized phantoms in modelling the breast to study and optimize the performance of digital mammography systems. The voxelized phantoms simulated adipose and fibroglandular tissue, ligaments, ducts and skin in three dimensions. The performance of digital mammography systems, specifically the radiation dose, noise in each voxel and scatter-to-primary ratios were studied. The scatter-to-primary ratios for different points in the image were estimated using a Monte-Carlo based simulation process. The simulation results were compared with measurements on images of sheets of polymethyl methacrylate obtained using a GE Senographe digital mammography unit. Based on the comparisons of contrast, scatter-to-primary ratio, radiation dose and signal-to-noise ratio it was concluded that the phantoms were an improvement on the homogeneous phantom used in the past for modelling breast tissue.

Reiser et al. [29] developed an analytic breast phantom consisting of simple shapes, aiming to capture the main features of the breast, to study DBT reconstruction. Projection data at various angles were computed analytically and the three dimensional volumes were reconstructed using filtered backprojection, expectation maximization and total variation algorithms. The total variation algorithm is a denoising technique that removes unwanted detail from the signal while preserving important details such as edges. The reconstruction results showed that the total variation denoising techniques used with filtered backprojection achieved highest contrast for mass lesions and best in-depth resolution.

Sechopoulos and Ghetti [33] evaluated different combinations of angular range and number of projections using DBT computer simulations in order to characterize their effects on the reconstructed image. The study was restricted to the effects of these two parameters and the reconstructed images were compared using contrast and artefact spread function. Although the angular intervals were not part of the study, it was concluded that increasing the angular range always increases the vertical resolution of the reconstruction. Vertical resolution decreased when the number of projections in the reconstruction set went over a certain threshold, which depended on the angular range. It was noted that the angular range was more important than other DBT variables to attain good reconstructed images.

Choi et al. [10] used a prototype DBT system to study the effect of the tomosynthesis parameters on the projections and the resulting reconstruction. The contrast-to-noise ratio and artefact spread function were used as metrics to compare images acquired with different sets of parameters. The CIRS breast phantom (Computerized Imaging Reference Systems, Norfolk, VA, USA) was used as a test object and the parameters were optimized using 32 acquisition sets with six angular ranges and eight projection views. The study concluded that a wide angular range improves image quality, especially in the vertical direction, but a large number of projections increases the electronic noise and decreases the contrast-to-noise ratio. A wider angular range also increases scan time, leading to higher probability of patient motion which degrades the image quality. This prototype system only tested the effects of tomosynthesis parameters and did not vary parameters related to mammography like the X-ray source, binning, system geometry and system resolution.

Baptista et al. [9] used the Monte Carlo N-Particle eXtended (MCNPX) code [27] to generate DBT projection data of a homogeneous breast phantom. Dose and image quality measures were used to evaluate DBT reconstruction using different combinations of parameters and at different X-ray source energies to optimize image quality and dose.

3 Methodology

This section describes the mathematical tissue-equivalent breast phantom, the components of a generic digital mammography system, their simulation using PENELOPE and extension of the latter for the generation of an X-ray image. The slot-scanning simulation is based on that described in [24] for a general-purpose slot-scanning digital X-ray imaging system. The parameters used for simulating DBT and re-constructing the projections are provided.

4 Results

The raw data from the simulated X-ray image is shown in Fig. 4. The phantom has been simulated including all its components and their properties. The top and bottom edges of the simulated image appear slightly low in intensity because this area is not covered by as many slots as other regions in the image. This is a common slot-scanning artefact on the edges of the image.

The simulation of the Nylon fibres is shown in Fig. 5. The differences in fibre thickness are clearly visible and their position and orientation are correct. The edges perpendicular to the slot direction are sharper than the ones perpendicular to the scan direction. This is due to the addition of the slots and is a characteristic of slot-scanning X-ray imaging. Post-processing techniques are used to correct this effect. The centres of the hemispheric masses shown in Fig. 6a are equally spaced but their varying diameters give the impression of varying spacing between them. These hemispheric masses, their location and thicknesses, have been correctly simulated. The CaCO₃ specks are clearly visible as shown in Fig. 6b, although their edges perpendicular to the slot direction are sharper than those perpendicular to the scan direction. Their location and orientation in the shape of a  + sign is simulated correctly. The line pair tool, shown in Fig. 7, clearly illustrates how edges perpendicular to the slot direction are sharper than the ones perpendicular to the scan direction. This is a characteristic of slot-scanning imaging that is simulated. The 3 mm tungsten block was simulated and imaged using Soteria to compare the two systems using MTF and DQE. The contrast, spatial resolution and noise characteristics can be compared using these image quality parameters. Figure 8a and b show the ESF and LSF functions respectively. The ESF of the real X-ray image (continuous line) is very similar to that of the simulated image (dashed line) in Fig. 8a, with the edge in the simulation being slightly noisier. The LSF, which is the derivative of the ESF, is also very similar in both the real X-ray (continuous line) and simulated X-ray (dashed line) images, with the simulation being slightly noisier. The normalized MTFs of the real X-ray image (continuous line) and simulated image (dashed line), shown in Fig. 9a, are comparable. These results show that the transfer functions of the simulation process are similar to those of the real mammography machine and their contrast and spatial resolution characteristics match. The normalized DQEs of the real X-ray image (continuous line) and the simulated image (dashed line) overlap each other indicating that the noise modelling of the two systems is also consistent andsimilar.

A selection of the angular projections of the tungsten block are shown in Fig. 10. As the angle of the projection increases from 90 to 180 degrees the size of the projection changes. The reconstructions of vertical slices through the tungsten block are shown in Fig. 11. The quality of the reconstructions decreases as the angular interval increases from 2 degrees in Fig. 11a to 10 degrees in Fig. 11c due to the decrease in the number of angular projections used to reconstruct the planar slice. In Fig. 11a using an angular interval of 2 degrees, 91 projections are used to reconstruct the planar slice compared to Fig. 11b in which 37 projections are used and Fig. 11c in which 19 projections are used. The reconstruction artefacts seen on the edges of the reconstructed slices are periodic.

Figure 12 compares the horizontal planar reconstruction of the tungsten block and its 90 degree projection. The reconstruction is performed using projections with angular range of 0 to 180 degrees and angular intervals of 2 degrees between projections. The intensity variations and artefacts visible in Fig. 12 are the result of the addition of many images. However, the crucial property of the images is that one can distinguish background and foreground. Changing parameters like number of images used and type of filters will give varying results. The results suggest that the simulation process can be modified to study the effects of DBT variables like angular intervals, number of projections and angular range on the resulting horizontal and vertical reconstructions.

Figure 13 shows four reconstructions of the same horizontal planar slice to illustrate the improvement of image contrast with increasing angular range and decreasing angular interval. Figure 13a and c have an angular range from 70 to 110 degrees with respect to the detector and Fig. 13b and d have an angular range from 80 to 100 degrees. The angular interval for Fig. 13a and b is 2 degrees and for Fig. 13c and d is 4 degrees. The number of projections is 21 for Fig. 13a, 11 for Fig. 13b and c and 5 for Fig. 13d.

The image contrast measures for Fig. 13b and d, which both have an angular range of 80 to 100 degrees but respective angular intervals of 2 and 4 degrees, are 204.6 and 199.5 respectively. For the angular range of 70 to 110 degrees (Fig. 13a and c), the contrast measures are 290.5 and 269.7 for angular intervals of 2 and 4, respectively. The best reconstruction result in terms of contrast is shown in Fig. 13a with 21 projections.

5 Discussion

Simulations of phantoms and digital breast tomosynthesis (DBT) systems for breast imaging are useful stages in the design and optimization of such systems [22, 33]. The aim of this study was to build on previous work [22, 24] to develop a tissue-equivalent 3D breast phantom and validate its application for slot-scanning digital mammography and use it in simulating digital breast tomosynthesis.The work done in [22] describes a two-dimensional phantom based on the CIRS phantom. In this paper, the mathematical, tissue-equivalent phantom [22] has been extended to three dimensions and its application for DBT has been shown. This extension allows the phantom to be rotated at any angle to facilitate DBT simulation. Also, unlike in [22], in this study individual slots have been simulated and this would enable the effects of different procedures used to combine the slots to be studied.

Results have shown that the transfer function and noise modelling of the Monte-Carlo simulation method correspond to those of a real mammography machine and that the Monte-Carlo method can simulate the phantom including all its components and their properties. These results agree well with those obtained when the Monte-Carlo approach was used to model a general-purpose slot-scanning digital X-ray system [24]. The top and bottom edges of the simulated images have low intensity due to a common artefact of the slot-scanning technique which results from fewer slots covering the image edges than the inner regions of the images. Another consequence of the slot-scanning technique on the quality of the simulations is difference in resolution between the slot and scan directions caused by the addition of individual slot data to produce the output image. The result is that edges parallel to the scan direction are sharper than those in the slot direction. However, post-processing techniques can be used to correct this effect. Also, the resolution of the raw images has been shown to be good enough to discern thin structures as seen in the simulation of the nylon fibres.

Different DBT reconstruction algorithms will yield planar slices of different quality and some may be better for the detection of calcifications than others [19]. The study done in [29] simulated a breast phantom to analyse the performance of different DBT algorithms. Contrast measures and vertical resolution were used as metrics for the analysis. The data required for such a study can be simulated using the system described in this paper.

The effects of parameters specific to DBT, like the angular range, angular interval and number of projection views, on the reconstructed images in a prototype system were considered in [10]. Angular projections were acquired at different angles and intervals to study limited-angle tomosynthesis and sparse-angle tomosynthesis. Only tomosynthesis parameters were studied and mammography parameters were kept constant. In this study, the contrast of the reconstructions improves with an increase in the number of projections but the noise in the reconstructed images also increases. These observations are similar to those of [10], which studied the effect of DBT parameters on the reconstructed images. The work described in this paper removes the need to construct a physical prototype and extends the work in [10] to enable the effects of both mammography and tomosynthesis parameters to be studied.

6 Conclusion

This paper has developed a general and flexible simulation tool for slot-scanning X-ray mammography and tomosynthesis. The tool can be used to simulate different phantoms and X-ray imaging systems. It has been validated using a slot-scanning digital mammography machine using standard image quality metrics. The simulation tool can be used to study the effect of tomosynthesis parameters on the projections and reconstructions of planar slices through the test object in order to aid system design.