Three-dimensional (3D) digital breast tomosynthesis (DBT) is being considered as a prospective replacement for digital mammography (DM) as the imaging modality of choice in breast cancer screening programs. It is crucial, however, to recognize that the overarching objectives of breast cancer screening programs encompass not only elevating the cancer detection rate, which can be achieved using DBT, but also diminishing breast cancer-related mortality by identifying advanced cancers at earlier stages. Therefore, it becomes imperative to establish whether DBT systems have the capacity to supplant DM in the endeavor to identify aggressive cancers in their incipient stages. Traditional clinical trials designed for such investigations usually are very lengthy, expensive and may involve exposing women to extra radiation. Such trials can significantly hinder regulatory evaluation, delaying patient access to safe and effective advanced technologies. In contrast, the expeditious and cost-efficient alternative lies in the application of sophisticated computational techniques, referred to as in silico trials (Badano et al 2018, Bakic et al 2018) 1 , which serve as a pragmatic source of regulatory evidence. The goal of our work, described in this paper, is to model the growth of breast lesions for longitudinal in silico trials studying the effect of breast cancer progression on imaging performance.

To facilitate the modeling of cancer progression within computational trial frameworks, we developed a 3D voxelized breast lesion model capable of generating distinct growth stages. Camp et al (2023) summarizes several approaches that have been reported for simulating in silico versions of breast lesions. As shown in table 1, these approaches can be categorized as patient-derived and parameterized models. Patient-derived models, while offering a heightened sense of realism (Bliznakova et al 2019, Dukov et al 2019), by incorporating features from actual patient data, exhibit certain limitations, including restricted tumor model resolution, a limited array of lesion shapes, and the requisite approval from institutional review boards (IRBs) to acquire patient data. Conversely, parameterized mathematical models enable the generation of a broader spectrum of lesion shapes with variable resolution (Bliznakova et al 2003). Notably, none of the models cataloged in table 1 are grounded in biological phenomena. In contrast to previous endeavors, the research detailed in this article takes into account the biological factors that exert an influence on tumor growth and morphology, including interstitial tumor pressure, metabolite diffusion, and the surrounding tissue's stiffness.

Table 1. Methods for breast lesion models from the literature.

Patient-derived models Segmentation of 2D slices from DBT or breast computed tomography (bCT) volumes and stacked to create the final 3D tumor morphology (Bliznakova et al 2019, Dukov et al 2019) Mathematical models Simple shapes (e.g. ellipsoids/spheres) (Bliznakova et al 2003)Random walk model (Hintsala et al 2009)Diffusion limited aggregation model (Rashidnasab et al 2013)Gaussian sphere model with iterative fractal branching (Sisternes et al 2015)

Previous computational and theoretical efforts have attempted to predict pattern formation in nature. Among them, the reaction-diffusion model developed by Turing is one of the most popular analytical frameworks for the explanation of spatial biological pattern formation (Turing 1952, Kondo and Asal 1995). This foundational model has subsequently evolved into the advection-reaction-diffusion model, which incorporates the influence of surrounding interstitial tumor pressure on the flow rates of various metabolites and growth factors. This expanded model has been instrumental in simulating the growth of diverse tumor types (Chaplain 1996, Frieboes et al 2010, Tang et al 2014, Yan et al 2017, Lujan et al 2019, Sauer et al 2020), including aggressive glioblastomas (Frieboes et al 2010, Yan et al 2017). Additionally, Weis et al (2013) have pioneered a similar approach, utilizing it to simulate the growth of breast lesions in MRI studies. It is worth noting that while Weis et al's work explored the impact of the lesion's microenvironment on the diffusion of drugs and nutrients, it did not address the growth of anisotropic lesions, a notable limitation of prior investigations.

The main objective of developing the lesion growth model outlined in this study is to facilitate longitudinal research studies (Osoba et al 2000, Steele et al 2000, Langendijk et al 2001, Lago et al 2023) investigating the progression of cancer. We intend to integrate the growth model into an in silico rendition of the National Cancer Institute-funded Tomosynthesis Mammographic Imaging Screening Trial (TMIST) (Pisano 2018) to investigate whether DBT can reduce the incidence of advanced breast cancer (Lago et al 2023). Our model would be used to grow lesions at different stages into breast phantoms of varying densities. The results of such a study could help inform the choice of imaging modality for breast cancer screening programs which aim to detect cancers at earlier stages, thereby reducing the incidence of advanced breast cancer in the population.

The methodology introduced in this present article models the growth of anisotropic breast cancer lesions specific to their breast environments. Furthermore, we have devised a technique for incorporating spiculations at various stages of lesion development. Lastly, we present digital mammography (DM) images depicting the evolving lesions within anatomical breast digital models, utilizing the openly accessible Virtual Imaging Clinical Trials for Regulatory Evaluation (VICTRE) in silico imaging pipeline (Badano et al 2018) 1 . These images showcase several stages of tumor growth, shedding light on the intricate dynamics of this process.

2.1. Tumor dynamics

Tumor growth encompasses two distinct phases, namely the avascular and vascular phases. During the avascular phase, the tumor grows up to a few mm (Chaplain 1996), depending on the diffusion of vital nutrients into the tumor from the surrounding tissue. Due to rapid cell division within a confined region, cells located toward the center experience very high interstitial pressure. Vasculature and insufficient lymphatic systems leading to fluid accumulation also play a vital role in the increase in the interstitial pressure in tumors. The interstitial pressure attains its peak at the tumor's core and diminishes progressively towards the periphery. This gradient in interstitial pressure engenders an outward convection of essential nutrients, leading to oxygen (O2) deprivation in the core, resulting in necrosis. In response to hypoxia, tumor cells begin the process of secreting chemicals known as tumor angiogenesis factors (TAFs), inducing new blood vessels to sprout from existing ones. The new vessels gradually penetrate the tumor, providing an adequate blood and nutrient supply required to maintain its growth, now in the vascular phase.

2.2. Generic Growth Model by Tang et al (2014)

In this work, we have modeled the dynamics of tumor growth, as elucidated above, through the utilization of advection-reaction-diffusion equations, a framework adapted from the work of Tang et al (2014). The model commences with the computation of interstitial tumor pressure, which is the summation of pressures induced by cellular and vascular perfusion. For the sake of simplicity, we have omitted the consideration of the lymphatic system's contribution to tumor pressure. Furthermore, the pressure attributed to vascular perfusion is factored into the interstitial tumor pressure exclusively during the vascular phase, when tumors become palpable (Gavaghan et al 2002). Since our focus in this study is primarily on tumors in their initial growth stages, specifically the avascular phase, our model begins by calculating interstitial tumor pressure stemming from rapid cell proliferation within a confined space. This pressure is calculated as a Gaussian-like function, with the highest pressure being exerted at locations closer to the tumor cell, as expected. In other words, the cell- induced pressure P experienced by a tumor cell located at X0, due to N neighboring cells, located at Xi can be expressed as (Tang et al 2014):

Here, l(X0, Xi ) is the euclidean distance between two locations in the tumor, and pi is a constant related to interstitial tumor pressure. Please note that in the original publication by Tang et al pi is expressed as a function of the tumor cell density, θ. However, in the actual computational implementation, pi was defined as a constant. The value of was determined empirically (Tang et al 2014):

where θ is the tumor cell density calculated as the ratio of the number of grid points in the surroundings occupied by tumor cells to the total number of neighboring locations (26 in our implementation). The value of the constant σ0 (0.15) was determined empirically (Tang et al 2014).

This model focuses on O2 and carbon dioxide (CO2) as the main macromolecules that affect the growth of the tumor. Other nutrients such as glucose, amino acids, fatty acids, vitamins, and micronutrients also play an important role in the growth of tumors. However, modeling their transport and permeability across the cell's plasma membrane is a complex process (Tang et al 2014). So for simplicity, only the effect of oxygen O2 is considered in this model. Although O2 is the only nutrient considered, it is crucial for aerobic metabolism in cells. Oxygen plays an important role in the tumor growth process by aiding cell metabolism, angiogenesis and metastasis. The gradients of pressure affecting the flow rate of the macromolecules such as O2 and carbon dioxide (CO2) and their final concentrations are calculated using advection-reaction-diffusion equations, such as (Tang et al 2014):

where n (mol/m3) represents O2 concentration, Dn is the diffusion constant in m2/s and u is the flow rate, ρ0 is the rate of oxygen supply in mol/(m3s), Ri is the radius of a blood vessel and Pw is the pressure gradient through the vessel wall. Ai is the cell activity and λ0 with units of mol/(m3s) 2 is an empirically-determined constant. The flow rate u is calculated as (Tang et al 2014):

where k is the hydraulic conductivity.

The cell activity is calculated as (Tang et al 2014):

where w represents CO2 concentration.

The first and second terms on the right side of the equation (3) represent diffusion and advection, respectively, while the third term is a source of O2 and the fourth represents cell activity consuming O2. The concentration of CO2 also can be determined using similar equations, but with the cell activities acting as a source of CO2. Please note that equation (3) is equivalent to the combination of equations (6)–(8) in the original publication. All parameter values used in the growth model including diffusion coefficients, supply and consumption rates of the various metabolites, and initial values of the nutrient and waste concentrations are summarized in table 1 by Tang et al (2014). References have also been provided to justify some of these choices, while the others have been optimized to cover the transition of the tumor from its avascular to vascular phase with reasonable computational efficiency.

The framework developed so far by Tang et al (2014) is summarized as part of figure 1. It begins with the calculations of interstitial tumor pressure (equation (1)), used to calculate the metabolite O2 and CO2 flowrates (equation (4)) and finally their concentrations (equation (3)). Based on these metabolite concentrations and cell activity levels (equation (5)), tumor cells are classified as active, quiescent or necrotic. Active cells, localized at the periphery of the tumor, are primarily responsible for engaging in cell proliferation. In contrast, quiescent cells occupy the intermediate layer between the active cells and the necrotic core. When an active cell enters the proliferation algorithm, it is treated as though it occupies the center of a cube with 26 neighboring positions, each representing a potential and equally probable site for proliferation. Our algorithm is, thus far, grounded in the work developed by Tang et al. However, a significant enhancement introduced in this study relates to cell proliferation, which is now informed by the stiffness characteristics of the various tissues within the tumor's local microenvironment (described in section 2.3), as opposed to random cell assignment described by Tang.

Figure 1. Algorithm for anisotropic tumor growth based on pressure fields.

Download figure:

Standard image High-resolution image

The simulation time scale for tumor growth is recorded in cell life cycles (CLCs). One CLC corresponds to approximately 25 iterations of the code, after which cells acquire sufficient energy to undergo division based on oxygen (O2) and carbon dioxide (CO2) concentrations.

2.3.  In situ growth

Our approach relies on the fundamental assumption that there exists a certain degree of correlation between tumor morphology and the mechanical properties of the surrounding tissue types. We conducted our tumor growth simulations within volumes extracted from a breast model comprising adipose, glandular, and ligament tissues. Additionally, the breast volume incorporated other anatomical structures such as ducts, arteries, and veins. Our hypothesis posits that modulating the elasticity of these various anatomical structures exerts an influence on the evolving shape of the growing lesion. This hypothesis aligns with earlier findings (Helmlinger et al 1997, Stein et al 2007) which have demonstrated the sensitivity of tumors to their microenvironments, with their growth being hindered under conditions of elevated mechanical stress. Prior models depicting growth of breast lesions (Weis et al 2013) and gliomas (Garg and Miga 2008) have incorporated the effect of mechanical stress by coupling the stress factor to the cell diffusion coefficient through an exponential function. However, a notable limitation of the previous approach by Weis et al (2013), is that despite accounting for tissue mechanics, the growth model generated predominantly isotropic tumors with minimal perturbations at their edges. The sole instances of anisotropic lesions in Weis' work were derived from patient data. Furthermore, Weis et al's model was primarily designed to predict the number of tumor cells rather than their final morphology. Additionally, their model was restricted to describing 2D lesion growth, rendering it unsuitable for generating the 3D voxelized lesion models essential for conducting in silico trials.

In our work, a more direct approach was adopted by integrating the mechanical constraints imposed by the local tissue structures into the cell proliferation algorithm. Essentially, stiffer tissues inhibit tumor cells from deforming or infiltrating them, causing tumor cells to preferentially grow into softer regions within the breast. As illustrated in figure 1, the algorithm adjusts the value of k0 in direct proportion to the elastic constants characterizing the tissues surrounding the proliferating tumor cells. The local pressures are manipulated by k0, thus facilitating the correlation between tissue stiffness and tumor growth dynamics:

Here, is the total interstitial tumor pressure calculated at X0 using equations (1) and (2), and is the updated pressure modified by k0 at X0. Following the algorithm shown in figure 1, a cell proliferation probability was assigned to each of the neighboring 26 locations, such that locations that had softer tissues (and lower pressures) had higher probabilities. Finally, a location for cell proliferation was randomly sampled, and a new cell was added in that location.

The elastic modulus for adipose, glandular and ligament tissues (Gefen and Dilmoney 2007) and corresponding k0 values are shown in table 2. One of the most frequently-detected types of cancer is ductal carcinoma in situ (DCIS) (Sakorafas and Tsiotou 2000), which starts out from within the ducts (Coleman 2019). If left untreated, most cases develop into invasive breast cancer (IBC) (Risom et al 2021). DCIS originates from the epithelial cells in the ducts and sometimes grows along the ducts for a prolonged (up to several years) period of time, after which it invades the surrounding tissues (Erbas et al 2006). To simulate lesions representing DCIS, ductal fluid was assigned a k0 value of 1. Since fat may not be as elastic as ductal fluid, a k0 value of 10 was assigned to the voxels representing fat. So far we found no evidence to indicate that tumors originating in other tissues of the breast can penetrate and continue growing inside blood vessels. Although angiosarcomas originate in the blood vessels, this form of cancer only accounts for 1 % of all breast tumors (Arora et al 2014). Thus this type of cancer has been ignored for the current study. So in this work, other structures included in the breast model, such as arteries, veins and skin were considered to be statistically impenetrable, with a k0 value of 106.

Table 2. Pressure multiplicative factor based on Elastic Modulus for different tissue types.

Tissue typeElastic modulusPressure multiplication factor [kPa](k0)Fat0.5–2510Glandular7.5–6645Ligaments80 000–400 0001062.4. Local anatomy modeling

For our simulations, breast models developed by Graff (2016) were used to simulate lesion growth in situ. As compared to previous approaches described in literature (Bakic et al 2001, 2002, 2003, Bliznakova et al 2003, Bliznakova 2016), the breast model developed by Graff (2016) incorporates a comprehensive set of anatomical features such as vasculature, ductal network and a realistic breast shape, in addition to the adipose and glandular tissues. These features make the breast model by Graff et al an attractive option for our specific application.

To create anisotropic pressure maps based on the local anatomy, 1 cm ×1 cm ×1 cm volumes were extracted from fatty breast models (Graff 2016) generated using the VICTRE toolbox (Badano et al 2018, Sharma et al 2019)1. Locations for volume extraction were chosen randomly from a list of potential candidates generated for each breast model. These volumes were selected after confirming that they lay within the breast boundary and had no overlap with air, muscle, nipple or skin.

The extracted anatomies were 3D voxelized volumes, with each voxel representing a different tissue type. As illustrated in figure 2(a), a central slice from one of the extracted volumes shows adipose and glandular tissues represented by black and blue, while the white voxels indicate ligaments. Notably, the ligaments are interspersed with many black voxels of fat, akin to 'holes'. While the presence of such minute holes in the breast model may not affect the appearance of the much larger DM or DBT images, they can affect the shape of tumors growing in their vicinity. Since these holes are represented as voxels of fat, a significant number of tumor cells can penetrate the ligaments through these holes resulting in isotropic shapes and unrealistic morphology.

Figure 2. Postprocessing of extracted anatomy volumes to reduce the number of holes. (a) A central slice from a 1 cm ×1 cm ×1 cm anatomy volume extracted from a fatty breast model. The black regions are fat, blue is glandular tissues and white are ligaments. An example of large holes in the ligaments is encircled. (b) Flowchart for postprocessing of the extracted volumes to reduce holes and (c) Process for reducing the number of holes in the ligaments. Fat and glandular tissues are represented by black and dark blue. Ligaments are indicated in cyan.

Download figure:

Standard image High-resolution image

To address this issue, an algorithm was devised (see figure 2(b)) to control the penetrability of the ligaments. First, the ligaments were segmented and a MATLAB-based Gaussian smoothing function 3 was used to fill the holes and enhance the thickness of the ligaments. Subsequently, the thicker ligaments were reinserted into the original volume by voxel substitution. By adjusting the standard deviation of the Gaussian function 3 , the number of holes in the ligaments could be controlled, and consequently, their impact on tumor shape. An illustrative example of the process for reducing the number of holes in the ligaments is depicted in figure 2(c).

Once the processing of the extracted volume was completed, it was translated into a 3D map of corresponding to the voxel tissue types. Subsequently, this map was used to manipulate the local tumor pressures during the cell proliferation process, following the methodology in figure 1.

2.5. Addition of spiculations

One of the most common indicators of breast cancer discernible in DM or DBT images is the presence of a mass with thin, needle-like projections called spicules, radiating from its center (Seto and Mardiyana 2019). Spicules are formed either through the tumor's direct invasion into the surrounding tissue or as a consequence of a desmoplastic reaction with the neighboring tissues (Griff and Dershaw 2002). To date, there has been no documented report addressing the modeling of spicule growth predicated on biological processes.

The spicules were modeled following a technique by Elangovan et al (2016). In their research, the authors detailed a method for generating spiculated lesions, which involved creating a central mass through a diffusion-limited algorithm and then attaching a set of spicules to the surface of this mass generated using patient data. Nevertheless, instead of relying on patient data, we opted for the utilization of our implementation of de Sisternes' model (Sisternes et al 2015) to construct a 3D map of spicules, as depicted in figure 3. The de Sisternes' model employs a stochastic Gaussian random sphere model to simulate the central mass and an iterative fractal branching algorithm to introduce intricate spicule structures emanating from the central mass (Sisternes et al 2015). Employing de Sisternes' model, a 3D voxelized spiculated mass was generated, with distinct voxel values assigned to the central mass and spicules. By virtue of this differentiation in voxel values of the central mass (represented in blue in figure 3) and spicules (represented in yellow in figure 3), the spicules could be effectively separated to create a 3D map. The number of spicules, their lengths, and curvatures were controlled with input parameters to the de Sisternes' algorithm. Finally, spicules of varying lengths were attached to the mass corresponding to different time points from our lesion growth model. This approach proved to be effective for modeling both isotropic and the majority of anisotropic masses. However, it is worth noting that, in instances where the central mass displays pronounced anisotropy, certain gaps may become apparent, as indicated by the encircled areas highlighted in red in figure 3). To mitigate the occurrence of gaps, we strategically generated spiculated de Sisternes' masses by incorporating an extremely small central mass accompanied by a set of elongated spicules. This deliberate approach minimized the likelihood of gaps when superimposed with the masses generated through the growth algorithm while enhancing the overall coherence of the generated structures.

Figure 3. Process of generating spiculated lesions for isotropic and anisotropic central masses. The spicules are indicated in yellow while the central mass is blue.

Download figure:

Standard image High-resolution image 2.6. Imaging

The lesion models representing various stages of growth were generated by saving the lesion morphology after sufficient CLCs (15, 17 and 20). These lesion models were subsequently incorporated into the same fatty breast models using voxel replacement techniques. The DM projections were generated using the GPU-accelerated Monte Carlo (MC-GPU) x-ray transport simulation code (Badal et al 2020), developed as part of the VICTRE pipeline (Badano et al 2018, Sharma et al 2019) 1 . The MC-GPU codes replicated the technical specifications and imaging characteristics of the Siemens Mammomat Inspiration system, currently in clinical use for DM and DBT. The MC-GPU codes were specifically used to model a direct 200 μm amorphous selenium-based detector with 2816 × 3584, 85 μm pixels. The radiographic spectrum suitable for the fatty breast models was simulated using a tungsten anode with a peak tube voltage of 30 kVp, filtered with 50 μm rhodium and 1 mm beryllium. To generate the DM images for each breast model containing the lesions, 3.7 × 1011 x-ray photons were simulated and tracked, which took 38 minutes on a system with 2 GeForce GTX 1080 NVIDIA GPUs and 32 GB of memory. The lesions were realized as denser (up to 1.5 times) glandular tissue for visualization purposes. All post-processing of the raw images generated by the MC-GPU codes was carried out using ImageJ. The process involved an initial conversion of the images to 16-bit format, followed by contrast enhancement, using window level and contrast limited adaptive histogram equalization tools on ImageJ, to facilitate accurate visualization and analysis.

2.7. Validation

While rigorous validation of the growth model is important, for our specific applications, it may be sufficient to demonstrate that our model can generate realistic mammograms. To gain some insight into the realism of our model, we conducted a reader study in which we presented 25 experienced readers in breast imaging with 25 digital mammography (DM) regions of interest (ROIs) from the DDSM dataset (Bowyer et al 1996), which contained real manifestations of breast cancer, and 25 simulated DM ROIs with inserted lesions from our models. In this study, readers were asked to individually score each ROI on a scale of 1–8, where 1 represented the least realistic and 8 represented the most realistic. The primary objective of the study was to ascertain whether the readers assigned comparable scores to both the simulated and real ROIs.

3.1. Modeling local anatomy

Figure 4 provides a visual representation of the impact of ligament hole preprocessing. The first column on the left, shows slices extracted from volumes originating from a fatty breast model, which have undergone post-processing involving the application of a Gaussian filter-based convolution kernel with varying standard deviations 3 . In this representation, adipose and glandular tissues are depicted in black, while ligaments are shown in green. The lesions themselves are indicated in red. As anticipated, the number of holes in the ligaments decreases with an increase in the standard deviation of the smoothing function. This reduction in the number of holes enables tumor cells (depicted in red) to infiltrate the fatty regions (in black, see figure 4) of the breast. This results in a lack of expected growth inhibition due to the presence of ligaments, ultimately leading to isotropic growth, as demonstrated in the top row of figure 4. By further increasing the standard deviation of the Gaussian filter-based convolution kernel 3 , the ligaments were effectively thickened yielding a spectrum of breast lesion morphologies at various locations within the breast, as showcased in the second column of figure 4. The third column of the figure illustrates simulated digital mammograms (DMs) that contain these same lesions. Additional images demonstrating these effects at different locations within the breast can be found in the supplementary material.

Figure 4. Effect of postprocessing of extracted anatomy volumes on final lesion morphology. Column 1 shows slices of a volume extracted from a fatty breast model, postprocessed with a Gaussian filter based convolution kernel of varying standard deviation 3 . The number of holes in the ligaments decrease with an increase in standard deviation. Column 2 shows the 3D voxelized lesion alone while column 3 shows ROIs from simulated DMs containing the same lesion in column 2.

Download figure:

Standard image High-resolution image 3.2. Effect of local anatomy on morphology

Figure 5 presents a range of lesion morphologies generated using our model. Column 1 illustrates the various shapes of the lesions, while Column 2 provides additional context related to the local anatomy. Specifically, the green and blue voxels signify ligaments and ducts, respectively, while the lesions are represented in red. The first two rows illustrate the influence of ligaments on lesion morphology. The ligaments in the second row exhibit fewer holes compared to those in the first row, resulting in significantly constrained lesion growth. The third row showcases the growth of DCIS, where tumors originate within ducts, progress along them, and eventually extend into the surrounding fatty and glandular regions. In Column 3, one can observe ROIs extracted from simulated digital mammograms (DMs) generated through the VICTRE pipeline. The lesions, highlighted in red circles, were simulated with attenuation coefficients similar to those of normal glandular tissues, featuring a higher (1.3 times) tissue density.

Figure 5. Lesions and their appearance in renderings and in the DM images. Lesions (column 1) grown in situ. The second column includes local anatomy information. The third column are DM ROIs containing the lesions in the first column. All lesions were imaged as glandular tissues, but with higher density. For better visibility, the lesions in all rows were imaged with 1.3 times higher density, except for the one in the third row, which was imaged with 1.8 times higher density.

Download figure:

Standard image High-resolution image

The lower portion of figure 5 (rows 4–6) presents a variety of models that result in isotropic growth. In the first two rows, the lesions exhibit minimal interaction with the ligaments, leading to predominantly isotropic growth. The lesion in the fifth row encounters some glandular tissue, resulting in a more irregular shape compared to the growth seen in the fourth row. The sixth row illustrates a scenario in which the tumor encounters a ligament with a significant number of holes. These holes facilitate cell migration to the other side of the ligaments, ultimately culminating in isotropic growth. Figure 6 provides insight into the generation of both anisotropic and isotropic masses, as well as the subsequent attachment of a 3D spicule map. The third column showcases projections of the spiculated masses alone, offering a detailed view of the fine appearance of the spicules in the projections. The last column presents ROIs containing the spiculated masses, captured alongside their respective breast environments. To enhance the contrast and visibility of the spicules in the ROIs, the lesions were rendered as calcium oxalate masses. In the x-ray projections, areas with a higher density of lesion cells per unit volume appear brighter. The lesion core, densely packed with cells, appears the brightest in comparison to the periphery and the spicules.

Figure 6. Examples of spiculated masses. Column 1 shows the central mass alone, while the second column shows the mass overlapped with spicules. The third and fourth columns show projections of the lesions alone and with appropriate tissue background respectively.

Download figure:

Standard image High-resolution image

All ROIs depicted in figures 4 and 5 (Column 3) encompass the central masses shown in Column 1, with 3D spicule maps attached to them. However, since these lesions were imaged as dense glandular tissues, the fine spicules were no longer visible due to tissue overlap. Figure 7 illustrates ROIs extracted from digital mammography (DM) images of breast models hosting lesions at various locations and stages of growth. To enhance the visibility of the lesions, their attenuation was adjusted to be similar to that of glandular tissue but with a 1.5 times higher density. This adjustment allowed even the smallest lesion, at 10 CLCs, to be distinguishable from the normal anatomical structures.

Figure 7. Lesion growth from an imaging point of view. 250 × 250 pixel ROIs selected from simulated DM images generated using the VICTRE pipeline and breast model containing a disc-like lesion at different locations in the breast model and at different stages of growth (CLCs).

Download figure:

Standard image High-resolution image 3.3. Comparison with patient images

As an initial step towards validating our model, we conducted a comparison between our simulated digital mammograms (DMs) of fatty breast models hosting lesions grown at different locations and real mammograms containing masses with similar shapes. Figure 8 illustrates this comparison. The left column shows real mammograms, featuring malignant tumors characterized by ill-defined margins. The middle column (Column B) displays simulated DMs generated using the VICTRE pipeline, featuring tumors of similar shapes. To enhance visibility, the inset ROIs provide a magnified view of the same masses but with increased contrast. Column C offers further examples of ROIs extracted from both real and simulated DMs, encompassing microlobulated and ill-defined masses. Upon comparing the masses in the left and right columns, we observed that, qualitatively, our model exhibited the capacity to generate masses that bear similarity to those encountered in actual clinical cases.

Figure 8. Comparison with patient images. (a) Real mammograms from the DDSM dataset (Bowyer et al 1996) with a malignant tumor with ill defined margins, (b) simulated DM of a fatty breast model with similar lesions growing inside. The inset 250 × 250 pixel ROI shows a magnified version of the same mass, but with better contrast, (c) pairs of ROIs from real as well as simulated DMs containing malignant tumors with ill-defined and microlobulated margins.

Download figure:

Standard image High-resolution image 3.4. Reader study

The reader study revealed that the mean scores for the simulated lesions and the DDSM lesions were 4.5 (standard deviation = 1.7) and 5.4 (standard deviation = 1.6), respectively, with a median score of 5 for both sets. These results suggest that the simulated ROIs containing the mass models obtained using the method described in this article were perceived as similarly realistic as lesions collected from the DDSM dataset. In other words, readers were not able to distinguish between the two, providing initial evidence of the realism of our simulation results.

In this article, we introduced an innovative technique that integrates tumor growth dynamics with tissue mechanics to produce anisotropic and realis