Improved dosimetry protocol for DBT

$D_{\text{gN}}^{\text{DBT}}$, which is the most essential quantity for estimating the $D_{\text{g}}$ of patients, was calculated for each combination of equipment geometry, X-ray spectrum, and breast model using Monte Carlo (MC) simulations. These values, expressed in mGy/mGy, represent the ratio of $D_{\text{g}}^{\text{DBT}}$ (in mGy) to the incident air kerma ($K_{\text{air}}$) (in mGy) at the reference point without considering backscatter. The reference point was positioned along the central axis of the upper surface of the breast, 4 cm from the edge of the detector closest to the chest wall. For each spectrum and target/filter combination, the half-value layer (HVL) was calculated using MC simulations [33] and expressed as the aluminum thickness with equivalent attenuation properties.

D_g is significantly affected by the depth of the major energy deposited in the glandular tissue. MC simulations have demonstrated that the dose distribution during mammography exhibits a high degree of heterogeneity [13, 34]. The deposited dose in nearly 40% of glandular voxels was found to exceed the $D_{\text{g}}$ within the breast model with a 4-cm compressed breast thickness (CBT) and 50% glandularity [34]. An alternative dose metric known as glandular depth dose ($D_{\text{g}}^{\text{dep}}\left(z\right)$) has been suggested to investigate the differences in $D_{\text{g}}^{\text{DBT}}$ between breast models with different glandular distributions [13, 34]. $D_{\text{g}}^{\text{dep}}\left(z\right)$ provides more detailed information on the distribution of the radiation dose along the depth direction, especially in areas with dense glandular tissue. Specifically, $D_{\text{g}}^{\text{dep}}\left(z\right)$ is defined as the glandular dose deposited in a specific slice at depth z.

To better understand the contribution of each voxel slice of the detailed breast model to $D_{\text{g}}$, we used a physical quantity: the normalized deposited energy in the glandular tissue at slice depth z, $E_{\text{gN}}\left(z\right)$, which can be calculated using the following equation: $E_{\text{gN}}\left(z\right)=\frac{D_{\text{g}}^{\text{dep}}\left(z\right)}{D_{\text{g}}}\cdot \frac{m_{\text{g}}\left(z\right)}{m_{\text{g}}}=\frac{D_{\text{g}}^{\text{dep}}\left(z\right)}{D_{\text{g}}}\cdot f_{\text{g}}\left(z\right),$ (1) where $m_{\text{g}}\left(z\right)$ is the glandular mass at slice depth z, and $f_{\text{g}}\left(z\right)$ is a function that represents the distribution of the glandular mass along the depth direction.

In traditional dosimetry, the $D_{\text{gN}}$ coefficients measure the $D_{\text{g}}$ undergoing mammography equivalently as multiplication results of the factors g, c, and s in Dance’s equation [8, 9, 35]. The g-factor is a dose conversion factor calculated under different spectra with reference to the breast model with 50% glandularity at different CBTs. The c- and s-factors are adjusted for an X-ray spectrum different from that of the Mo/Mo target/filter and glandularity different from 50%, respectively. The relative differences $\text{Δ}$ of $D_{\text{gN}}^{\text{DBT}}$ between the improved and traditional breast dosimetry can be calculated via $\text{Δ}=\frac{gcsT-D_{\text{g}}^{\text{DBT}}}{gcsT}\times 100\%,$ (2) where g, c, s, T, and $D_{\text{gN}}^{\text{DBT}}$ are the conversion factors under the same breast characteristics and beam conditions, corresponding to the traditional dosimetry standard and this study.

Detailed breast model and the improved combination with CRAF

In a previous study, we identified typical parameters of the breast, including external parameters such as the base diameter of the breast, the distance from the nipple to chest wall, skin thickness, and subcutaneous adipose thickness[28, 36]. Based on the assumptions of Bakic et al. [37] and Mahr et al. [38] regarding the anatomical growth of the lactiferous ducts and glandular tissue in the fibroglandular region, we altered the glandular distribution by randomly sampling adipose lobules to replace the glandular tissue. This allowed us to approximate the glandular distribution of our target breast model as close to that observed in clinical settings. The percentage of glandular tissue and the CBT in our target models were typical parameters selected from those reported for Chinese women in clinical statistics [28].

Based on the representative breast parameters, we developed the first set of detailed breast models for Chinese females. Table 1 lists the voxel sizes of models with different glandularities and CBTs. Each voxel represented a unique tissue. As the CBT decreased, the voxel size in the depth direction decreased to ensure the completeness of detail in the breast models. These models included four breast regions: skin, the adipose tissue region (subcutaneous adipose, posterior adipose, and Cooper’s ligaments), fibroglandular region (intraglandular adipose, glandular, lactiferous ducts, and lobules), and nipple region (lactiferous sinus and adipose). The skin thickness was set at 1.45 mm. The thickness of the subcutaneous adipose layer was set at 4 mm near the nipple and slightly greater than 4 mm near the chest wall. The fibroglandular region was the central region in the breast model, excluding the skin, subcutaneous adipose layer, and posterior adipose layer, similar to the central region in the simple breast model of the traditional standard. Detailed breast models with glandularity levels of 5%, 25%, 50%, 75%, and 100% were constructed. Glandularity is a physical quantity that indicates the proportion of glandular tissue and can be expressed as percentage volume or percentage mass glandularity. The percentage mass glandularity ($P_{\text{g}}^{\text{M}}$) can be calculated using $P_{\text{g}}^{\text{M}}=\frac{m_{\text{g}}}{M_{\text{fg}}}\times 100\%,$ (3) where $m_{\text{g}}$ is the mass of the glandular tissue, and $M_{\text{fg}}$ is the total mass of the fibroglandular region. The percentage volume glandularity ($P_{\text{g}}^{\text{V}}$) can be calculated via $P_{\text{g}}^{\text{V}}=\frac{V_{\text{g}}}{V_{\text{br}}}\times 100\%,$ (4) where $V_{\text{g}}$ is the volume of the glandular tissue, and $V_{\text{br}}$ is the volume of the breast model. The corresponding $P_{\text{g}}^{\text{V}}$ values for models with 5%, 25%, 50%, 75%, and 100% $P_{\text{g}}^{\text{M}}$ were 1.6%, 8.2%, 16.6%, 25.4%, and 34.4%. In the fibroglandular region, the 5% $P_{\text{g}}^{\text{M}}$ model contained only lactiferous ducts and intraglandular adipose tissue, whereas the 100% $P_{\text{g}}^{\text{M}}$ model contained only glandular tissue. Models with different glandularities were constructed by altering the number of adipose globules whose central points were uniformly sampled within the fibroglandular region.

Vertical slices of each breast model with different glandularities were deformed to compress these models to 2, 3, 4, 5, 6, and 7 cm in the CC view. This resulted in 30 compressed breast models for dose estimation in DBT. The compression algorithm divides the breast model into skin, the adipose tissue region, and fibroglandular tissue region and calculates the elasticity parameters of the breast tissues based on ultrasound velocity measurements and tissue densities. The algorithm separately compresses each vertical slice of breast tissue using the elastic modulus of the corresponding tissue region and determines the amount of compression based on the strain, which is calculated using the change in the thickness of the breast tissue. The algorithm preserves the volume of the breast during compression, as biological materials are considered incompressible, and considers the original dimensions and thicknesses of the breast model to calculate the final compressed dimensions [28]. Detailed breast models, in which the glandular tissue tends to be more concentrated in the central part of the sagittal plane, have glandular distributions that are different from those of the simple model. Previous research has indicated that, based on the simple breast model, the beam condition (target/filter combination, beam HVL) and breast characteristics (CBT and glandularity) are the main influencing factors of $D_{\text{gN}}^{\text{DBT}}$, in addition to the tomosynthesis system specifications [39].

To investigate the difference in $D_{\text{gN}}^{\text{DBT}}$ between the detailed and simple breast models, we constructed a series of homogeneous breast models with 4 cm CBT and 25%, 50%, 75%, and 100% $P_{\text{g}}^{\text{M}}$ by uniformly sampling the glandular voxels in the central part.

In our previous study, we combined the compressed breast model with the Chinese reference adult female whole-body voxel phantom (CRAF) to account for backscatter radiation from the female body. However, this reduced the computational efficiency of the MC simulation owing to the millions of voxels in CRAF. To enhance the efficiency and accuracy of the MC simulation, we improved the geometric construction of CRAF. According to the position of the breast in the female body, CRAF was rotated and cropped, and the major organs and tissues that provided backscattered particles were present. During construction of the imaging geometry, the nipple of the breast model was aligned with that of CRAF. The CRAF and cropped CRAF models are shown in Fig. 1.

Fig. 1Full-screen View

(Color online) (a) 3D display of the CRAF model, where the pink region represents the breast organ in CRAF. (b) Improved method to combine the breast model and CRAF. The major organs, including the heart, lung, and posterior muscle of the breast, were reserved, and the body tissue away from the breast was removed. The breast organ in CRAF was replaced via the detailed breast model

Irradiation geometries

The simulation geometry used to calculate two quantities, $D_{\text{g}}$ and $D_{\text{g}}^{\text{dep}}\left(z\right)$, is shown in Fig. 2(a). The geometry consisted of an X-ray tube and a series of geometric components. The X-ray tube had an X-ray source and a filter below it, which was employed for energy spectrum hardening. During the scanning process, the source and filter rotated together along the orbit around the center of rotation. The geometric components included a compression paddle, breast model, support paddle, and image receptor, which were aligned with the chest wall plane. The compression paddle was constructed from polycarbonate, whereas the support paddle was composed of carbon fiber. An air gap existed between the support paddle and the image receptor. The position of the compression paddle placed on top of the breast model changed with the CBT. The isotropic particle emission was set up at a cone angle that covered the image receptor edge. Four commercial scanning geometries were constructed. The device parameters are listed in Table 2.

Fig. 2Full-screen View

(Color online) Irradiation geometry to measure (a) glandular dose and (b) K_air

To simulate $K_{\text{air}}$, the breast model and support paddle were eliminated and a cylindrical air-filled ionization chamber (15 mm radius, 10 mm height) was positioned at the breast entrance plane [8, 40]. The compression paddle was elevated 40 cm above the upper plane of the ionization chamber to reduce scattered radiation from the paddle [8]. Air filled the world in all the simulations. The simulation geometry is illustrated in Fig. 2(b).

Monte Carlo simulation

A previous MC computer program for mammography dosimetry was modified to calculate $D_{\text{g}}$ for DBT using the Geant4 (version 10.6) simulation toolkit [41]. The irradiation geometries are described in detail in Section 2.3. A detailed validation of the mammography dosimetric program can be found in [28]. A simple breast model with 4-cm CBT and 50% $P_{\text{g}}^{\text{M}}$ was constructed using voxels to validate the modification of the X-ray tube rotation angle. This model had a uniform distribution of voxels labeled as glandular and adipose in the central part. Four different spectra (Mo/Mo 25 kV, Rh/Rh 29 kV, W/Rh 34 kV, and W/Ag 40 kV) were used to irradiate the model for projections from 0° to 30° in 5° increments.

Dosimetry quantities, glandular doses, and $K_{\text{air}}$ were estimated using independent MC simulations. Only photons were considered in the simulations (cutoff energy of 10 eV and cutoff length of 0.1 mm). All the secondary electrons were regarded as being locally deposited. The number of photons simulated for $D_{\text{g}}$ was between ${10}^8$ and ${10}^9$, depending on the parameters of the input models (breast glandularity and CBT). The statistical uncertainties for glandular dose and $K_{\text{air}}$ were less than 1%. The simulations were conducted using an Intel® Xeon® CPU E5-2680@2.3 GHz. The computational time required to process ${10}^8$ photons using a single core was several hours. The calculation speed was significantly improved to ensure simulation accuracy. To ensure that the statistical uncertainty for $D_{\text{g}}^{\text{dep}}\left(z\right)$ at all depths was less than 1%, the number of photons simulated for $D_{\text{g}}^{\text{dep}}\left(z\right)$ was set to $1.3\times {10}^{10}$. Python (v. 3.9) [42] scripts were written to generate input files and shell scripts for simulation automation.

The properties of adipose, glandular, and skin tissues, including their composition and density, were obtained from ICRU-Report 46 [43]. Values for Cooper’s ligament tissue in the detailed breast model were not specifically provided by the ICRU; therefore, the density and element composition of Cooper’s ligament tissue were substituted with those of muscular fibrous tissue [38].

The polychromatic X-ray spectra acquired from the spectral models of Boone et al. [44] and Hernandez et al. [45] were simulated. Simulations were conducted to cover the typical parameters available in clinical systems using various tube potentials (25, 28, 30, 32, 35, and 49 kV) and target/filter combinations (W/Rh, W/Al, Mo/Mo, Rh/Rh, and Rh/Ag). The polychromatic spectra exhibited an energy resolution of 0.5 keV. The heeling effect was not implemented.

Validation for simulation methods

As shown in Fig. 3, the t-factors reported by Dance et al. [9], which capture the variation in $D_{\text{gN}}$ due to non-zero projection angles of an X-ray source, were compared with those generated using the MC code employed in this study for four spectra and projection angles $\alpha$ ranging from 0° to 30°. These results are consistent with data obtained by Dance et al. All data differed within 3%, with the largest difference occurring when the angle was largest. The main reason for this variation is the method used to construct the simple model. The energy deposited in the adipose and glandular tissues in the central region of the simple breast model adopted by Dance was calculated based on the probability of interaction between the two tissues. However, the energy deposited in the glandular tissue could be obtained directly because of the separation of the adipose and glandular tissues in the voxel model used in this study.

Fig. 3Full-screen View

Comparison between the validated results and Dance’s data

Improved computational efficiency and accuracy considering backscatter

Table 3 lists the MC simulation computation times and the differences between the three combinations under the same irradiation parameters and breast model. These combinations included breasts without CRAF, CRAF+breast, and CRAF+breast (improved). The results obtained without considering backscattering were 2.4% lower than those obtained when considering body backscattering. This was because posterior adipose tissue existed in the model itself, which provided a portion of the backscatter dose for the glandular tissue. However, the presence of CRAF still resulted in more accurate simulation results. Compared with the combination method of the female body, which wasted considerable computational time in unnecessary organ voxels for particle transport [28], the improved method adopted in this study only retained the main body and organs that produced backscatter particles for glandular tissue. This resulted in a significant improvement (approximately 40 times) in computation speed while ensuring computational accuracy. This allowed us to simulate a large number of particles in a short period (several hours) and reduce statistical errors in the dose.

Normalized glandular dose - $D_{\text{gN}}^{\text{DBT}}$

3.3.1

Parameter dependence of breast characteristics and beam conditions on $D_{\text{gN}}^{\text{DBT}}$

The dependence of various model parameters on $D_{\text{gN}}^{\text{DBT}}$ must be investigated based on a detailed breast model. In the simulation, four breast models with different CBTs (4 cm and 5 cm) and $P_{\text{g}}^{\text{M}}$ (25% and 50%) were used for irradiation. We calculated the HVL for the energy spectra at different tube voltages and two target/filter combinations of Mo/Mo and W/Rh and found a good linear relationship between the HVL and $D_{\text{gN}}^{\text{DBT}}$. Linear regression analysis showed high consistency ($R^2>0.999$) and low residual error (within $1\times {10}^{-3}$) for Mo/Mo (Fig. 4(a)), and a good consistency ($R^2>0.95$) and a moderate residual error (within $4\times {10}^{-3}$) for W/Rh (Fig. 4(b)).

Fig. 4Full-screen View

$D_{\text{gN}}^{\text{DBT}}$ values calculated for the (a) Mo/Mo and (b) W/Rh Target/Filter combination, and their linear fit (continuous lines) with varying beam HVL

Figure 5 shows the relationship between $D_{\text{gN}}^{\text{DBT}}$ and the glandularity for different CBTs and the target/filter combination of W/Rh at 28 kVp. $D_{\text{gN}}^{\text{DBT}}$ values generally decreased as the $P_{\text{g}}^{\text{M}}$ concentration increased from 25 to 100%. However, we also observed an unexpected result: the model with 5% glandularity had a lower $D_{\text{gN}}^{\text{DBT}}$ than that with 25% $P_{\text{g}}^{\text{M}}$. This result contradicts the results obtained from literature based on simple breast models [9, 39]. The possible reasons for this discrepancy are discussed in the following section.

Fig. 5Full-screen View

Varying $D_{\text{gN}}^{\text{DBT}}$ values with $P_{\text{g}}^{\text{M}}$ for different CBTs and the target/filter combination of W/Rh at 28 kVp.

We examined the relationship between $D_{\text{gN}}^{\text{DBT}}$ and CBT at different tube voltages by setting the $P_{\text{g}}^{\text{M}}$ of the breast models to 25%. These data were fitted to the bi-exponential function proposed by Sobol and Wu (1997), which describes the relationship between $D_{\text{gN}}^{\text{DBT}}$ and CBT. The regression analysis revealed that this function accurately modeled our data (Figure 6) with a high goodness of fit ($R^2>0.99$) and low residual error (within $4\times {10}^{-3}$).

Fig. 6Full-screen View

(Color online) $D_{\text{gN}}^{\text{DBT}}$ values calculated for the W/Rh target/filter combination and bi-exponential fit (continuous lines) with varying CBTs

Glandular depth dose - $D_g^{dep}\left(z\right)$

We performed simulations on detailed and simple breast models with different $P_{\text{g}}^{\text{M}}$: 25%, 50%, 75%, and 100%. Each model had the same CBT value of 4 cm. The glandular tissue in the detailed models was not uniformly distributed along the depth direction, as shown in Fig. 7(a). We used a W/Rh 28 kV beam to irradiate these breast models and calculated the normalized $D_{\text{g}}^{\text{dep}}\left(z\right)$ for each model. Figure 7(b) compares the normalized $D_{\text{g}}^{\text{dep}}\left(z\right)$ curves for the detailed models and shows the relative errors of these four models during the MC simulation in the top-right corner. The linear attenuation coefficients of the models depended on the glandularity. As glandularity increased, so did the linear attenuation coefficient, because the density and equivalent attenuation coefficient of glandular tissue are larger than those of adipose tissue. Figure 7(c) shows a comparison of $E_{\text{gN}}\left(z\right)$ between the detailed and simple models. For the simple model, most of the deposited energy in the glandular tissue was concentrated in the upper surface of the central region, but for detailed models, the main energy-deposited area was concentrated between breast depths of approximately 1 cm and 2 cm. The greatest discrepancy in $E_{\text{gN}}\left(z\right)$ between the detailed and simple breast models was found near the breast surface and was primarily associated with depth variations of $f_{\text{g}}\left(z\right)$. When considering detailed breast models with varying glandularities, the main energy deposition area depicted in Fig. 7(c) was marginally shallower than the depth associated with the highest concentration of $f_{\text{g}}\left(z\right)$ because of the exponential attenuation of the X-ray beam.

Fig. 7Full-screen View

(Color online) (a) Glandular fraction along the z direction of detailed and simple breast models. (b) Normalized $D_{\text{g}}^{\text{dep}}\left(z\right)$ curves of detailed breast models with different $P_{\text{g}}^{\text{M}}$: 25%, 50%, 75%, and 100%. (c) $E_{\text{gN}}\left(z\right)$ of detailed and simple models.

Literature comparison

As demonstrated in Sect. 3.3.1, a strong linear correlation exists between $D_{\text{gN}}^{\text{DBT}}$ and the HVL. By fitting the available data, we calculated the $D_{\text{gN}}^{\text{DBT}}$ values for other HVLs and compared them with the $D_{\text{gN}}^{\text{DBT}}$ values obtained by multiplying with gcsT, as provided by Dance. The relative differences ∆ ranged from 18.6 to 51.0%, as summarized by CBTs and X-ray spectra. When the CBT and glandularity remained constant, ∆ increased with HVL because a higher HVL beam is more penetrating. Figure 8 illustrates the $D_{\text{gN}}^{\text{DBT}}$ results for various CBT and glandularities. As the CBT increased and glandularity decreased, the corresponding points progressively deviated from the identity line. This indicates that ∆ decreases with increasing glandularity and increases with increasing CBT. The highest ∆ value was observed for the breast model with 7-cm CBT and 5% glandularity.

Fig. 8Full-screen View

(Color online) Comparison of the estimated $D_{\text{gN}}^{\text{DBT}}$ values for different CBTs and glandularities under the target/filter combination W/Rh 28 kVp of Siemens Mammomat Inspiration. The data of each group with the same CBT correspond to 100%, 75%, 50%, 25%, and 5% $P_{\text{g}}^{\text{M}}$ from left to right. The identity line indicates equal $D_{\text{gN}}^{\text{DBT}}$ values for the detailed and simple breast models.