Scatter reduction for grid-less mammography using the convolution-based image post-processing technique

X-ray Mammography examinations are highly affected by scattered radiation, as it degrades the quality of the image and complicates the diagnosis process. Anti-scatter grids are currently used in planar mammography examinations as the standard physical scattering reduction technique. This method has been found to be inefficient, as it increases the dose delivered to the patient, does not remove all the scattered radiation and increases the price of the equipment. Alternative scattering reduction methods, based on post-processing algorithms, are being investigated to substitute anti-scatter grids. Methods such as the convolution-based scatter estimation have lately become attractive as they are quicker and more flexible than pure Monte Carlo (MC) simulations. In this study we make use of this specific method, which is based on the premise that the scatter in the system is spatially diffuse, thus it can be approximated by a two-dimensional low-pass convolution filter of the primary image. This algorithm uses the narrow pencil beam method to obtain the scatter kernel used to convolve an image, acquired without anti-scatter grid. The results obtained show an image quality comparable, in the worst case, to the grid image, in terms of uniformity and contrast to noise ratio. Further improvement is expected when using clinically-representative phantoms.


INTRODUCTION
One of the main challenges remaining in digital X-ray mammography is scattered radiation. Scattering reduces the quantitative usefulness of the mammogram, affecting the diagnosis of malignant lesions 1, 2, 3 .
Most of the planar mammography systems currently available in the market make use of anti-scatter grids as a technique to physically reduce the scattered radiation. Although this method helps to improve the image quality, it is not a complete solution. The grids do not block all the scattered radiation and attenuate a proportion of primary radiation, leading to an increase in the patient dose delivered, and add complexity and cost to the overall system 1,4,5 .
Image post-processing scatter reduction methods have emerged as a consequence of the limitations of the anti-scatter grids. Methods such as the convolution-based scatter estimation, have lately become attractive as they are faster and more flexible than direct Monte Carlo (MC) simulations 4,6,7 . This specific method is based on the idea that the scatter in the system is spatially diffuse, thus it can be approximated by a two-dimensional low-pass convolution filter of the primary image 3 .
In this paper, we compare the quality of clinical images acquired without anti-scatter grids that have been post-processed with the convolution-based scatter estimation method, and images acquired with anti-scatter grid.

METHODOLOGY
We assume that the image acquired from a digital mammography detector, I(x,y), is the linear combination of the energy deposited by the primary and the scattered X-ray photons, P(x,y) and S(x,y) respectively: ( , ) = ( , ) + ( , ). (1) The two dimensional (2D) description of the scatter can be seen as a low-pass filtered primary image 8 . However, the primary image is unknown and its determination is the principal objective of this study. To address this problem, we are following an approximation first introduced by Love and Kruger 8 : where, h(x,y) is a low-pass filter kernel that has to be determined.
With this approximation it is possible to calculate a new primary image, P'(x,y):

Mammography system and software used
Phantom images with and without anti-scatter grid were taken at Barts Health NHS Trust with a Hologic Lorad Selenia mammography system. The phantom used was a CIRS 010A test object. Both, grid and grid-less images, were obtained at constant exposure, 28kV and 100mAs.
The GEANT4 toolkit (version 10.01.p02) was the software chosen to study the scattered radiation in the simulated mammography geometry 9,10 . The analysis was done with the aid of MATLAB, R2013b with Image Processing toolbox, and ImageJ, 1.47V.

Definition of the Scatter Point Spread Function (SPSF)
SPSF kernels were obtained from the MC simulation using a narrow X-ray beam which was represented by a normallyincident 2D spatial delta function, i.e. narrow pencil beam method 7 . Energy from scattered radiation (S) was binned into 35µm radius (r) concentric annuli (up to 175mm). These were then normalized to the sum of the primary image and the scatter radiation, P(0)+S(r), and the area of the annuli, A(r), as described in 11 .
The final kernels were calculated from the scatter point spread function, see equation 4. As the kernel was considered isotropic and the function was normalized by the area, the value of every pixel of the kernel could be defined by calculating its distance to the centre, interpolating the SPSF(r) function and multiplying by the area of the pixel, 70x70µm 2 for the Hologic Lorad Selenia. The estimated primary image (P') was obtained following equations (2) and (3).

Geometry used in this study
The experimental geometry, i.e. Pencil Beam (PB) geometry, was defined to simulate the geometry of the Hologic Lorad Selenia mammography system, where the CIRS 010A phantom images were acquired. For the pencil beam methodology, the D-shaped CIRS phantom was changed to a cylindrical-shaped, simulating the scatter distribution at the centre of the phantom. From top to bottom, the simulated geometry included: X-ray source, compression paddle, breast phantom, support paddle, air gap and the detector. The Source to Image Distance (SID) was 660mm. The detector structure, up to the amorphous selenium (aSe) surface, was included in the geometry 12 , the data was obtained via private communication with Hologic Inc. Figure 1 shows both the realistic and pencil beam geometries and Table 1 gives more details about the composition and thicknesses.
In this report, the SPSF kernels were simulated using a narrow pencil beam and a spectrum of 28kVp W/Rh (HVL=0.538 mmAl). The photon cross-section "Electromagnetic physics option 4 package" was used both for the validation and the experiment. 20-60 runs, of 10 9 X-ray photons each, were simulated, ensuring uncertainties lower than 5% and calculated as suggested by Sempau et al. 13 . Figure 1. Diagram of a realistic geometry, used for validation purposes, (A) and of the experiment geometry (B). Figure A shows an example using a cone beam aligned with the chest wall of the breast phantom. Figure B shows the pencil beam example and a detector added after the air gap.

Prediction of the scatter-reduced image
To account for the difference in scattering between the breast phantom and the background areas, two different kernels were simulated, KPhantom and KBackground. For the background kernel, the pencil beam geometry was simulated with the absence of the phantom. The analysis only needed of these two kernels because of the characteristics of the chosen phantom. More kernels would need to be used if, for example, the thickness of the phantom was variable.
During the processing step, a segmentation algorithm based on thresholding was performed in order to separate both phantom and background areas, before the two convolutions were performed. The final scattering image was the sum of both convolutions, see Figure 2.
The predicted primary image (P') was calculated by subtracting the predicted scattering image, S'(x,y), from the grid-less image.

RESULTS
The processed grid-less image was compared to the image acquired with the use of a grid using several imaging metrics. Both images were acquired under the same acquisition conditions and dose in order to facilitate the comparison.
Before the scatter kernels were calculated, the geometry of the system and the kernel curves were validated with previous publications.

Geometry Validation
The values given in the report of the American Association of Physicists in Medicine (AAPM) task group 195, case 3 -Mammography and breast tomosynthesis, were used as a benchmark for validating the PSF simulations 14,15 . The results obtained showed an average discrepancy lower than 1% when compared with the results from the AAPM report.
The geometry followed in the validation is shown in Figure 1-A. Table 1 gives details about the materials, thicknesses and positions. A combination of two X-ray beams (cone and pencil beam) and two X-ray spectra (monoenergetic and polyenergetic) were used (see Table 2). The spectrum values were obtained from 14 . For each of them, the primary radiation, and scattering radiation (single and multiple Compton and Rayleigh events) were measured in 7 different regions of interest (ROIs), as described by 14 . Enough particles were simulated sufficient to produce Standard Error of the Mean (SEM) values equal to or lower than 1%.

SPSF Validation
The SPSFs have been validated against the data published by Diaz et al. 7 and Sechopoulos et al 16 . Table 3 shows the comparison with the previous published data, the results show excellent agreement with both Diaz and Sechopoulos. In this geometry the compression and support paddles were removed, leaving only a cylindrical phantom and an air gap between the phantom and an ideal detector. For all the cases, an energy spectrum of 26keV Mo/Mo and an air gap of 10 mm was used.

CIRS 010A Phantom: Grid vs. post processed grid-less images
The post-processed grid-less image, P' equation (3), and the grid image were compared in terms of uniformity and Contrast to Noise Ratio (CNR).

Uniformity
A comparison of the unprocessed grid-less image (I), the processed grid-less image (P') and the grid image is shown in Figure 3.
A horizontal profile, as highlighted in Figure 4, has been plotted for the three images. Grid-less Grid-less Processed Grid

Contrast to noise ratio
The CNR was calculated following equation (5). The area analysed corresponds to the area surrounding the step phantom, the equation was applied to the 5 different regions of interest (ROIs), see Figure 5-left.
where, ̅ and ̅ are the mean pixel value of the object and background ROI, respectively.
̅ is the standard deviation of the background ROI.

DISCUSSION
The analysis with this convolution-based technique shows an overall improvement of the processed grid-less image over the unprocessed grid-less image, and comparable results when compared with the grid image.
In terms of uniformity, the unprocessed grid-less image presents a non-uniformity higher than 100 grey values. This nonuniformity is reduced by more than 85% in the processed grid-less image. When compared to the grid image, it is possible to see that the image is flatter in the central region (5 vs. 10 grey values).
Towards the edges of the image, the non-uniformity of the processed image increases. In future measurements, further corrections should be applied around these areas, as suggested by Diaz et al. 7 . Towards the chest-wall area, the discrepancy is probably introduced by the effect that the plastic holder of the phantom (outside the collimated area) introduces in the scatter.
The CNR measurements give comparable results to those obtained with the grid image. For the first ROI of the step phantom the CNR is slightly lower, while ROIs 2-5 show an improvement.
An improvement in the overall result is expected for real mammographic images. The objects that can be found inside the phantom present an extreme case, covering different ranges of materials and introducing differences in the scattering. Taking the step phantom area, for example, the steps go from 0 to 100% glandular tissue. These differences are not taken into account in the convolution, as only one kernel that matches the phantom general material is used.

ROI
Processed grid-less image Grid image

CONCLUSIONS
Measurements made with a commercial planar mammography system were used with convolution-based image scatter reduction techniques in an attempt to generate images with high contrast to noise ratio, without the use of an anti-scatter grid.
Analysis of the images showed that the quality of the predicted primary image is comparable to the grid image in terms of uniformity and presents general improvements in the contrast to noise ratio. The results are expected to show further improvement when using more clinically-representative phantoms.
Further studies will include more realistic phantoms and a more complex analysis of the breast edge areas to gauge the level of complexity required to build the scatter point spread functions in a real system. The results obtained so far indicate that future grid-less mammography systems will be able to reduce the dose delivered to the patient without a compromise in the image quality.