Monte Carlo simulation of spatial frequency domain imaging for breast tumors during compression

Constance M. Robbins; Kuanren Qian; Yongjie Jessica Zhang; Jana M. Kainerstorfer

doi:10.1117/1.JBO.29.9.096001

14 September 2024 Monte Carlo simulation of spatial frequency domain imaging for breast tumors during compression

Constance M. Robbins, Kuanren Qian, Yongjie Jessica Zhang, Jana M. Kainerstorfer

Author Affiliations +

Journal of Biomedical Optics, Vol. 29, Issue 9, 096001 (September 2024). https://doi.org/10.1117/1.JBO.29.9.096001

Abstract

Significance

Near-infrared optical imaging methods have shown promise for monitoring response to neoadjuvant chemotherapy (NAC) for breast cancer, with endogenous contrast coming from oxy- and deoxyhemoglobin. Spatial frequency domain imaging (SFDI) could be used to detect this contrast in a low-cost and portable format, but it has limited imaging depth. It is possible that local tissue compression could be used to reduce the effective tumor depth.

Aim

To evaluate the potential of SFDI for therapy response prediction, we aim to predict how changes to tumor size, stiffness, and hemoglobin concentration would be reflected in contrast measured by SFDI under tissue compression.

Approach

Finite element analysis of compression on an inclusion-containing soft material is combined with Monte Carlo simulation to predict the measured optical contrast.

Results

When the effect of compression on blood volume is not considered, contrast gain from compression increases with the size and stiffness of the inclusion and decreases with the inclusion depth. With a model of reduction of blood volume from compression, compression reduces imaging contrast, an effect that is greater for larger inclusions and stiffer inclusions at shallower depths.

Conclusions

This computational modeling study represents a first step toward tracking tumor changes induced by NAC using SFDI and local compression.

1. Introduction

Neoadjuvant chemotherapy (NAC) is commonly prescribed for locally advanced breast cancer, often shrinking the tumor prior to surgery and thereby permitting more breast-conserving surgeries or rendering operable some tumors initially too large to be operable. Pathologic complete response (assessed histologically after surgery) is strongly associated with prognosis, but methods to predict this response early in the course of therapy are lacking. Structural imaging modalities and clinical palpation can assess changes in tumor size but can fail to distinguish between fibrosis and residual disease.¹^,² More success in outcome prediction has been demonstrated for dynamic contrast MRI,³^,⁴ at the disadvantage of high cost (limiting frequent measurement) and the need for exogenous contrast agents.

Near-infrared (NIR) optical imaging has emerged as a promising method for NAC monitoring due to its relatively low cost, non-invasiveness, and sensitivity to endogenous contrast from hemoglobin concentrations.⁴^–¹⁶ Increased total hemoglobin (tHb) relative to the background is thought to be related to angiogenesis within the tumor, and changes in vascularization in response to NAC are shown to occur earlier than structural changes,⁴^,⁵ including the “oxy-hemoglobin flare” observed as early as 1 day after the start of NAC.⁸^,¹¹ Over the course of therapy, multiple groups have reported a decrease in tHb in responding tumors but not in non-responding tumors¹⁰^,¹⁷^,¹⁸ or a greater tHb decrease in responding tumors.¹²^,¹⁹

Measuring the reaction of the breast to perturbations such as compression has also been used to provide additional contrast beyond that of baseline optical properties. A transient reduction in tHb of both tumor and healthy breast has been well documented in response to compression,²⁰^–²³ with the change in tHb being related to pressure distribution within the breast.²⁴ Carp et al.²³ reported a gradual recovery of tHb in normal tissue while displacement is held constant (and force is reduced due to viscoelastic relaxation) while tHb of the tumor region remains persistently low. This hemodynamic response has also been demonstrated to normalize in patients who achieve response to NAC.⁷

To elicit this type of hemodynamic response as well as to longitudinally track tumor and background optical properties in a low-cost format, we are pursuing a handheld device for breast imaging with localized compression based on spatial frequency domain imaging (SFDI), a diffuse optical imaging method that uses structured illumination to quantify tissue absorption and reduced scattering coefficients.²⁵^,²⁶ Because of its ability to quantify optical properties at a low cost and small form factor, the technique holds the potential for frequent monitoring during therapy, including in the patient’s home. As a reflectance-based 2D projection imaging method, SFDI suffers from partial volume effects, with superficial tissue contributing more strongly to the imaging result than deeper tissue.²⁷ Previously, we showed in polydimethylsiloxane (PDMS) phantoms that compression can be used to reduce the depth of stiff highly absorbing inclusions and improve optical contrast measured with SFDI.²⁸

We demonstrated that despite its high sensitivity to superficial tissue, SFDI can be used to detect tumor-mimicking inclusions at an initial depth comparable to palpable breast lesions and that localized compression can be used to increase the measured contrast. However, if longitudinal changes in contrast obtained with SFDI and local compression are to have utility as a biomarker for monitoring NAC, it is necessary to know what changes in imaging contrast would result from likely changes to tumor structure and function such as size,²⁹^,³⁰ stiffness,³¹^,³² and tHb concentration.¹⁰^,¹²^,¹⁷^–¹⁹ We leverage a combination of mechanical and optical simulations to explore the relationship between these tumor parameters and optical contrast under compression.

One method for Monte Carlo (MC) simulation of SFDI involves the acquisition of a spatially resolved simulation of point source illumination (in the spatial domain) and the application of the discrete Hankel transform to convert it into the spatial frequency domain. By contrast, the Gardner method avoids the need for discrete transforms and performs the simulation natively in the spatial frequency domain by calculating a frequency-dependent photon weight.³³ Both methods provide diffuse reflectance over a range of spatial frequencies and are suitable for homogenous or layered input media. However, neither can generate spatially resolved optical property maps of a complex 3D geometry.

To simulate the result of SFDI on an embedded highly absorbing inclusion (mimicking breast tumor), we employed MC to directly simulate the projection of a sinusoidal modulated light source at three phase shifts. The results are demodulated to obtain 2D diffuse reflectance maps, which are converted to optical property maps using a two-frequency lookup table (LUT) in the conventional manner for SFDI. In this method, diffuse reflectance information is limited only to those spatial frequencies that are simulated individually at three phases each, representing a considerable simulation run time. In addition, many photons must be launched for each MC simulation to limit noise in the detected reflectance with too few resulting in error in demodulation of the sinusoidal patterns. For this reason, only a single spatial frequency, $0.1 {mm}^{- 1}$ , is simulated at three phases, from which DC reflectance of $0 {mm}^{- 1}$ is also extracted according to the demodulation equations given in Gioux et al.²⁶ This spatial frequency pair was chosen to be consistent with that used in our group’s handheld SFD breast imager³⁴ and has been established to be sufficient for decoupling absorption ( $μ_{a}$ ) and reduced scattering coefficients ( $μ_{s}^{'}$ ).³⁵

To predict the final imaging result, we must model not only the 2D SFDI reflectance image that would arise from each input but also the change in tumor shape and depth from compression and the effect of compression on final optical properties in the tissue input. Finite element analysis (FEA) has been adopted together with image processing techniques to study tumor growth,³⁶ tracking,³⁷ and treatment.³⁸^,³⁹ In this study, FEA via ANSYS mechanical is used to predict the change in tissue geometry resulting from compression. Subsequent spatially resolved MC simulations are performed, both with and without modeling tissue blanching due to compression. Tissue blanching is modeled as occurring in all tissue, termed the “full blanching” condition, or as occurring in the background but not inclusion, the “background blanching” condition, as the internal stresses within the enclosed environment of tumors may prevent blood volume reduction from occurring in the same manner as healthy tissue.

2. Methods

The effect of tissue compression on the output of SFDI was examined across various tumor sizes, depths, and stiffness ratios (SRs). A breast tumor and surrounding healthy tissue were modeled as a stiff, spherical inclusion embedded in a softer cylindrical matrix, with optical properties based on typical tumor and healthy breast values. The deformation of this model under compression was simulated using ANSYS Mechanical, with new optical properties calculated based on the stress output of the mechanical simulation. MC SFDI simulations were then performed using Monte Carlo eXtreme (MCX).⁴⁰ This workflow of mechanical and optical simulations is illustrated in Fig. 1(a).

Fig. 1

(a) Summary of simulations performed. (b) Following mechanical simulations, tHb concentrations are adjusted on a voxel-by-voxel basis based on hydrostatic stress experienced by the model. The three models for the dependence of hemoglobin concentration on stress are illustrated here. In the “no blanching” model, tHb is held constant in both the inclusion and background regions. In the “background blanching” model, tHb of inclusion remains constant while tHb of the background region decreases linearly before plateauing at 20 mmHg. In the “full blanching” model, both inclusion and background regions experience a decrease in tHb with experienced pressure.

MC simulation is performed four times for each tissue input: one on the original geometry before compression and three on the compressed geometry corresponding to three models of the change in tHb concentration with pressure. In the “no blanching” condition, optical properties for inclusion and background remain equal to pre-compression values. In the “background blanching” and “full blanching” conditions, new optical properties are calculated on a voxel-by-voxel basis based on the stress distribution in the model for the background only and for both regions, respectively. The dependence of tHb on stress experienced is shown in Fig. 1(b). Results of the MC simulations were compared to the initial inclusion radius, $r_{0}$ , and depth, $d_{0}$ , and to the mid-compression surface radius of curvature, $r$ , and mid-compression depth, $d$ .

2.1.

Mechanical Simulation

For all mechanical simulations, a Young’s modulus of 4 kPa was used for the background material, based on results reported by Samani et al.⁴¹ for healthy breast. In that study, a Poisson’s ratio value of 0.495 was used in their FEA simulations, and the same value was used here. To mimic the stiffness contrast of breast tumors to healthy breast,⁴¹ the inclusions were given elastic moduli of 4, 16, 28, and 40 kPa, resulting in SRs of 1, 4, 7, and 10. The mechanical model consisted of a cylinder of background material of height 60 mm and radius 30 mm with the spherical inclusion located along the center axis of the cylinder, a variable distance from the top surface of the cylinder. Inclusion depths simulated were 3, 5, 7, 10, and 15 mm (measured from cylinder surface to inclusion surface) and inclusion radii simulated were 7.5, 12.5, 15, and 17.5 mm. A fixed boundary condition was applied to the cylinder bottom (opposite of the face to be imaged). A rigid plate (modeled as structural steel: Young’s modulus of $2 \times 10^{8} kPa$ , Poisson’s ratio of 0.3, and density of $7850 kg \cdot m^{- 3}$ ) was placed on the top surface of the cylinder (face to be imaged) with a frictional boundary condition with the cylinder and a friction coefficient of 0.1. This friction coefficient was intended to approximate the friction between the glass imaging head of the handheld SFDI device and the patient’s skin and was assumed not to differ greatly from the friction of lubricated glass on glass.⁴² To reduce the computational load, symmetry was used to reduce the size of the model, as is common in the field.⁴³ Only a 1/8 slice was simulated instead of the entire cylinder, with a frictionless support boundary condition applied to the faces cut through the cylinder. An upward displacement of 12.5 mm was applied to the bottom surface of the cylinder. These simulations took between 10 and 35 min each, for a total computation time of $\sim 24 h$ .

2.2.

Processing of FEA Output

For all FEA outputs, as well as the undeformed input meshes, the 1/8 slice was mirrored to form a full cylinder. The boundary surface between the inclusion and background nodes was extracted, as was the top surface of the background cylinder. The “iso2mesh” MATLAB toolbox was used to create a labeled volumetric image from these boundaries. The volumetric image had dimensions of $80 \times 80 mm$ for the top face and a depth of 62.5 mm ( $192 \times 192 \times 150 voxels$ ). Background-labeled voxels were extended out to the edges of the rectangular slab regardless of the outer surface of the FEA cylinder. The inclusion depth was extracted from the volumetric images for both compressed and uncompressed (UC) tissue inputs as the voxelization of the image caused depth to vary slightly from the inclusion depth of the FEA mesh. Similarly, the radius of curvature of the top surface of the inclusion was extracted from the volumetric images. A cross-section of the image through the center of the inclusion (plane orthogonal to the imaging surface) was taken, and the boundary pixels were identified with the MATLAB’s “bwboundaries” function. A circle was fit limited to those points located on the upper one-third of the inclusion.⁴⁴

Initial (pre-compression) tHb concentrations, tHb, for the main simulations were defined to be $50.6 μ M$ for tumor voxels and $17.6 μ M$ for background voxels, within the ranges reported by Grosenick et al.⁵ Simulations were also performed for tumor tHb of $63.3 μ M$ (a 25% increase) and background voxels remaining at $17.6 μ M$ , to verify the effect of increasing intrinsic tumor absorption contrast. The absorption coefficient at 805 nm (close to the isosbestic point) was calculated for each based on the extinction spectrum of hemoglobin⁴⁵ and assuming 75% tissue oxygen saturation. Though other chromophores are present in the tissue, including water, collagen, and lipid, their contribution in this region of NIR is small,⁴⁶^,⁴⁷ and only light absorption for tHb was considered. The resulting $μ_{a}$ was $0.0329 {cm}^{- 1}$ for background, $0.0948 {cm}^{- 1}$ for tumor, and $0.1185 {cm}^{- 1}$ for tumors with increased tHb. For all simulations, $μ_{s}^{'}$ was defined to be $10.37 {cm}^{- 1}$ for tumor voxels and $8.72 {cm}^{- 1}$ for background voxels.⁵ MC simulations of UC models were run using the specified baseline values for inclusion and background voxels. MC simulations were also run on all compressed geometries maintaining these initial optical properties for tumor and background voxels, with these results deemed the “no blanching” condition. The purpose of this condition was to model the effect of compression on change in optical contrast that arises from the change of shape and depth of a highly absorbing inclusion.

MC simulations were then run for the two tissue blanching conditions, in which the change in tHb concentration of the tissue is modeled based on the spatially varying stress experienced by the tissue input. FEA output stress was also converted to a volumetric image of the same dimensions as the tissue input. Hydrostatic stress was obtained by averaging the three normal stresses at each node. Hydrostatic stress at the center of each voxel was interpolated via MATLAB’s “scatteredinterpolant” function. Because the stress was undefined outside the cylinder boundary, stress located $\sim 1.5 mm$ inside the cylinder boundary was extrapolated outward to the edge of the volumetric image. Stress in the bottom 15 mm of the image was defined to be zero. As these voxels are far from the surface to be imaged, they contribute minimally to the diffuse reflectance. Hydrostatic stress values were converted to mmHg for ease of comparison with typical pressure in blood vessels.

A piecewise function was defined to model tissue blanching in response to pressure. Blood volume in the venous compartment (including veins and venules) was assumed to comprise 76.2% of systemic blood volume and to contain a maximum blood pressure of 20 mmHg.⁴⁸ For pressures between 0 and 20 mmHg, blood volume in the venous compartment was modeled to decrease linearly. Above 20 mmHg, the venous compartment was considered fully collapsed due to the lower pressure within veins,⁴⁸ and thus, no further decrease in blood volume occurred after 20 mmHg. For pressures below zero, the venous compartment maintains its original volume. As the compression pressures studied were below arterial pressure, no effect on the arteries was modeled.

Mid-compression tHb was calculated on a voxel-by-voxel basis using the function defined above and the corresponding hydrostatic stress matrix, and $μ_{a}$ at 805 nm (close to the isosbestic point) was calculated based on the extinction spectrum of hemoglobin.⁴⁵ For the “background blanching” condition, new $μ_{a}$ was calculated only for background voxels with the inclusion voxel optical properties remaining at pre-compression values. In the “full blanching condition, new $μ_{a}$ was calculated for all voxels (both inclusion and background). The effect of compression on $μ_{s}^{'}$ was not modeled.

2.3.

Implementation of Monte Carlo Simulations

For each input, three MC simulations were performed corresponding to three evenly spaced phase shifts of the spatially modulated light pattern, with all other simulation parameters held constant. Outputs from the three simulations were demodulated on a pixel-by-pixel basis in a fashion similar to conventional SFDI and consistent with our group’s handheld SFD breast imager. In our in vivo imaging of breast hemodynamics, projection of only one spatial frequency, from which reflectance at $0 {mm}^{- 1}$ is also extracted, allows fast image acquisition and fast processing of results. Likewise, simulation of one non-zero spatial frequency provides the same benefit here, while allowing the most consistency of methods with our in vivo SFDI results. The pair of 0 and $0.1 {mm}^{- 1}$ was chosen for both as it has been established to be sufficient for decoupling $μ_{a}$ and $μ_{s}^{'}$ .³⁵ Although lower spatial frequencies have been shown to provide higher sensitivity to deeper tissue,²⁷ lower frequencies do not allow for as efficient separation of $μ_{a}$ and $μ_{s}^{'}$ when paired with $0 {mm}^{- 1}$ .³⁵

The sinusoidal spatially modulated source was defined using the MCX source type “fourier” with spatial frequency $f_{x} = 0.1 {mm}^{- 1}$ and phase shifts of 0, $2 π / 3$ , and $4 π / 3$ . Cyclic boundary conditions (photons exiting from one face renter from the opposite face) were used for the side faces to mimic an infinitely tiled medium. Default boundary conditions were used for the remaining two, with photons counted as diffuse reflectance when escaping through the source-incident face and lost when exiting through the opposite face. With a total depth of 60 mm, a negligible amount of photons are lost through the bottom surface. A refractive index, $n$ , of 1.4 was used for background and inclusion material as this is a commonly used estimate for biological tissue.⁴⁹^,⁵⁰ A refractive index value of 1 was used for the region outside the sample, representing air. For the anisotropy factor, $g$ , a value of 0.8 was used, within the range of values reported in the literature.²⁷^,⁵¹^,⁵²

The raw output flux at each pixel of the tissue/air interface $(x, y)$ was demodulated to yield 2D maps of AC and DC flux. Equations of demodulation are given as:

Eq. (1)

φ_{DC} = (φ_{1} + φ_{2} + φ_{3}) / 3,

Eq. (2)

φ_{AC} = \frac{\sqrt{2}}{3} \cdot \sqrt{(φ_{1} - φ_{2})^{2} + (φ_{2} - φ_{3})^{2} + (φ_{3} - φ_{2})^{2}},

where

φ_{1}

,

φ_{2}

, and

φ_{3}

are the raw output flux of the three source phase shifts. This demodulation is performed in the same manner as the demodulation of intensities

I_{1}

,

I_{2}

, and

I_{3}

in conventional SFDI processing.⁵³ The distinction is that those intensities are of arbitrary units and must be compared with intensities obtained from a calibration phantom to obtain

R_{d}

, whereas these simulation outputs represented photon flux with units of

\frac{1}{{mm}^{2} \cdot s}

, and for that reason, they are denoted

φ

rather than

I

.

φ_{DC}

and

φ_{AC}

were then multiplied by timestep and area of the source-incident surface to produce unitless diffuse reflectances,

R_{d, DC} (x, y)

and

R_{d, AC} (x, y)

. Each pixel value of the diffuse reflectance maps is within zero and one and represents the fraction of light that is diffusely reflected at that location.

An MC LUT was generated by running simulations at three phases of the sinusoidal light input with $f_{x} = 0.1 {mm}^{- 1}$ for every combination of $μ_{a}$ from 0.01 to $0.28 {cm}^{- 1}$ (increments of $0.01 {cm}^{- 1}$ ) and $μ_{s}^{'}$ from 4.0 to $20 {cm}^{- 1}$ (increments of $1.0 {cm}^{- 1}$ ). For the LUT simulations, a homogenous slab of $24 \times 24 \times 150 voxels$ ( $10 \times 10 \times 62.5 mm$ ) was used. A “cyclic” boundary condition was defined for the four edge faces such that when a photon escaped from one side, it would re-enter the slab from the opposite face. This boundary condition allowed the simulation to mimic an infinitely wide slab with a relatively short computational time. A “total absorption” boundary condition was used for the bottom face as it was determined that a negligible number of photons reached the bottom surface of the slab. For each pair of optical properties, $1 \times 10^{7} photons$ were launched for each of the three source phases. $R_{d, DC}$ and $R_{d, AC}$ were averaged across the entire image and stored for each optical property combination. This quantity of photons launche was determined to be sufficient from visual inspection of the isolines of the resulting LUT, which were observed to be smooth.

MC simulations of the compressed and UC tissue geometries differed from the LUT simulations only in the tissue input and the number of photons launched. $1.4 \times 10^{8} photons$ were simulated per phase. More than an order of magnitude greater number of photons was used due to the larger tissue input and because spatial maps of optical properties were to be calculated, necessitating a lower signal-to-noise ratio than the LUT simulations in which $R_{d, DC}$ and $R_{d, AC}$ were averaged over the image. The tissue input definition was $192 \times 192 \times 150 voxels$ ( $80 \times 80 \times 62.5 mm$ ) and had continuously varying optical properties as determined above. All MC simulations of the UC tissue input as well as the “no blanching” and “full blanching” conditions took $\sim 33 h$ to complete. Simulations for the “background blanching” condition were performed later with higher computational power and required $\sim 2 h$ , collectively. Spatially varying optical properties were calculated from diffuse reflectance pixel-by-pixel using the LUT created above and MATLAB’s “scatteredinterpolant” function, yielding $μ_{a} (x, y)$ and $μ_{s}^{'} (x, y)$ .

2.4.

Processing of MCX Output

Optical contrast was quantified by subtracting the baseline value from the $μ_{a}$ map and summing the resulting pixel values over the area formed by projecting the inclusion boundary onto the surface. Though this sum would result in units of ${cm}^{- 1}$ , optical contrast was instead treated as having arbitrary units. For all conditions, the baseline $μ_{a}$ value was extracted from the final $μ_{a}$ map by averaging all pixels located at the image margins, at least 45 mm from the inclusion epicenter, at which location the inclusion contributes minimally to the result. This region of interest (ROI) was determined empirically by comparing the image margins of the UC simulation with the highest contrast and a case with no detectable contrast (the smallest inclusion at the deepest depth), and the image margins were found not to differ.

Though contrast calculated with the above method was used for most analyses, the tumor/normal (T/N) ratio was also calculated as this is a common measure in breast cancer optical imaging. In addition, a threshold of 5% difference from baseline (T/N ratio below 0.95 or above 1.05) was used to determine which results were considered detectable, as this level agreed with a visual inspection of the $μ_{a}$ maps.

3. Results

3.1.

Results without Change in tHb (“No Blanching” Condition)

MC-derived $μ_{a}$ maps of UC geometry show the expected trends with depth and radius: larger and shallower inclusions having the highest contrast over the background. These maps are shown in Fig. 2, with the highest contrast in $μ_{a}$ evident in the upper left corner ( $r_{0} = 17.5 mm$ , $d_{0} = 3 mm$ ). Moving to the right (corresponding to greater inclusion depth) or down (corresponding to smaller inclusion radius) in Fig. 2 both yield lower levels of contrast in $μ_{a}$ . All figures refer to the results using $50.6 μ M$ tHb for tumor and $17.6 μ M$ tHb for healthy breast.

Fig. 2

Summary of MC simulation results for UC tissue input. Each square corresponds to the $μ_{a}$ output for the inclusion radius, $r_{0}$ , and depth, $d_{0}$ , indicated. Dotted circles represent the projection of the inclusion outline onto the imaging surface.

When optical contrast is computed by baseline subtraction and summation of $μ_{a}$ over the area, these relationships can be quantified. When depth is held constant, contrast is found to have a linear relationship with inclusion radius, shown in Fig. 3(a). Equations of best fit of the form $C = a \cdot r_{0} + b$ are obtained for these linear relationships, where $C$ represents contrast and $r_{0}$ represents the initial radius. Coefficients are given in Table 1.

Fig. 3

Relationship of simulation output contrast to initial inclusion radius, $r_{0}$ , and depth, $d_{0}$ . (a) When depth is held constant, contrast is found to increase linearly with increasing radius. (b) When the inclusion radius is held constant, contrast is found to decrease exponentially with increasing depth. Plot colors adapted from Ref. 54.

Table 1

Summary of the relationship between contrast and initial radius when depth is held constant.

Initial depth (d0) (mm)	a	b
2.9	0.407	−0.598
4.6	0.267	−0.110
6.7	0.178	0.183
9.6	0.099	0.433
14.6	0.037	0.581

Equations of best fit are of the form C=a·r0+b.

When the radius is held constant, contrast is found to have an exponential relationship with inclusion depth, shown in Fig. 3(b). Equations of best fit of the form $C = a \cdot e^{b \cdot d 0}$ are obtained for these exponential relationships, where $C$ represents contrast and $d_{0}$ represents initial depth. Coefficients are given in Table 2.

Table 2

Summary of the relationship between contrast and initial depth when the radius is held constant.

Initial radius (r0) (mm)	a	b
7.5	3.11	−0.099
12.5	6.53	−0.142
15	8.50	−0.158
17.5	10.1	−0.162

Equations of best fit are of the form C= a · eb·d0.

This exponential decrease with depth is consistent with previous findings for PDMS phantoms with stiff, highly absorbing inclusions.²⁸

The effect of compression in all cases is to decrease the final depth compared to the initial depth ( $d < d_{0}$ ) and to increase the top surface radius of curvature compared to the initial radius ( $r > r_{0}$ ). Thus, when the optical properties of the inclusion and background are held constant with compression, optical contrast is gained from compression in all but two cases. The two cases that did not experience an increase in contrast consisted of the smallest and deepest inclusions with negligible contrast (T/N ratio between 1.0 and 1.05) both before and after compression. A representative example of the effect of compression on the initial geometry and the resulting change in MC-derived $μ_{a}$ is shown in Fig. 4. The effect of compression on the inclusion depth and radius of curvature is shown in the cross-sections of Fig. 4(a), and the resulting simulated $μ_{a}$ maps are shown in Fig. 4(b), with increased contrast evident for the compressed simulation. In this example, the initial radius was 17.5 mm, and the initial depth was 6.67 mm, with $SR = 2$ . After compression, the radius of curvature increased to 21.6 mm, and depth decreased to 4.58 mm. Background voxels were $μ_{a} = 0.0329 {cm}^{- 1}$ , and inclusion voxels were $μ_{a} = 0.0948 {cm}^{- 1}$

Fig. 4

Compression increases the inclusion radius of curvature and decreases depth, resulting in increased contrast. (a) Illustration of the effect of compression on input geometry. The black and white striped area indicates the phantom height reduction from compression. (b) Results of MC simulations on UC and compressed conditions. Dotted lines mark the outline of the inclusion projected onto the phantom surface.

Just as contrast before compression bore a strong relationship to $r_{0}$ and $d_{0}$ , contrast during compression is strongly related to final geometric parameters, $r$ and $d$ . Though the value of these properties during compression is not expected to be known in a potential clinical application of this work, it is of interest to quantify how the values of $r$ and $d$ affect optical contrast in this study to characterize the effect of compression-induced changes in shape and depth on optical contrast. Because $r$ and $d$ do not occur at discrete intervals as $r_{0}$ and $d_{0}$ do, all points were plotted on a set of 3D axes (with contrast on the $z$ -axis) instead of examining the relationship when each parameter is held constant. When contrast is plotted as a function of $r$ and $d$ , all points are found to lie within the same surface in 3D space, as shown in Fig. 5.

Fig. 5

Relationship of the final radius of curvature, $r$ , final depth, $d$ , and contrast for the solid phantom condition. Green squares indicate simulations of the 20 UC geometries ( $r = r_{0}$ and $d = d_{0}$ ), and black circles represent the 80 compressed simulations (20 inputs at 4 SR values).

The equation

Eq. (3)

z = (c_{1} + c_{2} y) \cdot e^{c_{3} x} + c_{4},

describing this surface, was obtained via least squares fitting (MATLAB’s “fit” function), where

y

is the radius,

x

is the depth,

z

is the contrast, and

c_{1}

,

c_{2}

,

c_{3}

, and

c_{4}

are constants obtained by the fit. This form was chosen to account for the previously observed linear relationship with radius when depth is held constant and the exponential relationship with depth when radius is held constant. The equation of the best-fit surface was

z = (0.802 + 3.694 r) \cdot e^{- 0.2259 d} 0.8871

. This surface describes a high degree of variance in contrast, with

R^{2} = 0.992

. This goodness of fit is also apparent from the visual inspection of Fig. 5. Thus, the SR is not necessary for the prediction of final contrast if the final depth and radius of curvature are known. This result also demonstrates that the top surface radius of curvature is a sufficient description of the shape and size of a deformed spherical inclusion. For this type of imaging, a flattened ellipsoidal inclusion mimics a larger spherical inclusion as only the top surface contributes significantly to optical contrast.

We have demonstrated that final contrast is determined only from final geometry parameters $r$ and $d$ . However, in a potential clinical application of this work, mid-compression $r$ and $d$ are not expected to be known, though it is of interest to predict contrast change (difference between final and initial contrast) as a function of the input parameters $r_{0}$ , $d_{0}$ , and SR. All simulation cases in the “no blanching” condition exhibited an increase in contrast induced by compression, with the exception of two cases of small and deep inclusions, which exhibited negligibly low contrast values both before and during compression. These results are visualized below in Sec. 3.3 along with the two tissue blanching conditions.

When simulations are repeated with tHb concentration in the tumor voxels increased from 50.6 to $63.3 μ M$ tHb, the same trends as above are found. The initial radius still exerts the most influence on contrast change, followed by initial depth and a weak influence from the SR. When the $63.3 μ M$ simulations are compared directly with the $50.6 μ M$ simulations, those with greater hemoglobin concentration in the tumor voxels experience greater contrast gain from compression, other variables being equal (not shown). This finding is consistent with the intuition that, as tumor voxels are brought close to the surface via compression, more highly absorbing tumor voxels will cause a greater gain in imaging contrast.

3.2.

Results with Change in tHb (“Background Blanching” and “Full Blanching” Conditions)

In simulating tissue blanching with compression, the same compressed geometry is used as in the “no blanching” condition, with final optical properties calculated from hydrostatic stress according to the relation described in “processing of FEA output” above. Again, all figures refer to the results using $50.6 μ M$ tHb for tumor and $17.6 μ M$ tHb for healthy breast. Representative examples of the effect of compression in the “background blanching” and “full blanching” condition and the resulting changes to the MC-derived $μ_{a}$ contrast are shown in Fig. 6. Inclusion $r_{0}$ , $d_{0}$ , and SR are the same as in Fig. 4. Figure 6(a) shows the effect of compression on both the inclusion geometry cross-section and the optical property distribution of the inclusion and background.

Fig. 6

Summary of the two tissue blanching models, in which compression decreases $μ_{a}$ of the inclusion, resulting in a reduction of optical contrast for some conditions. (a) Illustration of the effect of compression on the input model with simulated tissue blanching. Color represents simulation input $μ_{a}$ . The black and white striped area indicates the phantom height reduction from compression. (b) Results of MC simulations on UC and compressed conditions. Dotted lines mark the outline of the inclusion projected onto the phantom surface.

In the solid phantom condition, optical contrast was found to be dependent on the final geometry alone. Although the SR in that condition influenced the shape change of the inclusion with compression (and thus influenced $r$ and $d$ ), its value was not independently associated with final optical contrast. Conversely, in the “tissue blanching” condition, the SR strongly influenced optical contrast through its action on tHb reduction, and $r$ and $d$ alone described a comparatively small portion of the variance in contrast. As was done for the solid phantom condition, a surface of best fit of the form given in Eq. (3) was fit for the prediction of contrast from $r$ and $d$ . The equation of the best-fit surface was $z = (0.160 - 1.738 r) \cdot e^{- 0.1815 d} - 1.722$ . This surface was a poor fit, failing to describe much of the variation in contrast, with $R^{2} = 0.182$ .

In the “background blanching” condition, compression has the effect not only to flatten the inclusion and decrease depth but also to reduce $μ_{a}$ of the background, particularly the region of background directly above the inclusion. Although reduction in depth and flattening of the inclusion would act to increase optical contrast, the reduction in $μ_{a}$ above the inclusion acts to decrease optical contrast. Figure 6(b) shows the combined influence of these factors on the simulation $μ_{a}$ results, showing a decrease in optical contrast with compression. In the “full blanching” condition, the same reduction of $μ_{a}$ above the inclusion occurs, and the inclusion itself also experiences a reduction in $μ_{a}$ . These factors result in a reduction of contrast with compression that is more pronounced than that of the “background blanching” condition.

3.3.

Compression-Induced Contrast Change for all Conditions

In the “no blanching” condition, the vast majority of simulations experienced an increase in optical contrast from compression, across all sizes, depths, and SRs. In the two tissue blanching conditions, some simulations exhibited a decrease in contrast induced by compression, as demonstrated by the examples in Fig. 6. For simulations experiencing a compression-induced decrease in contrast, some maintained positive, though reduced, contrast values mid-compression, whereas others exhibited negative contrast values mid compression, where $μ_{a}$ above the inclusion was lower than the background value. The percentage of simulations exhibiting contrast increase, contrast decrease that remained positive compression, and negative mid-compression contrast is reported in Table 3.

Table 3

Percentage of simulations falling into three categories based on the effect of compression on the contrast value.

SR	Contrast increased with compression (%)	Contrast reduced, remained positive (%)	Contrast became negative (%)
“No blanching” condition
1	95	5	0
4	100	0	0
7	100	0	0
10	95	5	0
“Background blanching” model
1	80	15	5
4	30	40	30
7	30	30	40
10	20	35	35
“Full blanching” model
1	65	30	5
4	0	40	60
7	0	30	70
10	0	25	75

Results are separated into sub-tables for each of the three blanching models. For each table row, n=20 (five inclusion depths and four inclusion sizes).

Optical contrast before and during compression are illustrated for the two extremes of stiffness ratio, $SR = 1$ and $SR = 10$ , in Fig. 7. For comparison with other work, these values of contrast have also been converted to the T/N ratio using the projection of the inclusion boundary onto the surface as the “tumor” ROI and image margins greater than 45 mm from inclusion epicenter as the “normal” ROI. For UC simulations, the T/N values range from $\sim 1$ to 1.2 (0% to 20% increase in $μ_{a}$ ) and after compression values range from $\sim 0.9$ to 1.3 (10% decrease to 30% increase) depending on the blanching model used.

Fig. 7

Measured contrast for uncompressed (UC) and compressed (C) cases, shown for the two extremes: $SR = 1$ and $SR = 10$ . Initial inclusion depth is indicated by color, and inclusion radius is indicated by marker size. Plot colors adapted from Ref. 54.

This reduction in contrast that occurs at high SRs is explained by the greatest hydrostatic stress (and thus greatest tHb reduction) occurring in the region directly above the inclusion, in which the soft background material is squeezed between the stiff inclusion and the compression plate. A greater reduction in tHb over the inclusion than in the surrounding areas serves to offset the contrast from the high tHb content of the inclusion itself. In the “full blanching” conditions, there can also be significant tHb reduction of the inclusion top surface as well, contributing to negative contrast in the compressed state.

To assess how the initial simulation parameters of SR, depth ( $d_{0}$ ), and radius ( $r_{0}$ ) affect the contrast change with compression, Pearson correlations were calculated for three variables and the difference between final and initial contrast for the three blanching models. These values are reported in Table 4. In the “no blanching” condition, the initial radius exerts the most effect on contrast change, with large inclusions gaining the most contrast from compression, followed by initial inclusion depth, with shallow inclusions gaining the most contrast from compression. In the “background blanching” model, SR and depth are negatively correlated with contrast change, that is, shallow inclusions with low stiffness are associated with a gain of contrast, whereas stiffer and deeper inclusions are associated with a reduction in contrast. In these cases, though inclusions are brought closer to the surface and blanching does not occur in the inclusion itself, the effect is outweighed by the blanching of the background material above the inclusion. Due to the high sensitivity of SFDI to superficial tissue, this blanching above the inclusion results in a decrease in contrast with compression for stiff and deep inclusions. When considering size in the “background blanching” condition, a large inclusion radius results in a wider range of outcomes in contrast change, whereas small radius inclusions experienced only small contrast changes, be they positive or negative. Finally, in the “full blanching” model, stiffness, as well as the inclusion radius, is again negatively correlated with change in contrast, with larger inclusions experiencing the most contrast reduction from compression. However, it was noted that, as for the “full blanching” condition, the correlation of radius with contrast change is reversed when only $SR = 1$ is considered, that is, for $SR = 1$ , larger inclusions experienced a more positive change in contrast, whereas for all other SR, large inclusions experienced a more negative change in contrast. It can be said for both “background blanching” and “full blanching” models that a larger radius is associated with a higher magnitude of contrast change, the sign of which is determined by SR and, to a lesser extent, depth.

Table 4

Pearson correlation of initial simulation properties of SR, inclusion depth, and inclusion radius on contrast change from compression.

	“No blanching” condition		“Background blanching” model		“Full blanching” model
	Correlation	p-value	Correlation	p-value	Correlation	p-value
SR	0.11	0.35	−0.61	<0.001	−0.68	<0.001
Depth ( $d_{0}$ )	−0.53	<0.001	−0.50	<0.001	−0.01	0.91
Radius ( $r_{0}$ )	0.74	<0.001	−0.12	0.30	−0.42	<0.001

Correlations with p-value<0.05 are indicated in bold.

When simulations are repeated with tHb concentration in the tumor voxels increased from 50.6 to $63.3 μ M$ tHb for the “full blanching” condition, the trends in contrast change are found to be the same as those in Fig. 7, with SR exerting the greatest effect on contrast change. When the $63.3 μ M$ simulations are compared directly with the $50.6 μ M$ simulations, the increase in tumor tHb appears to have little effect on the trend in contrast change from compression in the “full blanching” model (not shown).

4. Discussion

In the “no blanching” condition, we examine the effect of inclusion size and depth on optical contrast. We confirm that, with optical properties of the inclusion and background held consistent, contrast is determined primarily on the final depth and shape of the inclusion top edge. This indicates that, for this measure of optical contrast, a compressed ellipsoidal inclusion is effectively indistinguishable from a spherical inclusion of radius equal to the ellipsoid top surface radius of curvature. In characterizing change in contrast based on initial parameters, the most gain in contrast was found for large inclusions and for initially shallow inclusions. The high contrast gain for large inclusions is explained by the idea that, for a given depth reduction, a larger mass of inclusion voxels is brought closer to the surface, resulting in a larger gain. Though deep inclusions experienced a greater absolute reduction in depth than initially shallow inclusions, these depth changes amounted to smaller contrast gains than shallow inclusions due to the exponential relationship between depth and contrast. When tHb content within the tumor is increased, the expected responses of greater baseline contrast and greater increase of contrast are observed. The relationship with stiffness is more complex, as stiff inclusions experienced a greater depth reduction than softer inclusions but had less tendency for the top surface to deform and become flatter. These two tendencies worked to counteract each other, and the SR exerted very little effect on contrast change.

Though these results represent a greatly simplified model, we may still consider these findings in terms of expected trends in a tumor during response to NAC, and we can predict that a tumor that becomes smaller would exhibit lower initial contrast and lower contrast gain from compression than prior to therapy, and a reduction in tumor tHb would exhibit the same effect. A reduction in tumor stiffness is not expected to have a large effect on initial contrast or contrast gain, except in cases in which the tumor is still quite large and quite deep, which might result in a smaller contrast gain. Thus, we can hypothesize that lower initial contrast and reduced contrast gain from compression could indicate a combination of tumor shrinkage and tHb reduction without being inconsistent with loss of tumor stiffness.

When a change in tHb in response to compression is modeled in the “background blanching” and “full blanching” conditions, the final geometry (depth and radius of curvature) is not sufficient for the predication of compressed contrast due to the strong effect of SR on $Δ tHb$ . The relationship between initial inclusion properties and contrast change with compression is also more complex, with the SR exhibiting the greatest effect in both blanching models. The soft material superficial to stiff inclusions experienced significant hydrostatic stress, and thus a reduction of tHb, in this volume region. Due to the high sensitivity to superficial tissue of SFDI, this region of reduced tHb resulted in lower contrast between inclusion and background regions than seen before compression. Although this reduction in contrast from compression is a potential concern for clinical use of this imaging method, it is also possible that the information can be gleaned from this effect to its dependence on tumor size and stiffness.

Although $SR = 1$ experiences a small positive contrast gain from compression, all other SRs experience a reduction in contrast, with a higher magnitude associated with higher stiffness. This reduction is also more pronounced in the “full blanching” condition as the surface of the inclusion also undergoes a reduction in tHb. After the SR, the initial inclusion depth shows the greatest correlation with contrast change in the “background blanching” model, with deep inclusions experiencing the most contrast reduction. In the “full blanching” condition, large inclusions experience more contrast reduction for SR values above one. More sophisticated mechanical modeling as well as in vivo imaging of breast tumors under compression using SFDI will be necessary to validate these trends. An increase in baseline tHb concentration within the tumor was not observed to have a significant effect on contrast gain in the “full blanching” condition.

Limitations of this study include the very low compression pressure used. The reaction force to the applied displacement was $\sim 0.5 N$ , compared with an applied force of $\sim 17 N$ by Sajjadi et al.⁷ That study involved whole-breast compression with a larger contact area than that used here, but any reasonable estimate of contact area yields a significantly higher pressure than that used here. Stronger compression likely would have induced a larger change in lesion depth, as well as different behavior in tissue blanching. For stiff inclusions ( $SR > 1$ ), the highest stresses were always experienced directly above the inclusion, and thus, the highest reduction in tHb happened in this location. The surrounding background material experienced lower stresses and thus a lower reduction in tHb. Negative contrast in many cases arose from the fact that blanching of the background surface occurred primarily above the inclusion with only small changes in tHb occurring further from the inclusion. This localized blanching of the surface served to offset the high tHb concentration of the inclusion, resulting in near-zero or negative final contrast. As we modeled tissue blanching to plateau at 20 mmHg, we would expect greater compression pressure to increase the area of the background that is experiencing the maximum tHb reduction. Higher pressure would result in maximal blanching of much of the background material surface, rather than only the region on top of the inclusion. With optical properties of the background surface closer to uniform in $x$ and $y$ , the high tHb content of the inclusion may still produce a positive contrast under these conditions.

At least one study reports a decrease in breast tissue $μ_{s}^{'}$ resulting from compression,⁵⁵ but the effect of tissue compression on scattering properties is generally not well characterized in the literature. This work assumed that $μ_{s}^{'}$ remains constant during compression, and further research is needed to account for possible changes in scattering during compression.

Although this simplified model accounts for spatially varying $Δ tHb$ induced by compression, the effect on time-varying hemodynamic changes is not covered by this work, nor are any changes in tissue oxygen saturation in response to compression. The effect of applied pressure on breast hemodynamics has been modeled in a simulation study by Nioka et al.,⁵⁶ and further research in this area would be very useful. Utilizing diffuse optical tomography, differential changes in tissue oxygen saturation, as well as the time-course of changes in tHb, have been demonstrated between tumor and healthy breast in response to compression.⁷^,²³ A greater decrease of oxygen saturation was observed in the tumor, and as the pressure decreased due to viscoelastic relaxation (with displacement held constant), a persistent reduced tHb level was observed in the tumor in contrast to a gradual recovery of tHb in the normal tissue. This differential response was shown to normalize in patients responding to NAC.⁷

In addition, the viscoelastic character of breast tissue was not modeled. The tissue was modeled as a solid material with hydrostatic stress being used to determine tissue blanching; however, tissue consists of a solid porous structure containing interstitial fluid, and future work should directly consider the effect of applied pressure of interstitial fluid pressure.⁵⁷

Compression-induced hemodynamics represent a complex interplay between blood inflow and outflow and local metabolic oxygen demand. Differential hemodynamic responses likely arise from mechanical differences;⁵^,⁴¹^,⁵⁸ increased metabolic activity within the tumor;⁴^,⁵⁹ abnormalities in the vascular network, including disorganized and tortuous blood vessels;⁶⁰^,⁶¹ and increased interstitial fluid pressure within the tumor.⁶²^,⁶³ Though not directly comparable with compression-induced changes, hemodynamic changes in response to breath-holding have also been shown to differ between tumor and healthy breast,⁶⁴ and normalization of these differences can be used for NAC response prediction,⁶ further highlighting the value of tracking time-varying changes in hemoglobin concentrations for breast cancer therapy monitoring. A more complex dynamic model would be needed to predict how these hemodynamics will be evident in SFDI.

5. Conclusion

In this work, combined mechanical and optical simulations were used to show how attributes of a breast tumor-mimicking inclusion influenced measured contrast with SFDI and tissue compression. The two tissue blanching models used in this study predicted that compression of a stiff tumor would result in a reduction of contrast between tumor and normal tissue due to squeezing the healthy tissue superficial to the tumor, an effect of higher magnitude for large and very stiff tumors. By contrast, the “no blanching” model predicted that contrast will be gained from compression across all inclusion parameters, though this model is of limited applicability to in vivo imaging. Though no definitive conclusion may be reached at this time and much additional research is needed, this work represents a first step toward the potential to use SFDI with local tissue compression for tracking changes in breast tumors such as reduction in size and reduction in stiffness.

Disclosures

The authors have no relevant conflicts of interest to disclose.

Code and Data Availability

Data and code used to generate results and figures are available on Code Ocean, https://doi.org/10.24433/CO.4469356.v1.

Acknowledgments

This work was supported by the National Institutes of Health (Grant No. T32-EB003392) – BiRM (Biomechanics in Regenerative Medicine) Training Program.

References

1.

M. Kaufmann et al., “Recommendations from an international consensus conference on the current status and future of neoadjuvant systemic therapy in primary breast cancer,” Ann. Surg. Oncol., 19 (5), 1508 –1516 https://doi.org/10.1245/s10434-011-2108-2 (2012). Google Scholar

2.

K. H. Koga et al., “Monitoring the response to neoadjuvant chemotherapy in breast cancer,” Neoadjuvant Chemotherapy - Increasing Relevance in Cancer Management, 25 –48 Intechopen( (2013). Google Scholar

3.

R. Mclaughlin and N. Hylton, “MRI in breast cancer therapy monitoring,” NMR Biomed., 24 (6), 712 –720 https://doi.org/10.1002/nbm.1739 (2011). Google Scholar

4.

R. Choe and T. Durduran, “Diffuse optical monitoring of the neoadjuvant breast cancer therapy,” IEEE J. Sel. Top. Quantum Electron., 18 (4), 1367 –1386 https://doi.org/10.1109/JSTQE.2011.2177963 IJSQEN 1077-260X (2012). Google Scholar

5.

D. Grosenick et al., “Review of optical breast imaging and spectroscopy,” J. Biomed. Opt., 21 (9), 091311 https://doi.org/10.1117/1.JBO.21.9.091311 JBOPFO 1083-3668 (2016). Google Scholar

6.

J. E. Gunther et al., “Dynamic diffuse optical tomography for monitoring neoadjuvant chemotherapy in patients with breast cancer,” Radiology, 287 (3), 778 –786 https://doi.org/10.1148/radiol.2018161041 RADLAX 0033-8419 (2018). Google Scholar

7.

A. Y. Sajjadi et al., “Normalization of compression-induced hemodynamics in patients responding to neoadjuvant chemotherapy monitored by dynamic tomographic optical breast imaging (DTOBI),” Biomed. Opt. Express, 8 (2), 555 –569 https://doi.org/10.1364/BOE.8.000555 BOEICL 2156-7085 (2017). Google Scholar

8.

A. Tank et al., “Diffuse optical spectroscopic imaging reveals distinct early breast tumor hemodynamic responses to metronomic and maximum tolerated dose regimens,” Breast Cancer Res., 22 (1), 1 –10 https://doi.org/10.1186/s13058-020-01262-1 BCTRD6 (2020). Google Scholar

9.

B. J. Tromberg et al., “Predicting responses to neoadjuvant chemotherapy in breast cancer: ACRIN 6691 trial of diffuse optical spectroscopic imaging,” Cancer Res., 76 (20), 5933 –5944 https://doi.org/10.1158/0008-5472.CAN-16-0346 CNREA8 0008-5472 (2016). Google Scholar

10.

E. A. Lim et al., “Diffuse optical tomography changes correlate with residual cancer burden after neoadjuvant chemotherapy in breast cancer patients,” Breast Cancer Res. Treat., 162 (3), 533 –540 https://doi.org/10.1007/s10549-017-4150-7 BCTRD6 (2017). Google Scholar

11.

D. Roblyer et al., “Optical imaging of breast cancer oxyhemoglobin flare correlates with neoadjuvant chemotherapy response one day after starting treatment,” Proc. Natl. Acad. Sci. U. S. A., 108 (35), 14626 –14631 https://doi.org/10.1073/pnas.1013103108 (2011). Google Scholar

12.

P. G. Anderson et al., “Optical mammography in patients with breast cancer undergoing neoadjuvant chemotherapy: individual clinical response index,” Acad. Radiol., 24 (10), 1240 –1255 https://doi.org/10.1016/j.acra.2017.03.020 (2017). Google Scholar

13.

Q. Zhu et al., “Breast cancer: assessing response to neoadjuvant chemotherapy by using US-guided near-infrared tomography,” Radiology, 266 (2), 433 –442 https://doi.org/10.1148/radiol.12112415 RADLAX 0033-8419 (2013). Google Scholar

14.

D. R. Busch et al., “Optical malignancy parameters for monitoring progression of breast cancer neoadjuvant chemotherapy,” Biomed. Opt. Express, 4 (1), 105 –121 https://doi.org/10.1364/BOE.4.000105 BOEICL 2156-7085 (2013). Google Scholar

15.

J. M. Cochran et al., “Tissue oxygen saturation predicts response to breast cancer neoadjuvant chemotherapy within 10 days of treatment,” J. Biomed. Opt., 24 (2), 021202 https://doi.org/10.1117/1.JBO.24.2.021202 JBOPFO 1083-3668 (2018). Google Scholar

16.

W. Zhi et al., “Predicting treatment response of breast cancer to neoadjuvant chemotherapy using ultrasound-guided diffuse optical tomography,” Transl. Oncol., 11 (1), 56 –64 https://doi.org/10.1016/j.tranon.2017.10.011 (2018). Google Scholar

17.

Q. Zhu et al., “Noninvasive monitoring of breast cancer during neoadjuvant chemotherapy using optical tomography with ultrasound localization,” Neoplasia, 10 (10), 1028 –1040 https://doi.org/10.1593/neo.08602 (2008). Google Scholar

18.

M. G. Pakalniskis et al., “Tumor angiogenesis change estimated by using diffuse optical spectroscopic tomography: demonstrated correlation in women undergoing neoadjuvant chemotherapy for invasive breast cancer?,” Radiology, 259 (2), 365 –374 https://doi.org/10.1148/radiol.11100699 RADLAX 0033-8419 (2011). Google Scholar

19.

A. E. Cerussi et al., “Diffuse optical spectroscopic imaging correlates with final pathological response in breast cancer neoadjuvant chemotherapy,” Philos. Trans. Math. Phys. Eng. Sci., 369 (1955), 4512 –4530 https://doi.org/10.1098/rsta.2011.0279 (2011). Google Scholar

20.

S. Jiang et al., “Measurement of pressure-displacement kinetics of hemoglobin in normal breast tissue with near-infrared spectral imaging,” Appl. Opt., 48 (10), D130 –D136 https://doi.org/10.1364/AO.48.00D130 APOPAI 0003-6935 (2009). Google Scholar

21.

D. R. Busch et al., “Blood flow reduction in breast tissue due to mammographic compression,” Acad. Radiol., 21 (2), 151 –161 https://doi.org/10.1016/j.acra.2013.10.009 (2014). Google Scholar

22.

B. Wang et al., “Dynamic schema for near infrared detection of pressure-induced changes in solid tumors,” Appl. Opt., 47 (16), 3053 –3063 https://doi.org/10.1364/AO.47.003053 APOPAI 0003-6935 (2008). Google Scholar

23.

S. A. Carp et al., “Hemodynamic signature of breast cancer under fractional mammographic compression using a dynamic diffuse optical tomography system,” Biomed. Opt. Express, 4 (12), 2911 –2924 https://doi.org/10.1364/BOE.4.002911 BOEICL 2156-7085 (2013). Google Scholar

24.

R. Al Abdi et al., “Mechanical and hemodynamic responses of breast tissue during compression,” in Biophotonics Congr.: Biomed. Opt. 2020 (Transl., Microsc., OCT, OTS, BRAIN), JTu3A.4 (2020). https://doi.org/10.1364/TRANSLATIONAL.2020.JTu3A.4 Google Scholar

25.

D. J. Cuccia et al., “Quantitation and mapping of tissue optical properties using modulated imaging,” J. Biomed. Opt., 14 (2), 024012 https://doi.org/10.1117/1.3088140 JBOPFO 1083-3668 (2009). Google Scholar

26.

S. Gioux, A. Mazhar and D. J. Cuccia, “Spatial frequency domain imaging in 2019: principles, applications, and perspectives,” J. Biomed. Opt., 24 (7), 071613 https://doi.org/10.1117/1.JBO.24.7.071613 JBOPFO 1083-3668 (2019). Google Scholar

27.

C. K. Hayakawa et al., “Optical sampling depth in the spatial frequency domain,” J. Biomed. Opt., 24 (7), 071603 https://doi.org/10.1117/1.JBO.24.7.071603 JBOPFO 1083-3668 (2018). Google Scholar

28.

C. M. Robbins et al., “Feasibility of spatial frequency-domain imaging for monitoring palpable breast lesions,” J. Biomed. Opt., 22 (12), 121605 https://doi.org/10.1117/1.JBO.22.12.121605 JBOPFO 1083-3668 (2017). Google Scholar

29.

M. Choi et al., “Assessment of pathologic response and long-term outcome in locally advanced breast cancers after neoadjuvant chemotherapy: comparison of pathologic classification systems,” Breast Cancer Res. Treat., 160 (3), 475 –489 https://doi.org/10.1007/s10549-016-4008-4 BCTRD6 (2016). Google Scholar

30.

W. F. Symmans et al., “Measurement of residual breast cancer burden to predict survival after neoadjuvant chemotherapy,” J. Clin. Oncol., 25 (28), 4414 –4422 https://doi.org/10.1200/JCO.2007.10.6823 JCONDN 0732-183X (2007). Google Scholar

31.

H. Jing et al., “Early evaluation of relative changes in tumor stiffness by shear wave elastography predicts the response to neoadjuvant chemotherapy in patients with breast cancer,” J. Ultrasound Med., 35 (8), 1619 –1627 https://doi.org/10.7863/ultra.15.08052 JUMEDA 0278-4297 (2016). Google Scholar

32.

O. Falou et al., “Evaluation of neoadjuvant chemotherapy response in women with locally advanced breast cancer using ultrasound elastography,” Transl. Oncol., 6 (1), 17 –24 https://doi.org/10.1593/tlo.12412 (2013). Google Scholar

33.

A. R. Gardner and V. Venugopalan, “Accurate and efficient Monte Carlo solutions to the radiative transport equation in the spatial frequency domain,” Opt. Lett., 36 (12), 2269 –2271 https://doi.org/10.1364/OL.36.002269 OPLEDP 0146-9592 (2011). Google Scholar

34.

C. M. Robbins et al., “Hand-held multi-wavelength spatial frequency domain imaging for breast cancer imaging,” Proc. SPIE, 10874 108740U https://doi.org/10.1117/12.2510399 PSISDG 0277-786X (2019). Google Scholar

35.

S. Tabassum et al., “Feasibility of spatial frequency domain imaging (SFDI) for optically characterizing a preclinical oncology model,” Biomed. Opt. Express, 7 (10), 4154 –4170 https://doi.org/10.1364/BOE.7.004154 BOEICL 2156-7085 (2016). Google Scholar

36.

G. Lorenzo et al., “Tissue-scale, personalized modeling and simulation of prostate cancer growth,” Proc. Natl. Acad. Sci. U. S. A., 113 (48), E7663 –E7671 https://doi.org/10.1073/pnas.1615791113 (2016). Google Scholar

37.

Y. Zhang et al., “Dynamic lung modeling and tumor tracking using deformable image registration and geometric smoothing,” Mol. Cell. Biomech., 9 (3), 213 –226 https://doi.org/10.1201/b12753-46 (2012). Google Scholar

38.

J. T. Oden et al., “Dynamic data-driven finite element models for laser treatment of cancer,” Numer. Methods Partial Differ. Equ., 23 (4), 904 –922 https://doi.org/10.1002/num.20251 NMPDEB 1098-2426 (2007). Google Scholar

39.

M. N. Rylander et al., “Optimizing heat shock protein expression induced by prostate cancer laser therapy through predictive computational models,” J. Biomed. Opt., 11 (4), 041113 https://doi.org/10.1117/1.2241310 JBOPFO 1083-3668 (2006). Google Scholar

40.

Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express, 17 (22), 020178 https://doi.org/10.1364/OE.17.020178 OPEXFF 1094-4087 (2009). Google Scholar

41.

A. Samani, J. Zubovits and D. Plewes, “Elastic moduli of normal and pathological human breast tissues: an inversion-technique-based investigation of 169 samples,” Phys. Med. Biol., 52 (6), 1565 –1576 https://doi.org/10.1088/0031-9155/52/6/002 PHMBA7 0031-9155 (2007). Google Scholar

42.

, “Friction and friction coefficients,” https://www.engineeringtoolbox.com/friction-coefficients-d_778.html (2004). Google Scholar

43.

B. Szabó and I. Babuska, “Boundary value problems,” Finite Element Analysis: Method, Verification and Validation, 71 –110 2nd ed.John Wiley & Sons, Inc. ( (2021). Google Scholar

44.

N. Chernov, “Circle fit (Pratt method),” https://www.mathworks.com/matlabcentral/fileexchange/22643-circle-fit-pratt-method (2021). Google Scholar

45.

S. A. Prahl, “Tabulated molar extinction coefficient for hemoglobin in water,” https://omlc.org/spectra/hemoglobin/summary.html (1998). Google Scholar

46.

T. Durduran et al., “Diffuse optics for tissue monitoring and tomography,” Rep. Prog. Phys., 73 (7), 076701 https://doi.org/10.1088/0034-4885/73/7/076701 RPPHAG 0034-4885 (2010). Google Scholar

47.

L. V. Wang and H.-I. Wu, Biomedical Optics, John Wiley & Sons, Inc., Hoboken, New Jersey, USA (2009). Google Scholar

48.

J. G. Betts et al., “Anatomy and phsiology,” (2013). Google Scholar

49.

D. T. Delpy et al., “Estimation of optical pathlength through tissue from direct time of flight measurement,” Phys. Med. Biol., 33 (12), 1433 –1442 https://doi.org/10.1088/0031-9155/33/12/008 PHMBA7 0031-9155 (1988). Google Scholar

50.

C. K. Hayakawa, J. Spanier and V. Venugopalan, “Comparative analysis of discrete and continuous absorption weighting estimators used in Monte Carlo simulations of radiative transport in turbid media,” J. Opt. Soc. Am. A, 31 (2), 301 https://doi.org/10.1364/JOSAA.31.000301 JOAOD6 0740-3232 (2014). Google Scholar

51.

R. Cubeddu et al., “Effects of the menstrual cycle on the red and near-infrared optical properties of the human breast,” Photochem. Photobiol., 72 (3), 383 –391 https://doi.org/10.1562/0031-8655(2000)0720383EOTMCO2.0.CO2 PHCBAP 0031-8655 (2000). Google Scholar

52.

S. L. Jacques, “Optical properties of biological tissues: a review,” Phys. Med. Biol., 58 (11), R37 –R61 https://doi.org/10.1088/0031-9155/58/11/R37 PHMBA7 0031-9155 (2013). Google Scholar

53.

S. Gioux et al., “Special section guest editorial: special section on spatial frequency domain imaging,” J. Biomed. Opt., 24 (7), 071601 https://doi.org/10.1117/1.JBO.24.7.071601 JBOPFO 1083-3668 (2019). Google Scholar

54.

, “cbrewer: colorbrewer schemes for Matlab,” https://www.mathworks.com/matlabcentral/fileexchange/34087-cbrewer-colorbrewer-schemes-for-matlab (2021). Google Scholar

55.

S. A. Carp et al., “Compression-induced changes in the physiological state of the breast as observed through frequency domain photon migration measurements,” J. Biomed. Opt., 11 (6), 064016 https://doi.org/10.1117/1.2397572 JBOPFO 1083-3668 (2006). Google Scholar

56.

S. Nioka et al., “Simulation study of breast tissue hemodynamics during pressure perturbation,” Adv. Exp. Med. Biol., 566 17 –22 https://doi.org/10.1007/0-387-26206-7_3 AEMBAP 0065-2598 (2005). Google Scholar

57.

A. L. Darling et al., “Interstitial fluid pressure in soft tissue as a result of an externally applied contact pressure,” Phys. Med. Biol., 52 (14), 4121 –4136 https://doi.org/10.1088/0031-9155/52/14/007 PHMBA7 0031-9155 (2007). Google Scholar

58.

T. A. Krouskop et al., “Elastic moduli of breast and prostate tissue under compression,” Ultrason. Imaging, 20 (4), 260 –274 https://doi.org/10.1177/016173469802000403 ULIMD4 0161-7346 (1998). Google Scholar

59.

D. R. Busch et al., “Toward noninvasive characterization of breast cancer and cancer metabolism with diffuse optics,” PET Clin., 8 (3), 345 –365 https://doi.org/10.1016/j.cpet.2013.04.004 (2013). Google Scholar

60.

J. A. Nagy et al., “Why are tumour blood vessels abnormal and why is it important to know?,” Br. J. Cancer, 100 (6), 865 –869 https://doi.org/10.1038/sj.bjc.6604929 BJCAAI 0007-0920 (2009). Google Scholar

61.

M. D. Abeloff et al., “Vascular and interstitial biology of tumors,” Abeloff’s Clinical Oncology, 105 –124 4th ed.Churchill Livingstone( (2008). Google Scholar

62.

R. K. Jain, Y. Boucher and N. Wolmark, “Interstitial hypertension in human breast and colorectal tumors,” Cancer Res., 52 (22), 6371 –6374 https://doi.org/10.1007/s12195-020-00661-w CNREA8 0008-5472 (1992). Google Scholar

63.

S. D. Nathanson and L. Nelson, “Interstitial fluid pressure in breast cancer, benign breast conditions, and breast parenchyma,” Ann. Surg. Oncol., 1 (4), 333 –338 https://doi.org/10.1007/BF03187139 (1994). Google Scholar

64.

M. L. Flexman et al., “Optical biomarkers for breast cancer derived from dynamic diffuse optical tomography,” J. Biomed. Opt., 18 (9), 096012 https://doi.org/10.1117/1.JBO.18.9.096012 JBOPFO 1083-3668 (2013). Google Scholar

Biography

Constance M. Robbins is a postdoctoral associate in the radiology department of the University of Pittsburgh. He earned his PhD in biomedical engineering at Carnegie Mellon University, advised by Dr. Jana Kainerstorfer. He received his BS in bioengineering from the University of Pittsburgh in 2015. Their doctoral research focused on non-invasive optical imaging of breast hemodynamics, and their current work involves combined electrophysiological recording and wide-field calcium imaging in mouse brain cortex.

Kuanren Qian is a PhD student in the Computational Bio-Modeling Lab at the Department of Mechanical Engineering of Carnegie Mellon University. He received his BS in mechanical engineering from the University of California, San Diego, in 2018 and his MSR in mechanical engineering from Carnegie Mellon University in 2020. His current work includes isogeometric analysis, neuron growth simulation, phase field method, truncated hierarchical B-splines, machine learning, and physics-informed modeling.

Yongjie Jessica Zhang is the George Tallman Ladd and Florence Barrett Ladd Professor of Mechanical Engineering at Carnegie Mellon University with a courtesy appointment in biomedical engineering. She received her BEng in automotive engineering and MEng in engineering mechanics from Tsinghua University, China, and her MEng in aerospace engineering and engineering mechanics and PhD in computational engineering and sciences from the Institute for Computational Engineering and Sciences (now Oden Institute), the University of Texas at Austin. Her research interests include computational geometry, isogeometric analysis, finite element method, data-driven simulation, image processing, and their applications in computational biomedicine and engineering. She is a fellow of ASME, SIAM, IAMBE, AIMBE, USACM, SMA, and ELATES at Drexel. She is the editor-in-chief of Engineering with Computers.

Jana M. Kainerstorfer is an associate professor of biomedical engineering at Carnegie Mellon University. She earned her PhD in physics from the University of Vienna, Austria, in partnership with the National Institutes of Health and worked as a postdoctoral fellow at Tufts University. Her lab at Carnegie Mellon University focuses on optical imaging technology development for monitoring cerebral and tissue health. She serves on program committees at national and international conferences, is an associate editor for Journal of Biomedical Optics, and is a senior member of Optica.

CC BY: © The Authors. Published by SPIE under a Creative Commons Attribution 4.0 International License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.

Citation Download Citation

Constance M. Robbins, Kuanren Qian, Yongjie Jessica Zhang, and Jana M. Kainerstorfer "Monte Carlo simulation of spatial frequency domain imaging for breast tumors during compression," Journal of Biomedical Optics 29(9), 096001 (14 September 2024). https://doi.org/10.1117/1.JBO.29.9.096001

Received: 17 April 2024; Accepted: 19 August 2024; Published: 14 September 2024

JOURNAL ARTICLE
19 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY

GET CITATION

KEYWORDS

Monte Carlo methods

Tumors

Tissues

Breast

Voxels

Optical properties

Biomedical optics

Significance

Aim

Approach

Results

Conclusions

1.

Introduction

2.

Methods

Fig. 1

2.1.

Mechanical Simulation

2.2.

Processing of FEA Output

2.3.

Implementation of Monte Carlo Simulations

Eq. (1)

Eq. (2)

2.4.

Processing of MCX Output

3.

Results

3.1.

Results without Change in tHb (“No Blanching” Condition)

Fig. 2

Fig. 3

Table 1

Table 2

Fig. 4

Fig. 5

Eq. (3)

3.2.

Results with Change in tHb (“Background Blanching” and “Full Blanching” Conditions)

Fig. 6

3.3.

Compression-Induced Contrast Change for all Conditions

Table 3

Fig. 7

Table 4

4.

Discussion

5.

Conclusion

Disclosures

Code and Data Availability

Acknowledgments

References

Biography

Show All Keywords

Keywords/Phrases

Search In:

Publication Years