A method for estimating the oxygen consumption rate in multicellular tumour spheroids

1Medical imaging has advanced to the stage where it is now technologically feasible to incorporate positron emission tomography (PET)/CT information into radiotherapy planning. Strategies for using this information have been proposed that range from assisting definition of boost dose regions to the design of IMRT plans which spatially modulate dose voxel-by-voxel on the basis of the functional images, commonly referred to as dose painting. However, these strategies involve significant departures from the current practice of delivering uniform population-based dose to the effective CT-derived treatment volume, and thus away from the clinical evidence base that informs current treatments. While the incorporation of functional imaging holds much promise for radiotherapy planning, for example, by selectively boosting regions of hypoxia, there is still a large amount unknown about the biology underpinning it and how to relate this to the images obtained.

2Hypoxia is one of the major hallmarks of solid tumours; it occurs when the oxygen requirements of a tissue outstrip the ability of the local vasculature to supply oxygen. It has been known since the 1950s that hypoxic tissue is more radioresistant than well-oxygenated tissue, and this factor has a large impact on a treatment, requiring higher levels of radiation to elicit the same cell kill. In the context of radiotherapy, the oxygen effect is quantified by the oxygen enhancement ratio, which is the increased effect of ionizing radiation on cell kill in the presence of molecular oxygen. This is taken as unity for anoxic cells and has been shown to be around 2.5–3 in well-oxygenated cells in vitro. In vivo tumours typically have heterogeneous oxygen distributions throughout. Imaging oxygen distribution in tumours has consequently been the focus of much interest. PET with tracers, for example 18F-FMISO, are capable of mapping regions of hypoxia in an image. However, this technology has a resolution of about 3–5 mm, whereas oxygen diffusion limits are typically 100–200 μm and thus not possible to directly resolve. The question of how we translate an intensity distribution on a PET scan to a radiotherapy dose prescription remains largely unanswered.

3In order to effectively plan irradiation of a tumour based on macroscopic scale images, it is important to understand oxygen gradients within a tumour and what impact these have on the image formation. Since the early 1900s, many mathematical descriptions of oxygen distribution in tissue have been proposed. Early mathematical models focused on the diffusion of oxygen into a cylindrical mass of tissue. This was extended to estimate the diffusion from a cylindrical blood vessel. Validating such models presents significant obstacles, as directly measuring oxygen concentration in tissue is not straightforward. As a result, most model testing has used histological sections stained for both hypoxia and vasculature to estimate the diffusion radius and compare this to the predicted limit. One of the major limitations of this approach is that the presence of vasculature does not necessarily translate to a well-oxygenated blood supply, as tumour blood vessel functionality is temporally heterogeneous and tumour microvasculature is often chaotic and poorly perfused. For this reason, in order to begin to understand the transport and metabolic dynamics of oxygen in extravascular tissue, an initially avascular and ideally in vitro experimental model is desirable. For this work, we have studied oxygen diffusion and consumption using multicellular tumour spheroids (MTS), as these fit both criteria.

4Tumour spheroids are a bridge between standard two-dimensional monolayer cultures and in vivo models. As tumours, spheroids are three-dimensional aggregates of cancer cells which, beyond the diffusion distance of oxygen, naturally form regions of hypoxia. Additionally, their signalling and metabolic profiles are more similar to in vivo cells than monolayers. However, similar to monolayers, spheroids are relatively straightforward to culture and easier to manipulate than tumours. For these reasons, spheroids have been widely used to investigate the development and consequences of hypoxia. Measuring spheroid oxygen gradients in situ is possible using pO2-sensitive electrodes. However, the physical disruption of the three-dimensional environment by the electrode itself can perturb the internal oxygen gradient. Oxygen-sensitive markers, such as EF5 and pimonidazole, are widely employed to detect hypoxia both clinically and in preclinical models, for example spheroids, without perturbing the cells. In this work, we will use EF5 to determine the extent of hypoxia in MTS owing to its well-defined oxygen-dependent kinetics.

5Early interest in MTS began in the late 1970s. Spheroid growth has been examined by several authors; in the formulation by Conger & Ziskin, spheroids are modelled as growing exponentially and then approximately linearly, achieving an approximately constant viable rim thickness. Mueller-Klieser derived a semi-analytical piecewise model of oxygen distribution that considers the oxygen gradient subject to various boundary conditions for EMT6/Ro cells. This and its variations provided an effective way to estimate the O2 distributions using oxygen probes and a series of equations for the different spheroid regions. The model is useful where the boundaries are clearly known but is not in itself predictive of how and where these boundaries occur, nor of the limits of these regions. Extensions of this model have been used in photodynamic formulations to estimate the threshold dose for the production of reactive oxygen singlets.

6Of late, there has been renewed interest in tumour spheroids and the scope for their application has increased dramatically—spheroids have been used in radiation biology as a means to test fractionation and other parameters in a controllable environment, in chemotherapy to act as a model for drug delivery and even to investigate cancer stem cells. Cancer spheroids have also shown potential as a model for exploring FDG-PET dynamics in solid tumours. In all these applications, the diffusion of oxygen and the rate at which it is consumed plays a large role and better understanding of these dynamics would be of considerable benefit. To robustly investigate oxygen consumption in MTS, we derive an analytical model which describes oxygenation through a spheroid and explicitly predicts the consumption rate for a given spheroid, as well as a prediction of the extent of the hypoxic, necrotic and proliferating regions. This model is validated here using stained cross sections of DLD1 tumour (human colorectal cancer) spheroids.

7Consider the diffusion equation, which relates the change in oxygenation with time to the spatial distribution. This is written aswhere D is the diffusion constant of the tissue, C the volume of oxygen per unit mass and a(r) the oxygen consumption rate at a point r. In previous literature on this subject, the quantity C has been defined as a concentration but this can lead to dimensional inconsistencies. In this work, we derive a robust relationship between the partial pressure and C, by the manipulation of Henry's law. The oxygen supply is assumed to be in steady state, and the tissue consumes oxygen at a rate a(r). Thus, the equation can be re-written in terms of the spherical Laplacian asThe form of a(r) is of importance in solving this equation. Initially, it is assumed that a is approximately constant. At the edge of the spheroid, a radius of ro, the volume of oxygen per unit mass is Co. At a distance of rn from the centre, C = 0 and dC/dr = 0, as illustrated in figure 1. Necrosis may set in at partial pressures slightly above zero, as cells are severely hypoxic at partial pressure of 0.8 mmHg or below, but as this is close to zero it does not affect the generality of the model. With these boundary conditions, we can solve equation (2.2) and writewhere C expresses the volume of oxygen gas per unit tumour mass. This is essentially the semi-analytical solution presented by Mueller-Klieser. In order to manipulate this equation via Henry's law, we re-express this in terms of the mass of oxygen per unit volume. Assuming that the density of the tumour, ρT is similar to that of water and the density of oxygen gas is 1.331 kg m−3, then we can convert this by simply multiplying C by Henry's law constant K varies with temperature, and for oxygen gas at human body temperature, it may be expressed as 2.2779 × 10−4 m3 mmHg kg−1. Setting Ω = K = 3.0318 × 107 mmHg kg m−3, then equation (2.3) can be expressed in terms of partial pressure asWe can further extend the model to a fully analytical form by deriving an expression for rn, the anoxic limit. By definition, C(rn) = p(rn) = 0. As ro = rn + rc, equation (2.3) may be rearranged to yield a third-order polynomial equation that characterizes the spheroid.

8There must also be a certain spheroid size where the partial pressure of oxygen just reaches zero at the centre, which is the greatest radius the spheroid can obtain and still be entirely viable. If we define this radius as the ‘diffusion limit’ rl, then solving equation (2.6) for rc = ro = rl, this is simply

9The general solution to the identity in equation (2.5) can be solved analytically with trigonometric identities. If we define χ as a function of ro so thatthen the thickness of the viable region, rc and the anoxic region rn are given simply byand The novel result of this analysis is an expression that explicitly enables us to find the viable rim and anoxic regions in a spheroid at any stage of growth from first principles. This allows some interesting insights into spheroid growth. Consider rl as defined in equation (2.6); this can be thought of as the radius for a spheroid where the oxygen partial pressure at the centre is exactly zero. It is important to note that rl does not depend on the size of the spheroid. The limit of the viable rim, rc, can be found by taking the derivative of equation (2.5) with respect to the spheroid radius ro. This is simply

10The implication of this is that the viable rim thickness decreases as the spheroid grows in size, tending towards this limit if the consumption stays approximately constant. The maximum thickness of the viable region is then rc = rl and this decreases asymptotically towards the limit rc = rm as the spheroid grows. Some previous models for MTS have assumed a constant thickness for the viable rim. A growth analysis by Burton suggested that for a spherical tumour with central necrosis, the rim thickness should tend asymptotically towards of the diffusion limit that can be reached without necrosis, which is in agreement with the form we have derived. The quantity rm is similar to the one-dimensional diffusion limit derived for a cylindrical geometry by Secomb et al.. This new formulation makes three distinct and testable predictions, which are (1) the diffusion limit rl is inversely related to the consumption rate and is independent of spheroid radius;(2) the oxygen consumption rate and oxygen partial pressure at any point in the spheroid can be estimated once rl or spheroid boundaries are experimentally derived; and(3) the viable rim thickness, rc, decreases towards a theoretical minimum of rm with increasing spheroid size.

11(1) the diffusion limit rl is inversely related to the consumption rate and is independent of spheroid radius;

12(2) the oxygen consumption rate and oxygen partial pressure at any point in the spheroid can be estimated once rl or spheroid boundaries are experimentally derived; and

13(3) the viable rim thickness, rc, decreases towards a theoretical minimum of rm with increasing spheroid size.

14There is some evidence that the consumption rate may vary with spheroid size, although the exact relationship is not entirely clear. In such a case, a constant rl would hold only if the consumption rate for a spheroid of a given size can be treated as approximately constant; the current analysis can be used to investigate this, and the findings and implications of such are examined in the Results and Discussion sections of this study. A model considering Michaelis–Menten-type form for a(r) was also developed, but as the constant km (the oxygen concentration at which the consumption rate is half its maximal value) is typically below 1 mmHg, this model form was found to have only negligible differences from the constant consumption form shown here, and has been omitted from this work.

15DLD1 (human colorectal carcinoma) cells were used to produce multicellular spheroids using the liquid overlay technique as described previously. Briefly, tumour cells were seeded in 24-well plates coated with 1% (w/v) agarose, at a final concentration of 10 000 cells per well in 200 μl growth medium (10 mM glucose). After the initial 4-day coalescent period, the growth medium was changed every 2 days. At set timepoints between 4 and 17 days after initiation, spheroids were incubated with 300 μM EF5 at 37°C for 4 h. Spheroids were then fixed in 4% paraformaldehyde prior to cryopreservation. Frozen spheroids were cut into 10 μm sections via microtome, and the section corresponding to the central cut through the diameter was selected. Sections were then dual-stained for EF5, using a Cy3-conjugated antibody, Ki-67, a marker of proliferation, using an AlexaFluor 488-conjugated antibody. Slides were counterstained with the nuclear dye Hoechst 33342. Sections were imaged at 10× magnification using a Leica epifluorescent microscope. The oxygen partial pressure outside the spheroids in solution during growth was po = 100 mmHg. The pO2 was monitored adjacent to the spheroid in 10 wells using the Oxylite oxygen sensor probe (Oxford Optronix). Although spheroids were maintained in incubators at a controlled temperature of 37.2°C, the temperature in each well was also measured using a digital thermometer (RS electronics) which indicated 36.9°C.

16Spheroids were grown as outlined, and nine dual-stained DLD1 spheroids up to 17 days of growth inclusive were obtained for analysis. Artefacts irrelevant to the image were removed manually before automated analysis; these artefacts included parts of other spheroids in the image plane and orphaned staining molecules that had become dislodged during the sectioning process. A MATLAB script then performed a reproducible analysis on each spheroid. This consisted of the following steps: — the image was segmented and the centroid was found;— the radial distance from the centroid to the closest region above a certain threshold was determined for every 1° rotation. This gave the values for rn; and— by stepping back from the image edge towards the centroid, the radial distance from the outer edge of the spheroid to the centroid was determined at 1° intervals, allowing determination of ro.

17— the image was segmented and the centroid was found;

18— the radial distance from the centroid to the closest region above a certain threshold was determined for every 1° rotation. This gave the values for rn; and

19— by stepping back from the image edge towards the centroid, the radial distance from the outer edge of the spheroid to the centroid was determined at 1° intervals, allowing determination of ro.

20The same algorithm and process was performed on each image, and the results plotted over the image to ensure accuracy. An example is shown in figure 2. This analysis ensured consistency between spheroids and reproducibility of results. Three hundred and sixty values for ro, rc and rn were averaged and the standard deviation was calculated for each spheroid image. This is important as during sectioning the spheroids can become slightly deformed and many measurements ensured the reproducibility of results as well as enabling the estimation of uncertainty. Results of segmentation were not critically dependent on choice of threshold, owing to the steep gradient in the images.

21Images were separated into their green and red colour channels to find the boundary between Ki-67 and EF5 staining. The two colour channels were both examined to find the average point of transition. The distance from the centroid to the edge of this image was calculated in 1° steps around the centroid. If the detected pixel lay beyond the spheroid boundary ro or inside the anoxic region rn, it was counted as an artefact and excluded. An example of this is illustrated in figure 3.

22The diffusion distance rl can be estimated experimentally by combining and re-arranging equations (2.5) and (2.6), yielding an equivalent expression for this value ofUsing the image analysis technique previously outlined, it is possible to measure the spheroid radius ro and the thickness of the viable region rc for any spheroid and rl can be readily computed through equation (4.1). If the oxygen consumption rate and diffusion constant stay approximately constant at all stages of growth, then values for rl should also remain approximately constant regardless of spheroid dimensions. Thus, a value for rl can be determined with the data from each image to determine whether this relationship holds.

23The approximate consumption constant a can be estimated by manipulating equation (2.6) which would yieldFor this work, we have assumed the oxygen diffusion constant D is close to that of water (D = 2 × 10−9 m2 s−1). The viable regions, rc, comprise both proliferating and hypoxic cells stained, respectively, as Ki-67 staining (green) and EF5 staining (red). The point at which these boundaries meet is the threshold level for hypoxia binding. By manipulating equation (2.4) and finding the solution to the resulting equation, it is possible to find the distance rp that corresponds to any given oxygen partial pressure p. This iswhere ϕ(p) is given byEF5 binding occurs maximally in the presence of hypoxia at partial pressures below 10 mmHg. Thus, we may assume that the partial pressure at the interface between the proliferating and hypoxic regions is approximately 10 mmHg. The distance from the spheroid centre at which this transition occurs can be estimated from the images and compared to the model prediction from equations (4.3) and (4.4) for p = 10 mmHg, using the experimentally obtained value for consumption rate a. These can be compared in order to investigate and validate the model and estimate its predictive power.

24The model predicts that the viable rim thickness of the spheroid will tend to decrease with the growth of the spheroid towards the limit rm.

25Examples of stained spheroid sections are shown in figure 4 and in the electronic supplementary material. Three hundred and sixty measurements per spheroid over nine different spheroids were obtained, yielding values for ro, rc and rn. Using equation (4.2), a single value of rl was obtained for each spheroid. These values are shown in figure 5, with bars of ±1 s.d. shown for all spheroids. Despite the significant difference in spheroid size, quality and age, the derived value was reasonably constant at 233 ± 22 μm. The approximately constant rl gives rise to an approximately constant estimate for a for all spheroids studied through the use of equation (4.2), as given in table 1, with standard deviation of ±1.4 × 10−7 m3 kg−1 s−1, or approximately 1.82 ± 0.35 × 10−3 mol kg−1 min−1.

26The distance from the centre of the spheroid to the interface between the hypoxic and proliferating region was denoted rp = 10. Using the derived value for a given in table 1, solutions to the model (equation (2.4)) for an oxygen partial pressure of 10 mmHg were obtained. These are compared to the data obtained experimentally and depicted in figure 6. Error bars of ±1 s.d. are shown along the vertical (rp = 10) axis. For clarity, the uncertainties on the x-axis (ro) are not shown on this figure. These range between 16 and 30 μm for the various data points. A regression analysis on the solid line data shown and the measured points yields a co-efficient of determination of R2 = 0.9635 for rl = 233 μm, indicating a high level of agreement between the model and experimental data.

27The model predicts a small variation of the viable radius with increasing spheroid radius; experimental data and model predictions are shown in figure 7 with the vertical bars representing 1 s.d. Again, the uncertainty on the x-axis (ro) is not shown for clarity. All the data points lie within 1 s.d. of the predicted value, and almost all lie within the envelope, indicating agreement. However, it should be noted that for the range of values of ro available, the model predicts that the absolute maximum difference in rc would be only This is less than the uncertainty in the measurements for ro and rc, making conclusive validation impossible. However, all the measured data lie between the rl and rm, as predicted by the model and underlying theory, indicating that the model predictions are not inconsistent with the experimentally derived data.

28For the DLD1 spheroids studied, rl had a relatively constantly value of 233 ± 22 μm, with the low deviation indicating that the consumption rate for all spheroids studied was approximately constant at a = 7.29 ± 1.4 × 10−7 m3 kg−1 s−1, in good agreement with previous data from similar approaches. It was important to further verify whether the model could predict the oxygen partial pressure at any point within the spheroid: as the spheroid was dual-stained, the boundary between the Ki67 and EF5 staining implied that the partial pressure at this point was approximately 10 mmHg. This could be theoretically predicted by equations (2.3) and (2.4) and compared to the value found with analysis of sectioned spheroid images. The model predictions and measured data were found to be in good agreement, with a fit of R2 > 0.96. This also provides good evidence that the constant value of rl was not merely due to chance. Finally, the prediction that rc will tend towards the limit rm as ro increases was also tested. This was done by using the experimentally derived value of rl and comparing modelled model values of rc with the measured data. All of the data points investigated lay within 1 s.d. of the model predictions. However, greater resolution would be needed to confirm this trend, as the decrease in rc is difficult to observe, being less than the uncertainty in the measured values of rc. This is a challenge to overcome, because the majority of the uncertainty in rc is due to irregularities in spheroid growth or distortion induced during sectioning. In summary, the model predictions were observed experimentally, indicating strongly that the model provides a quantitative description of the spheroid growth characteristics for this cell line and that boundaries are influenced mainly by oxygen diffusion limits.

29One advantage of the proposed model is that it allows direct prediction of the spheroid dimensions and properties from first principles, and predicts these regions well at any stage of spheroid growth. Specifically, the characteristic equation of the spheroid allows the determination of vital parameters, such as anoxic radius rn and viable rim thickness rc from the physical properties, yielding a model that is fully analytical. It also enables the prediction of oxygen profiles without invasive measurement which may perturb the oxygen gradient within a spheroid; examples of which are illustrated in figure 8. These profiles are similar in form to previously measured results and can be matched by adjusting the parameter values.

30There has been some discussion in the literature regarding whether the variation in consumption rate changes as spheroids grow. Previous work used oxygen probes in situ to estimate the partial pressure gradient and from this estimated the value of a for 15 EMT6/Ro spheroids; while the standard deviation of this dataset was only 10% of the mean value, early growth spheroids appeared to have consumption rates markedly above average. Consumption has been estimated to vary up to 50% with cell size in some cell lines, while minimal variation is seen in other cell lines. For this analysis, the variation in consumption rate should have quite an effect on the spheroid properties; as, this indicates a slowly increasing rl as consumption decreases. Such a scenario is depicted in figure 9a with simulated drops in the consumption rate of 10, 20, 50 and 75% across the dataset contrasted with constant consumption. This trend was not observed in our data. Similarly, a large decrease in a would have an equally striking impact on the interface measurements. This is depicted in figure 9b. Simulated curves can be contrasted to the measured data to determine how good the fits are and the results of this are given in table 2. The fit to measured values become increasingly poor as the variation is increased, and markedly less linear, yielding a negative R2 at 75% variation. These fits and figures indicate strongly that appreciable variation in consumption rate was not seen for the cell line, growth conditions and size range of spheroids analysed here.

31There are several reasons that might explain this discrepancy; firstly, the EMT6/Ro and DLD1 cell lines have different growth kinetics. Owing to the rapid growth of DLD1 spheroids, all of our datasets were taken at points where ro > rl with central anoxia; by contrast, the analysis of the EMT6/Ro data suggests that many of the spheroids analysed in their work had ro < rl, without central anoxia, making direct comparison between the two approaches difficult. Early work on tumour spheroids suggested that at very early growth their consumption rate is high, rapidly dropping to a constant level as spheroids grow. Other authors treat the oxygen consumption as relatively constant at normal oxygen tensions, falling rapidly at extremely low partial pressures. In both of these cases, the value a in this work can be considered as the plateau oxygen consumption rate under normal oxygen pressures. This is an avenue worthy of further exploration, and more investigations between different cell lines would be required to answer this with certainty. Finally, the methodology of this and previous approaches is quite different, as the current method takes a first-principles approach and validates using sectioned stained images.

32The oxygen consumption of a two-dimensional array of cells was measured by using a Seahorse Extracellular Flux Analyser (Seahorse Biosciences, Billerica, CA, USA). This indicated the quantity of oxygen consumed by DLD1 cells (moles of oxygen per minute per 80 000 cells) was 250 pmol min−1 per 80 000 cells. This implies that each cell consumes 1.25 × 10−18 m3 of oxygen every second. To see how this compares to our values of a, we first assume that these consumptions are equivalent and that cells have the same density as water. From this, an estimate of cell volume can be obtained. Under the assumption of these cells being approximately spherical, this yields a cellular radius of 7.42 ± 0.47 μm. To investigate how accurate this estimate was, individual DLD1 cells were suspended and their area were measured under a microscope. This yielded a mean value for suspended cell radius of 6.97 ± 0.68 μm, which is in very good agreement with estimates. The estimate of D is also a potential confounding factor; the estimate of a is dependent on D by equation (2.6) and this likely varies among different tumour types. However, as other authors have noted, it is likely to be close to water and as such the value used here is a reasonable estimate, though better determination for a given spheroid type would improve accuracy. The value of po is another potentially confounding factor—previous investigations have indicated that a diffusion depleted zone can exist up to 100 μm surrounding the spheroid where the oxygen partial pressure is below that measured in the bulk medium. To circumvent this, oxygen probe measurements were taken as close as possible to the spheroid edge and this value is taken as po. The value for rp = 10 is dependent on this measurement as seen in equation (4.4), and the good agreement between measured and predicted rp = 10 indicates that this value was a realistic estimate. As, better determination of po yields a more accurate estimation of rl and consumption rate. While the literature indicates that EF5 binds between 0.8 and 10 mmHg, a sensitivity analysis was performed to estimate the resultant errors if the upper limit for EF5 binding was significantly under- or overestimated; even at 20% error (upper limit for binding occurring at 8–12 mmHg) the projected rp = 10 difference was only 4.5 μm, indicating that potential errors introduced by binding sensitivity are small and effectively negligible.

33The current model does not consider the glucose consumption, focusing instead on the oxygen diffusion. The evidence presented here indicates that for the case of DLD1 spheroids, oxygen diffusion is by far the dominating factor in spheroid growth. A variation of the model was derived for oxygen consumption that varied with Michaelis–Menten kinetics, but this was found not to differ markedly from the model presented here for the low values of km and is consequently not detailed here. The model is time independent, but in theory it could be extended to factor in growth rate. This is beyond the scope of the present work, but despite this time independence, the model does make useful predictions about growth and growth limits. Moreover, it allows the quantification of oxygen diffusion through avascular tissue, including oxygen distributions and likely patterns of hypoxia and proliferation. While the method outlined in this work allows calculation of the anoxic boundary, the relationship between hypoxia and necrosis is complicated. Literature reviews suggest that there are cases when hypoxia and necrosis coincide, cases where necrosis sets in after hypoxia and even cases where necrosis can precede hypoxia; although this behaviour is not seen with the cell line and experimental conditions in this work, and the scope of our modelling is limited to the proliferation, hypoxia and necrosis patterns shown. For cell lines that do not readily aggregate in spheroids, alternative modalities such as multicellular layers (MCL) have shown promise in investigating oxygen distribution with a one-dimensional diffusion model, but it is currently difficult to determine predictive boundaries with this method. The spheroid results presented serve as an important validation of the underlying theory, and future work will include applying similar models to histological measurements of oxygen in tissues with known microvessel vascular maps, and contrasting this with macroscopic imaging to better quantify likely hypoxia distribution in tumours, and therefore better inform strategies for dose painting in radiotherapy.

34We present a method for estimating oxygen diffusion and consumption rate in MTS, using an analytical solution of the spherical reaction–diffusion equation. The advantage of this model is that it provides a prediction of the extent of viable, hypoxic and anoxic regions in a spheroid. It also predicts analytically how these regions will change with increasing spheroid size. We have shown that model predictions of the spherical oxygen diffusion limit, the hypoxic boundary and the viable rim thickness are all in agreement with experimental measurement using stained sections of DLD1 spheroids, and the estimated oxygen consumption rate was in agreement with measured values.