Mouse insulinoma 6 (MIN6) pancreatic beta cells were kindly provided by Dr. James D. Johnson, University of British Columbia (originally acquired from Dr. Jun-ichi Miyazaki, Osaka University). MIN6 cells in complete medium were seeded in tissue culture treated polystyrene flasks at a density of 4–8 × 10⁵ cells/cm². Complete medium consisted of Dulbecco’s Modified Eagle medium (DMEM; Gibco #10313-021, Thermo Fisher Scientific) supplemented with 10% fetal bovine serum (FBS; HyClone #SH3039602, Fisher Scientific), 1% penicillin/streptomycin (Gibco #15140-122, Thermo Fisher Scientific), 1% L-glutamine (Gibco #25030-081, Thermo Fisher Scientific), and 0.1% 2-mercaptoethanol (Fisher Chemical #O3446I, Fisher Scientific). Cells were cultured in an aseptic incubator at 37°C and 5% CO₂. Medium changes were performed every 48 h. Cells were passaged at 90% confluency using 1X TrypLE Express Enzyme (Gibco #12605-028, Thermo Fisher Scientific).

Adherent MIN6 cells were passaged twice before seeding in AggreWell^™400 24-well plates (#34415, STEMCELL Technologies) to generate islet-like clusters (ILCs). The following seeding densities were used: 2.4 × 10⁵ cells/microwell (200 cells/ILC), 6 × 10⁵ cells/microwell (500 cells/ILC), 1.2 × 10⁶ cells/microwell (1,000 cells/ILC), 2.4 × 10⁶ cells/microwell (2,000 cells/ILC), 3.6 × 10⁶ cells/microwell (3,000 cells/ILC). Cells were cultured for 48 h in complete medium at 37°C and 5% CO₂. Cell seeding and ILC harvesting steps were performed according to the AggreWell^™400 manufacturer’s protocol.

ILCs harvested from each macrowell of an AggreWell^™400 plate were immediately transferred to a new well in a 48-well plate (∼1,200 ILCs/well). A live/dead staining solution was prepared with final concentrations of 5 µM calcein-AM (Invitrogen #C1430, Thermo Fisher Scientific) and 5 µM ethidium homodimer-1 (#ab145323, Abcam) in phosphate buffered saline (PBS; Gibco #21300025, Thermo Fisher Scientific). To each well, 150 µL of staining solution were added. Samples were incubated in the dark at room temperature for 30 min and then imaged using an Olympus IX81 inverted system microscope and MetaMorph^® acquisition and analysis software.

ILCs from each of the five seeding density conditions described above were harvested from one macrowell of an AggreWell^™400 plate and immediately transferred to a 15-mL centrifuge tube and allowed to settle (∼3,000 ILCs/tube). Once settled, the ILCs were resuspended in 600 µL of DMEM and divided evenly into 1.5 mL tubes, each containing 200 µL (∼1,000 ILCs/tube). The OCR was determined as per previously reported methods (Papas et al., 2007a; Papas et al., 2007b). Briefly, the OCR was measured on a calibrated Instech Fluorescence Lifetime Micro Oxygen Monitoring System (Instech Laboratories, Plymouth Meeting, PA; FOL/C2T175PFOL) using Amplifier NeoFox software. Immediately before measuring, the DMEM was removed from the ILCs and replaced with 200 µL of serum-free DMEM. The sample was mixed to ensure a homogeneous solution and 170 µL was transferred to each reading chamber where the oxygen partial pressure (pO₂) was monitored. Once the measurement was complete, the ILCs were carefully collected and stored at −80°C in preparation for DNA quantification.

ILCs were lysed by performing a freeze-thaw cycle and sonicating for 15 s on ice. DNA quantification was performed using the Quant-iT™ PicoGreen™ dsDNA Assay Kit (Life Technologies, Carlsbad, CA; Cat #P7589) according to the manufacturer’s protocol.

The OCR calculation was performed as previously reported by Papas et al. by determining a linear region of the measured pO₂ versus time plot and calculating the slope (Papas et al., 2007b). This was then normalized to the DNA content to give an OCR/DNA ratio.

A stock solution of 2.5% alginate (Manugel GHB, IFF Nutrition & Biosciences (formerly FMC Biopolymer); 60–65% guluronic acid; 1% solution viscosity 50–100 mPa s) and a stock suspension of 0.5 M calcium carbonate (CaCO₃; Avantor #1301-01, VWR) were both prepared in pH 7.4 HEPES-buffered saline (Fisher BioReagents #BP310 and #BP3581, Fisher Scientific) and autoclaved at 121°C for 30 min. Glucono-δ-lactone (GDL; Sigma-Aldrich #G4750, MilliporeSigma) was dissolved in HEPES-buffered saline and sterilized with a 0.2-µm nylon syringe filter (Fisherbrand™ #09-719C, Fisher Scientific) immediately before use. Using the above reagents, an alginate hydrogel mixture was prepared with final concentrations of 2% alginate, 30 mM CaCO₃, and 60 mM GDL. A 9:1 volume ratio of alginate-to-medium was maintained.

MIN6 ILCs were generated using AggreWell^™400 plates as previously described and added to the gel mixture. A seeding density of 1.2 × 10⁶ cells/macrowell was used to achieve 1,000 beta cells/ILC. GDL was mixed in immediately prior to casting to trigger internal gelation with a final concentration of 1,000 ILCs/mL gel. The gel mixture was loaded via syringe into a cylindrical perfusion device outfitted with a 2-mm diameter stainless steel rod as a vascular template. The device was incubated at room temperature for 60 min for adequate gelation. The vascular template was then removed, leaving a central hollow channel, and the device was connected to a tubing loop for perfusion culture. Device assembly is described in detail in our previous work (Fernandez et al., 2021). A general schematic of the device geometry is provided in Figure 3A.

Perfusion of culture medium was established using a peristaltic pump set up inside an incubator at 37°C. O₂ levels were set to 2% (hypoxia), 10% (arterial), or left at ambient 18.6% (in vitro normoxia); CO₂ was maintained at 5% and N₂ made up the balance (Wenger et al., 2015). Hypoxic and arterial conditions were established using a BioSpherix ProOx C21 system to control CO₂ and O₂. In vitro normoxic conditions were established using a regular cell culture incubator with CO₂ control. Culture medium was pre-incubated at the desired O₂ tension for 24 h before beginning perfusion. A culture medium flow rate of 3 mL/min was selected to mimic laminar blood flow through a medium to small artery 2 mm in diameter (Reynolds number Re < 50). The medium was drawn from a 30-ml reservoir of DMEM and perfusion was run for 48 h.

Fixed samples were sent to a core facility of McGill University for paraffin embedding, sectioning, deparaffinization and staining. Sections of 6-µm thickness were taken orthogonal to the perfusion channel. Per device, 8 stained sections were acquired: 2 × hematoxylin and eosin (H&E) for general tissue visualization, 2 × anti-cleaved caspase-3 (CC3; New England BioLabs, 9,661) to visualize apoptotic cells, 2 × anti-hypoxia-inducible factor 1-alpha (HIF-1α; Novus Biologicals, NB100-479) to visualize oxidative stress, and 2 × anti-Ki67 (clone B56) (#3168022D, Fluidigm Corporation) to visualize proliferative cells. Ki67-stained samples were omitted from analysis due to patchy staining results. Imaging was performed using a Zeiss Axio 200M automated inverted microscope and Zen Pro acquisition and analysis software. Images were taken at 20X magnification and stitched to assemble a complete cross-section of the sample.

COMSOL Multiphysics 5.5 and MATLAB™ 2021b were used to develop and render models simulating the perfusion device. The model uses two sets of physics options: 1) laminar flow and 2) transport of dilute species (oxygen). The boundary conditions prescribed for 1) are: laminar flow at the inlet with a known flow rate; zero pressure at outlet (results in gauge pressure profile). For (2), the boundary conditions are: a fixed inlet concentration of oxygen; advection dominant transport at the outlet (pre-programmed in COMSOL). All models were solved in COMSOL and post-processed in MATLAB.

The perfusion device geometry was replicated in the COMSOL interface and assigned to different materials (cell culture medium or cell-laden alginate). The physical properties of the cell culture medium were assumed to be the same as water at 37°C. The cell-laden alginate was modelled as a porous medium domain with modified Michaelis-Menten oxygen consumption kinetics supported and developed in current literature as shown in Table 1 (Buchwald, 2011).

The system was solved as one fully coupled model in COMSOL using a convection-diffusion model (Eq. 1) with a reaction term (Eq. 2) representing oxygen consumption. This approach is consistent with analytical models without advection (Dulong and Legallais, 2007) or using oxygen partial pressure instead of oxygen concentration in the convection-diffusion equation (Komatsu et al., 2017). − D e f f ∇ 2 c O 2 + u ⋅ ∇ c O 2 = R (1) R = [ − R m a x , O 2   c O 2 c O 2 + K M M , O 2     δ ( c O 2 > c c r ) ] X (2) Where D _eff : effective diffusivity of cell-laden alginate; c _O2 : concentration of oxygen; c _cr : critical concentration of oxygen below which cell necrosis occurs; u: velocity field; R: oxygen consumption; R _max,O2 : maximum oxygen consumption rate; K _MM,O2 : Michaelis-Menten coefficient for oxygen diffusion; δ: Heaviside function; and X: cell fraction.

All single channel models were solved in a 2D axisymmetric system while 2-channel systems were solved in 3D. The main controlled parameters of the simulation are the flow rate Q, cell fraction X, and channel geometry. The model assumes rigid device walls, Newtonian fluid, fully developed laminar flow, no thermal diffusion, no pressure loss, no backflow, negligible body forces, and a homogeneous cell fraction throughout the hydrogel. A “normal” mesh size was used for an optimal trade-off between computation time and visual quality. A relative tolerance of 1 × 10⁻⁶ (0.0001% error) was selected as the convergence criterion for all oxygen maps. A relative tolerance of 1 × 10⁻⁵ (0.001% error) was selected for the parametric analyses to reduce computation time. The solving method employed was a tolerance-based direct solver. The Heaviside function (δ) serves to set the oxygen consumption rate to zero in regions below the necrotic oxygen partial pressure (c _cr ; see Table 1), it also helps in avoiding negative concentrations and divergence (Komatsu et al., 2017).

All assumed model parameters are listed in Table 1. The oxygen consumption rate R _max was either found experimentally for MIN6 or taken from literature for human islets. The cell fraction X was calculated using a volume fraction V _cells /V _total assuming islets are spheres with radius r. The effective diffusivity of cell-laden alginate D _eff was calculated as a weighted average of the diffusivity in cell tissue D _tissue and in alginate D _alg (Eq. 3). The Henry’s constant used is for water at 37°C. D e f f = D t i s s u e X + D a l g ( 1 − X ) (3)

An analysis of the model parameters listed in Table 1 was performed using a parametric sweep in COMSOL Multiphysics. Default parameter values were assigned to reflect the experimental conditions used during the perfusion studies described above. Individual parameter values were varied over a selected range and plotted at six points—three within 1 mm of the device entry and three within 1 mm of the exit. These points were located at different radial positions between the perfusion channel and the gel periphery (0, 250 μm, 500 µm). A “coarser” mesh size was applied near the gel periphery and a “normal-fine” mesh was applied closer to the perfusion channel to resolve the kinetics. These points also facilitate an assessment of the effects of channel length on theoretical oxygenation.

MIN6 ILC sizes were measured using the manual freeform selection tool in ImageJ software. Size distributions, volume mean moment, and standard deviations were calculated in Microsoft Excel and plotted in GraphPad Prism 5.01.

MIN6 viability was analyzed by measuring fluorescence intensity with Fiji software. Mean values and standard deviations were calculated in Microsoft Excel and plotted in GraphPad Prism. One-way analysis of variance (ANOVA) and the Tukey multiple comparisons test were performed in GraphPad Prism to compare the means of the sample conditions (α = 0.05).

Oxygen consumption rates normalized to DNA content were calculated in Microsoft Excel. One-way ANOVA was performed in GraphPad Prism to compare the means of the sample conditions (α = 0.05).

Due to large file sizes, stitched histology images were resized to 50% and exported from Zen Pro software. ILC size, distance, and mean grey value measurements were acquired using Fiji software and analyzed in Microsoft Excel and GraphPad Prism. Statistical analysis of mean grey values was performed using a one-way ANOVA followed by the Tukey multiple comparisons test (α = 0.05).

To determine appropriate ILC seeding conditions for perfusion device studies, MIN6 cells were seeded at varying densities in 24-well AggreWell™400 plates and visualized using live/dead staining and fluorescence microscopy. Size distributions and representative staining images are shown in Figure 1. The bar graphs in Figure 1 indicate that the ILC size distribution narrows as the seeding density decreases. The volume mean moment (De Brouckere mean diameter) D [4,3] calculated for each condition ranges from 208 to 306 µm. Conditions A (3,000 cells/ILC) and B (2,000 cells/ILC), exhibit broad size distributions, with >14% of total ILCs measuring <100 µm in diameter. ILCs in Conditions A and B are also visibly irregular in shape when contrasted to imaging results from Conditions C–E. This is likely due to an overflow of excess cells in the well plates, causing inadequate aggregation and irregular clumping. D [4,3] decreases from Condition B through Condition E, as does the corresponding standard deviation.

Overall cell survival for each condition was quantified using live/dead staining images. Adherent MIN6 controls were also stained and analyzed (images not shown), resulting in a mean viability of 96.8%. For MIN6 ILCs, the mean viability ranges from 84.1 to 96.6%. No significant difference is observed between the viability of the control and the ILCs at Conditions A, B, and E (p > 0.05). Significant differences are observed between the control and Conditions C and D (p < 0.001). Viability is observed to be higher in Conditions A and B likely because of the broad size distribution and relatively high proportion of small ILCs (<150 µm diameter) experiencing high mass transport. Similarly, viability is high in Condition E and is attributed to the uniformly smaller ILCs (D [4,3] = 208 µm).

To further assess appropriate ILC seeding conditions, the oxygen consumption rate (OCR) was measured for Conditions A–E and normalized to DNA content as presented in Figure 2. A breakdown of the numerical results is also available in Supplementary Table S1. Single-cell MIN6 controls were also analyzed. Using a one-way ANOVA, no significant difference was observed between the means of the conditions, although it was noted that the standard deviations were smallest for Conditions B–D. From the literature, human research islets cultured with 10% fetal bovine serum exhibit an OCR/DNA of 202 (SD 87) nmol/min/mg (Papas et al., 2007a). The OCR values measured for MIN6 ILCs are therefore 4–6 times higher than in human research islets.

Given the size distribution, viability, and OCR data collected, Condition C was selected as the most appropriate ILC seeding density for future studies. The density of 1,000 beta cells per aggregate best approximates the beta cell content of native human islets (Pisania et al., 2010) while the relatively narrow size distribution facilitates computational modelling parameters. Moreover, the indication of central necrosis occurring in the ILCs (Figure 1C, right) increases the likelihood of observing a measurable response to hypoxia during oxygen studies. Mean OCR normalized to DNA content for Condition C is 780 (SD 261) nmol/min/mg. At a low concentration of 1,000 ILCs/mL gel, this corresponds to an equivalent of approximately 70,000 human islets per device, with each device having a capacity of 18 mL cell-laden gel.

To predict the oxygen profiles within the perfusion device, a computational model was developed simulating the device geometry, which is a cylinder with length 110 mm, diameter 12.7 mm, and central perfusion channel diameter 2 mm (Figure 3A). Parameters were obtained both from literature and from experimental results as summarized in Table 1.

To assess the effect of the major parameters on the final oxygen concentration profiles, parametric analyses were performed as shown in Figures 3B–E. Points were selected at the device entry and exit at radial distances of 0, 250, and 500 µm from the perfusion channel to observe any notable effects at different positions. Vertical dashed lines on the plots indicate ranges or points of interest with respect to each parameter. Both cell fraction X (Figure 3C) and maximum oxygen consumption rate R _max (Figure 3D) have a significant impact on device oxygen tensions. This occurs within the range of interest for islet transplantation applications. For example, changing from a low cell dose to a therapeutic dose may shift the oxygen tension by several orders of magnitude. Similarly, a change in cell type from human islets to MIN6 causes an overall decrease in oxygen tension due to an increase in oxygen consumption.

Large changes in diffusivity D may cause large changes in oxygen tension (data not shown). However, in practice, hydrogels such as alginate tend to have similar diffusivity values at the concentrations used for cell encapsulation (Mehmetoglu et al., 1996; Buchwald, 2011). Therefore, experimental changes in alginate concentration are not expected to significantly affect device oxygen profiles. For changes in X and R _max, notable differences are observed between the trends at the device entry and exit, especially with increasing total oxygen consumption in the device. A major difference between entry and exit is observed for changes in flow rate Q at low values (Figure 3E). At the experimental flow rate of 3 mL/min, we would anticipate a slight decrease in oxygen tensions along the length of the device. Larger decreases are anticipated at flow rates below 1 ml/min. A similar plot representing the parametric analysis of the dimensionless Reynolds number Re is found in Supplementary Figure S1.

Given these results, immediate applications of this model would focus on the effects of cell type (i.e., OCR) and cell concentration, which significantly influence the oxygen distribution. Future work should further investigate longitudinal profiles within the device.

Based on the results shown in Figure 1 and Figure 2, MIN6 cells were seeded at Condition C (1,000 cells/ILC), encapsulated in alginate at a concentration of 1,000 ILCs/mL gel, and cast in perfusion devices. Devices were cultured under DMEM perfusion for 48 h at normoxic (140 mmHg or 18.6% O₂), arterial (76 mmHg or 10% O₂), and hypoxic (15 mmHg or 2% O₂) environmental conditions. Samples were then fixed and stained for hematoxylin and eosin (H&E; general tissue visualization), anti-cleaved caspase-3 (CC3; apoptosis), and anti-hypoxia inducible factor-1α (HIF-1α; response to hypoxia). Using H&E-stained sections, ILC morphology was categorized as Intact, Early Disruption, or Advanced Disruption, as illustrated in Figure 4A. Intact (Int) ILCs are smoothly rounded and lack any visible breakage. Early Disruption (ED) samples exhibit limited breakage, typically visible at the ILC periphery or indicated by small holes in the ILC interior. Advanced Disruption (AD) is characterized by pervasive ILC breakage including large interior holes, elongated cells, and peripheral cell debris. Intact ILCs were analyzed for CC3 and HIF-1α staining intensity (Figure 4B). While no statistical differences were observed between conditions, the HIF-1α results suggest a possible increasing trend in positive staining with decreasing oxygen tension. Increased expression of the HIF-1α transcription factor would indicate a response to hypoxia comprising an upregulation in processes such as angiogenesis to promote oxygen delivery.

Intact (Int) and Early Disruption (ED) ILCs were analyzed with respect to distance from the perfusion channel (Figure 4C). Advanced Disruption (AD) samples were excluded as the extreme breakage prevented reliable measurements. Under normoxic conditions, Int and ED ILCs are distributed over a notably larger distance (>5 mm) than observed under arterial and hypoxic conditions (both <2 mm). Under hypoxia, Intact ILCs were rarely observed and severely limited to <0.6 mm from the channel. This aligns with the computational simulation in Figure 4D, which illustrates the radial oxygen distribution under each oxygen condition. In each simulation, a white contour indicates an oxygen tension of 7 mmHg (∼1% O₂). As impairments in beta cell function have been observed below 7 mmHg (Papas et al., 1996), this value has been used as a lower threshold in previous models (Fernandez et al., 2014). Contrasting normoxia, arterial, and hypoxia simulations, the approximate position of this contour shifts from 0.45 mm to 0.25 mm to ∼0 mm from the perfusion channel, respectively. Similarly, a red contour representing 0.1 mmHg O₂ shifts from 0.66 mm to 0.5 mm to 0.2 mm from the perfusion channel, respectively. Current literature suggests a threshold of 0.1 mmHg O₂ for necrosis, implying islets beyond this contour would be dead and no longer consuming oxygen (Johnson et al., 2009; Buchwald, 2011). For all experimental conditions in Figure 4C, beyond the range of ED ILCs, the remaining distance to the device periphery was populated solely with AD ILCs (not quantified due to pervasive breakage).