Dynamic functional and mechanical response of breast tissue to compression

Physiological tissue dynamics following breast compression offer new contrast mechanisms for evaluating breast health and disease with near infrared spectroscopy. We monitored the total hemoglobin concentration and hemoglobin oxygen saturation in 28 healthy female volunteers subject to repeated fractional mammographic compression. The compression induces a reduction in blood flow, in turn causing a reduction in hemoglobin oxygen saturation. At the same time, a two phase tissue viscoelastic relaxation results in a reduction and redistribution of pressure within the tissue and correspondingly modulates the tissue total hemoglobin concentration and oxygen saturation. We observed a strong correlation between the relaxing pressure and changes in the total hemoglobin concentration bearing evidence of the involvement of different vascular compartments. Consequently, we have developed a model that enables us to disentangle these effects and obtain robust estimates of the tissue oxygen consumption and blood flow. We obtain estimates of 1.9±1.3 µmol/100mL/min for OC and 2.8±1.7 mL/100mL/min for blood flow, consistent with other published values.


Introduction
Near infrared spectroscopy and imaging of breast cancer has seen rapid progress over the past decade as a non-invasive, low-cost method to probe the tissue physiological state through the quantification of tissue chromophores, such as oxy and deoxy-hemoglobin, water and lipids. Various research groups have reported success in characterizing breast tumor properties for diagnosis [1][2][3][4][5][6][7][8][9][10][11][12][13][14], as well as treatment monitoring purposes [5,15,16]. A common characteristic of the vast majority of published optical studies is their focus on static, single-shot images. Dynamic optical imaging, made possible by technological advancement in recent years has allowed the exploration of a new dimension -time-resolved changes in tissue properties. Largely, these methods have been applied to functional brain imaging [17][18][19][20][21][22]. However, there is a growing trend to build optical spectroscopy systems that measure the dynamic features of breast tissue, either intrinsic [10,22] or in response to external mechanical stimulation (compression) [23][24][25].
Tissue dynamics following compression are governed by the interplay of tissue biomechanics and metabolic activity, and thus offer access to a rich information content related to these mechanical and metabolic characteristics. We have recently reported [24] preliminary data describing an opportunity to estimate tissue metabolic parameters (oxygen consumption and blood flow) by analyzing the time resolved changes in hemoglobin oxygen saturation in response to a step change in external compression force. The values of these metabolic parameters are known from other imaging modalities to be elevated in malignant breast lesions [26][27][28]. More fundamentally, tissue mechanical properties (which influence the viscoelastic response to compression) are known to vary in the presence of abnormalities, a problem widely investigated in the field of breast elastography [29][30][31][32]. Dynamic optical measurements following tissue compression permit the (indirect) characterization of both of these processes and thus may lead to the identification of novel optical disease markers to complement those derived from steady state measurements.
In this paper we present measurements of breast tissue viscoelastic relaxation after step compression, together with the associated hemodynamics. Here we extend our previous results [24] using a new set of measurements on 28 healthy subjects, by employing a longer duration compression to unveil the effects of viscoelastic relaxation on the observed hemodynamics, together with continuous monitoring of the force necessary to maintain the compression plates in place. We note evidence of the involvement of different vascular compartments during the relaxation process, which enables us to propose a model for disentangling the tissue relaxation effects from those due to tissue metabolism. Finally, we use this model to obtain robust estimates of oxygen consumption and blood flow from an objectively chosen subset of the measurements and compare these with our earlier results.

Experimental Setup
The mammography-like compression system has been previously described by our group [24]. Briefly, the system (shown in Fig. 1) consists of a set of parallel plastic compression plates with integrated optical fibers for light delivery and collection. A computer controlled translation stage mounted on a platform is used to vary the distance between the plates to achieve various levels of compression. The top plate is attached to the translation stage by means of two strain gauges that allow the monitoring of the compression force. The bottom plate is mounted on the same platform as the translation stage and it contains the optical fibers.
A frequency domain near infrared spectrometer was used to acquire the tissue measurements. The ISS Imagent system (model 96208, ISS, Champaign, IL) uses four photomultiplier tube detectors and eight time-multiplexed laser diodes emitting at 635, 670, 690, 752, 758, 782, 811 and 831 nm, respectively. The modulation frequency is set to 110 MHz and a full acquisition cycle over all wavelenghts is completed 12.5 times a second. The eight 400 µm illumination fibers are assembled into a 2 mm fiber bundle. Three 3 mm diameter fiber bundles are used to collect diffusely reflected light. The tips of the fibers are arranged along a line parallel to the edge of the lower compression plate and in direct contact with the breast tissue, about 1.5 cm from the edge of the plate. The detection fiber bundles are placed at 1.5, 2.0 and 2.6 cm, respectively from the illumination bundle.

Human subject protocol
Measurements were conducted under the Massachusetts General Hospital Institutional Review Board (IRB) protocol #2003-P-000135 on 28 healthy volunteers aged 23 to 56, with no open wounds on the breast, no breast implants, and no breast biopsy within the previous three months. The following procedure was followed: first, the volunteer's breast was placed in contact with the lower compression plate, with her chest wall resting against the plate's edge. The upper plate was then brought down until it just made contact with the breast. Baseline optical properties were recorded for 90 seconds, followed by application of 6 lbs of compressive force, after which the breast was kept under compression for another 90 seconds with continuous optical and force recording. Then the compression was released and the entire cycle was repeated a total of 3 times. We used 6 lbs (approximately one third of typical mammographic compression) of compression because it was comparatively comfortable for the subjects and our previous study [24] indicated that 6 lbs of force were sufficient to induce the tissue dynamics that we are interested in investigating.

Frequency domain measurements
Tissue optical properties (i.e, the absorption (μ a ) and reduced scattering coefficients) were derived from measurements at multiple source-detector separations using a multi-distance fitting algorithm as described in our previous publication [24] as well as in Refs. [33,34]. Briefly: (1) (2) where ω is the modulation frequency of the light sources and c t is the speed of light in tissue. I AC is the measured AC amplitude of the diffusely reflected light, Φ the corresponding phase shift, S Φ and S AC are the slopes characterizing the assumed linear relationship between the source detector separation r, and Φ and ln(r 2 I AC ), respectively. Oxy-hemoglobin ([HbO 2 ]) and deoxy-hemoglobin ([HbR]) concentrations were then derived from a spectral fit to the measured optical properties assuming a water fraction value of 30% [35][36][37].

Partial breast occlusion model
Following a similar approach to our previous publication [24] we write a mass balance for oxygenated hemoglobin (HbO 2 ) within the measurement volume (V t ): (3) where [HbO 2 ] is the molar concentration of oxygenated hemoglobin in tissue, OC is the tissue oxygen consumption (mol O 2 .s −1 ), F in and F out are the amounts of blood flow into and out of the measurement volume (in L.s −1 ), respectively, [HbT] b is the blood total hemoglobin concentration, and S a ,O 2 and S v ,O 2 are the arterial and venous oxygen saturation, respectively.
We note that unlike in our previous publication [24], we have not assumed F in = F out . Instead, we will take advantage of having a real-time measurement of the tissue total hemoglobin concentration ([HbT] t ) to obtain the relationship between inflow and outflow, given by the following equation: (4) Also, we note that the near-infrared measurement of tissue hemoglobin saturation (SO 2 ) represents a weighted average of the arterial and venous saturations: Hence: (5) Finally, since [HbO 2 ] = SO 2 [HbT] t , by chain rule: (6) Rewriting equation (3) gives: (7) Dividing by [HbT] t and grouping terms by SO 2 (t) gives: (8) Because we have [HbT] t (t) as a measurement, we can also compute , thus we will recast equation (8) as an ordinary differential equation with time-variable coefficients in the following form: (9) where (10) (11) Equation (9) is then solved numerically subject to the initial condition SO 2(t=0) = SO 2,init (the measured initial value of SO 2 ) using the ode45 function in Matlab (Mathworks Inc., Natick, MA). Low pass filtering of the [HbT] t trace (0.1 Hz) is necessary to obtain reasonable derivative values. When results for flow are presented later in the paper, they refer to F in . Metabolic parameters are derived volumetrically, i.e. , in units of µmol/100mL/min and mL/ 100mL/min, respectively. Note that oxygen consumption and blood flow are assumed constant during the measurements, due to the need for a sufficient time interval to obtain stable fitting results. While this is likely a very good assumption for oxygen consumption, blood flow is likely recovering during tissue relaxation, and the fit value represents a measure of the average blood flow during the measurement period.
In our data analysis we assumed S a ,O 2 =0.98 (normal physiological value), f =0.2, and [HbT] b =0.72 mM. The choices of f and [HbT] b deserve further explanation. While we were not able to find any published data on the arterial blood volume fraction in breast tissue, such measurements exist for brain tissue [38,39], leading us to choose 20% as a reasonable first approximation.With respect to the choice of [HbT] b , the normal hemoglobin concentration is in the range of 120-160 g/L for females [40]. However, the blood in the region probed by our spectrometer is dominantly contained in capillary vessels, whose hematocrit values is known to be 2-5 times lower than the systemic hematocrit (due to the Fahraeus effect) [41]. For the purposes of this study, we assumed a value of 140 g/L for systemic hemoglobin concentration and a 3 times dilution factor in the capillaries, resulting in a value of 0.72 mM for [HbT] b . We note that the choice of [HbT] b does influence the flow results, but not the calculated oxygen consumption, because the inflow and outflow terms appear as part of a product with [HbT] b in the mass balance equations. Figure 2 shows example recordings from three different subjects. The left column subplots display total hemoglobin ([HbT] t ), the middle column subplots display hemoglobin oxygen saturation (SO 2 ) and right side subplots display the amount of force necessary to hold the compression plates in place, respectively. The shading indicates the periods when the upper compression plate is stationary, after having been brought down until at least 6 lbs of force were applied to the tissue. Consistent, repeatable features are observed in these plots for all the parameters measured.

Compression induced dynamics
An initial fast decrease in [HbT] t is observed when compression is applied followed by a slow recovery after the compression plates stop moving. SO 2 also decreases somewhat as compression is applied and the change appears to correlate with the concomitant reduction in [HbT] t . This initial fast decrease is followed by a continued slower decrease that generally follows an approximately exponential decay profile. As described in our previous publication [24] we believe this is due to the tissue oxygen consumption exceeding oxygen delivery as compression induces a reduction in blood flow. An important feature can be noted at the end of the first compression cycle in Fig. 2 Subject 1, and to a lesser extent in Fig. 2 Subject 2 -the SO 2 trace curves up towards the end of the compression plateau. This appears to be a result of SO 2 modulation by the reduction in compression force leading to blood inflow, reflected in the increasing [HbT] t . A possible mechanism that explains this coupling is proposed in §4.2.
Several features are distinguished in the force dynamics. First, the force ramps up as the compression stage applies pressure to the tissue. Then, once the control system has determined that a force level of 6 lbs has been exceeded, the compression plate stops and the recorded force (and hence pressure) rapidly decreases over the first few seconds, followed by continued gradual reduction with an exponential-like profile during the rest of the compression plateau. The time of the transition to the exponential decay profile was automatically determined using a custom Matlab script (described in the Appendix), and the pressure drop during the two periods is analyzed separately. Table 2 summarizes the average pressure decrease from the entire volunteer set. The fast drop is largely constant in relative terms while the amount of slow relaxation decreases in later compression cycles. Table 3 presents the average temporal parameters of the pressure relaxation. While the fast pressure drop occurs over only a few seconds and is relatively constant in length, the slow drop exhibits progressively longer time constants during successive compression cycles.
A notable feature of the results shown in Fig. 2 is the relationship between blood volume and pressure during the relaxation process. Fig. 3 shows a plot of [HbT] t vs. pressure for three of our subjects. The rate of change of [HbT] t with respect to pressure is lower at higher pressure (low apparent vascular compliance) and higher at lower pressures (high apparent vascular compliance). The apparent increased compliance at lower pressures may be further amplified by reactive hyperemia, that can occur in occluded blood vessels and has a vasodilatory effect.
It is important to note some potential sources of systematic errors in the hemoglobin measurements presented. We assumed a fixed 30% water volume fraction and considered the contribution to absorption of other chromophores such as lipids and collagen to be negligible. This choice was made because our PMT based instrumentation is not well suited for quantification of these chromophores, which is best done using measurements beyond 850 nm. If water and lipid concentration deviate from our assumptions, errors of up to 15-20% and 8-10% may affect the estimated values of oxy and deoxy-hemoglobin, respectively. The largest errors would be encountered in subjects with high water volume fractions. Fortunately, the volume fractions of water and lipids generally counter-vary [47], thus limiting the estimation errors. At the same time, neglecting collagen absorption may result in over-estimation of deoxy-hemoglobin [48], with an associate underestimation of SO 2 . Finally, a further source of errors may arise from the shallower penetration depths of 635 and 670 nm light compared to the 690 to 830 nm range -while multi-wavelength fitting generally improves the accuracy of hemoglobin measurements, if significant tissue heterogeneity is present in the optically probed volume, a small compression dependent variation may affect the apparent hemoglobin concentration. Nevertheless, while these error sources make subject to subject comparisons rather difficult, they have only limited impact on the hemodynamic features described in this section, as they are observed on an individual subject basis.

Metabolic parameter estimation
As mentioned in the previous section, we noted an apparent correlation between SO 2 and [HbT] t , clearly evidenced for example in Fig. 2 Subject 1 where the tissue hemoglobin saturation begins to increase towards the end of the first compression cycle. Such behavior is not predicted by the partial breast occlusion model because of the constant flow assumption ( §2.4), hence pre-processing of the experimental data is necessary for obtaining reliable estimates of the tissue metabolic parameters.
The inter-dependence of SO 2 and [HbT] t is likely due to the influx of oxygenated arterial blood during tissue relaxation, as well as potential time-varying partial volume effects due to tissue mechanical rearrangement in the field of view of the optical measurement. As a first order approximation we believe it is reasonable to assume the variation in blood volume is confined to the arterial vascular compartment and we propose to "subtract" the influence of this arterial blood inflow on the tissue SO 2 . To this end, assuming the arterial saturation S a,O 2 =0.98, we calculate a "corrected" SO 2 in the following fashion: (13) where [HbT] t,begin is the total hemoglobin concentration at the beginning of the current compression period (just after the compression plate has stop moving down). Further justification for this correction is given in §4.2. The effect of this method is exemplified in Fig.  4, where the adjusted SO 2 is fit well by the partial breast occlusion model.
A number of criteria had to be met in order to obtain meaningful estimates of tissue metabolic parameters. Measurements were excluded from analysis if: 1. Evidence of fiber coupling or instrument saturation artifacts were present. The detector gains were set in the uncompressed state of the breast, and in a number of cases the reduction in absorption during compression lead to signal intensity changes that exceeded the dynamic range of the ISS instrument.

2.
Excessive noise was present in the SO 2 (t) measurement (peak to peak noise more than 30% of the total saturation change). Measured saturation changes were small (on the order of 1-4 %) and the non-linear relationship between SO 2 and hemoglobin concentrations led to enhancement of measurement uncertainty. Averaging multiple measurements together with a constant force measurement protocol are possible solutions.

3.
The second derivative of a second order polynomial fit to SO 2,corrected was negative. Such behavior is not predicted by the partial breast occlusion model and is likely caused by time-varying partial volume effects which are not accounted for in our model. True tomographic reconstruction of the dynamic tissue optical properties is likely necessary to overcome this phenomenon.

4.
A fast increase in [HbR] was noted to occur as compression was applied. This is again likely due to a partial volume effect as a result of tissue rearrangement (tissue with higher HbR (e.g. glandular tissue) comes into view). Again, dynamic tomographic imaging is needed to illuminate this effect.
Approximately 21% of the 44 successfully completed 3 cycle breast measurements met all four criteria (as detailed in Table 4) and were used to estimate tissue blood flow and oxygen consumption. Moving to tomographic imaging as well as optimizing the measurement protocol is hoped to significantly reduce the percentage of rejected measurements in future studies. To quantify the improvement brought by the inclusion of [HbT] variation and SO 2 correction in the metabolic model we also processed these data sets using the method in Ref. [24], which was based on the assumption of constant blood volume (d[HbT]/dt = 0) during compression plateaus. Table 5 displays the average volumetric blood flow and oxygen consumption estimated from the current volunteer data set using both the methods. The oxygen consumption (OC) and blood flow (BF) values from the new model are fairly consistent from cycle-to-cycle, unlike the ones estimated using the previously published model [24] (despite being applied to the same data). Also, the coefficient of determination, R 2 indicates a much closer fit to experimental data for the model in this publication and the corresponding metabolic estimates also exhibit lower relative standard deviations across subjects. The poor performance of the constant [HbT] model (Ref. [24]) can be traced to its failure to fit measurements with uptrending SO 2 such as the one exemplified in Figure 4. For numerical reasons, the failed cases resulted in high estimates for the metabolic parameters, hence the much higher average blood flow and oxygen consumption. In contrast, the enhanced partial breast occlusion model achieves a good fit to all the selected measurements. The simpler model in Ref. [24] can match this level of performance only on datasets pre-selected to have nearly flat [HbT] (less than 4% variation), as was done in Ref. [24]. Of the 28 individual (single cycle) compression recordings that meet all the criteria listed in Table 4, only 8 exhibit less than 4% [HbT] variation, highlighting the improved practicality of the method described in this paper.

Discussion
Monitoring tissue dynamics presents a great opportunity for enriching the information content of optical measurements if a near-infrared spectroscopy (NIRS) system with sufficient time resolution is used. The total hemoglobin concentration (a surrogate measure of blood volume) and hemoglobin oxygen saturation are the physiologically relevant parameters. During compression, blood volume dynamics (and hence [HbT] t ) are affected solely by the tissue mechanical properties (in particular, the tissue viscoelasticity), while the evolution of hemoglobin saturation is governed by the interplay of tissue biomechanics and tissue metabolism (in particular, oxygen consumption), as illustrated by our results.

Pressure induced hemodynamics
As shown in Fig. 2 the relative magnitude of the initial [HbT] t reduction increases while the amount of recovery decreases with repeated compression. A possible explanation is that during the first compression easily movable liquids are pushed out of the tissue (this includes blood through the circulatory system and interstitial fluid through the capillaries and lymphatic system) and the tissue collagen matrix begins to stretch. While the tissue recovers part of its volume when compression is released due to its elasticity, the breast thickness does not return to its initial uncompressed level, indicating the collagen matrix remains stretched. Subsequent compression encounters less of an elastic resistance and a larger relative amount of blood is removed and at the same time less recovery occurs during the compression plateau.
These characteristics of the [HbT] t variation correlate well with the pressure measurements. The fast pressure drop phase reflects the inertia of fluid flow, while the slow pressure relaxation occurs as the tissue collagen matrix stretches [30]. Constant load breast tissue creep measurements have shown these phenomena to occur over the scale of 10 seconds and minutes, respectively [30]. While the time constants measured for the slow relaxation phase of our data agree well with the literature, the fast pressure drop appears shorter that described in the literature. However, this discrepancy can be easily explained if we note that our compression is applied progressively, during a period of 10-15 seconds, and thus fluid removal from tissue substantially occurs before the compression process is completed. Hence, after the compression plate has stop moving we are only able to observe the tail end of the fluid flow.
Since [HbT] t is representative of the total blood volume, (in the absence of significant hematocrit changes) the data in Fig. 3 represents evidence that vascular compartments of varied compliance participate in the blood return process. At the initially higher pressure values compliance is low, likely representing volume changes in the arterial compartment, followed by increased compliance as the pressure relaxes, likely representing volume changes in the capillary and venous compartments.A possible interpretation of these phenomena is the involvement of the tissue collagen matrix as a support structure during the compression/ relaxation process. The stress-strain relationship for breast tissues is non-linear, with the tissue stiffness increasing rapidly beyond a few percent strain [42,43]. As vessel size reduces due to compression the tissue matrix begins to support a larger percentage of the external pressure, in effect shielding the vasculature from further compression. The more compliant capillary and venous compartments experience a faster reduction in size, and begin to be shielded by tissue at lower pressures than the arterial compartment. Hence at the larger pressure values associated with full compression we expect volume changes to occur predominantly in the arterial compartment. Such an effect, coupled with the presence of venous valves that prevent backflow suggests that increases in blood volume due to pressure relaxation are likely occuring predominantly in arteries and arterioles. This forms the basis of our model for "correcting" the SO 2 measurements to disentangle metabolic changes from relaxation driven effects.
A downside of the point-measurement approach is the possibility that in moderately dense breasts a boundary between adipose and glandular tissue might occur near the field of view of our diffuse reflection measurement. In such a case, compression may lead the more stiff glandular tissue to "push away" the adipose parenchyma near the skin resulting in a net increase in the glandular tissue volume fraction in the optically probed volume. Such a process may explain why the enhanced metabolic model still cannot fit the measurements in a number of cases (Inclusion Criteria 3 and 4 in §3.2): the coming into view of more metabolically active glandular tissue with decreased saturation and increased blood volume will cause unexpected behavior including an accelerating decrease in SO 2 and an increase in HbR during compression even as the total blood volume may be decreasing.

Metabolic parameter estimation
A first opportunity for validation of the oxygen consumption and blood flow estimates is the comparison with the findings of Ref. 24 A more stringent test consists of a comparison with MRI and PET studies using tracerbased methods [26,27,49,50]. These have found breast tissue blood flow to be in the 4-6 mL/100mL/ min range, while one study ( [26]) quotes an oxygen consumption rate of 19.8 µmol O 2 /100mL/ s. Our estimate of blood flow, 2.8 ml/min/100 ml is between 1.3 and 2 times lower than the MR and PET values. When taking into account the compression induced flow reduction we believe our estimate can be considered to be in fairly good agreement with literature values, especially given the inherent variation due to the different modeling approach and average blood vessel size sensitivity (our method is particularly sensitive to capillary vessels, while MR and PET tracer studies have a significant signal contribution from larger vessels). With respect to oxygen consumption, our estimate is much lower than that offered by Beaney et al. However, it is possible to argue that Beaney's estimate is on the high side, since given the quoted blood flow rate of 3.96 mL/100 mL/min, an oxygen extraction fraction of over 50% is implied, much higher than the typical 25% in a healthy person. Another reason for the difference may be that we measure a likely less metabolically active region of the breast (near the skin). To elucidate this difference we believe a dynamic diffuse optical tomography study is required, in order to measure regional differences in metabolic parameters.
While the values of the optically derived metabolic parameters appear to only have moderate accuracy in absolute terms (stemming from propagation of the HbO/HbR estimation errors described in §3.1, as well as fairly simple modeling approach), the ability to obtain consistent estimates from a moderately size subject population gives us confidence to attempt to use the partial breast occlusion model to image breast cancer patients. While the malignancy threshold may vary from person to person, we expect a repeatable contrast in oxygen consumption and blood flow between healthy and diseased tissue.

Conclusion
Breast tissue compression exposes a rich variety of dynamic features that offer additional opportunities for identifying disease markers, especially those associated with malignancies. This could be achieved by assessing tissue mechanical properties through [HbT] t dynamics, as well as metabolic differences gleaned from the evolution of tissue oxygen saturation following compression. Both mechanical and metabolic transients occur over similar timescales, and hence must be considered together for good accuracy. Clearly, blood volume variations are driven by pressure variations, which in turn are governed by the tissue viscoelasticity. This paper shows that blood volume variations convey information about the viscoelastic properties of breast tissue, and that this information can be used to increase the accuracy of metabolic modeling. Future work will include developing a Windkessel model [44] to more accurately describe the response of different vascular compartments to pressure variations and its linkage to the observed optical properties, using the foundation offered by previous work on non-linear elastic [45] or poroelastic [46] representations of biological tissue. Clinically, tomographic optical imaging may permit within-patient region of interest comparisons of these new biomarkers to obtain diagnostic information, without the need for absolute quantification. Experimental setup. A computer controlled translation stage applies compression to a subjects' breasts while diffuse reflection measurements are acquired at several source detector separations. Strain gauges are integrated into the compression frame for force monitoring. Representative hemodynamic measurement. An initial decrease in blood volume and a small decrease in saturation are observed during the initial compression, followed by a slow increase in blood volume and a continued decrease in saturation. Data from three subjects is presented, with the total hemoglobin [HbT] trace on the left, followed by the hemoglobin saturation time course, and finally, the compression force needed to maintain the position of the compression plates on the right. Total hemoglobin concentration vs pressure plots from three subjects, showing lower apparent compliance at high pressure and higher apparent compliance at lower pressures. Data is taken from measurements performed on the same volunteers shown in Figure 2.