Characterization of three‐dimensional bone‐like tissue growth and organization under influence of directional fluid flow

The transition in the field of bone tissue engineering from bone regeneration to in vitro models has come with the challenge of recreating a dense and anisotropic bone‐like extracellular matrix (ECM). Although the mechanism by which bone ECM gains its structure is not fully understood, mechanical loading and curvature have been identified as potential contributors. Here, guided by computational simulations, we evaluated cell and bone‐like tissue growth and organization in a concave channel with and without directional fluid flow stimulation. Human mesenchymal stromal cells were seeded on donut‐shaped silk fibroin scaffolds and osteogenically stimulated for 42 days statically or in a flow perfusion bioreactor. After 14, 28, and 42 days, constructs were investigated for cell and tissue growth and organization. As a result, directional fluid flow was able to improve organic tissue growth but not organization. Cells tended to orient in the tangential direction of the channel, possibly attributed to its curvature. Based on our results, we suggest that organic ECM production but not anisotropy can be stimulated through the application of fluid flow. With this study, an initial attempt in three‐dimensions was made to improve the resemblance of in vitro produced bone‐like ECM to the physiological bone ECM.

complex bone ECM structure. The recapitulation of the physiological bone ECM in vitro has received too little attention from researchers (de Wildt et al., 2019).
Tissue engineering strategies are nowadays increasingly applied for the creation of in vitro models of healthy or pathological bone, aiming at improving preclinical treatment development while addressing the principle of reduction, refinement, and replacement of animal experiments (3Rs) (Balls et al., 1995;Caddeo et al., 2017;Holmes et al., 2009;Owen & Reilly, 2018). Changes in bone's ECM composition and organization are characteristic for bone pathologies like osteoporosis, osteogenesis imperfecta, and bone metastasis (Bala & Seeman, 2015;Bishop, 2016;Karunaratne et al., 2016;Matsugaki et al., 2021).
Therefore, in vitro models that aim at studying changes in bone ECM under the influence of treatments would benefit from improved control over organic matrix formation and subsequent mineralization (de Wildt et al., 2019;Hadida & Marchat, 2020). As the organic bone ECM with mainly collagen type 1 functions as a mineralization template (Wang et al., 2012), the improvement of collagen network organization and density may enhance the biomimicry of in vitro produced bone ECM. However, the mechanism by which collagen forms a dense anisotropic or lamellar network in vivo is poorly understood. It is well accepted that in vivo bone morphology and mass is regulated by osteocytes which sense interstitial fluid flow through their lacuna-canalicular networks (Burger & Klein-Nulend, 1999). In vitro, studies mainly focus on manipulating the preferred orientation of collagen producing osteoblasts to induce tissue anisotropy. This might be accomplished by mechanical loading like cyclic stretch or directional fluid flow (Matsugaki et al., 2013;Xu et al., 2018).
These studies are, however, mainly performed in a controlled but simplified two-dimensional (2D) environment for a short period of time, not representative for the in vivo situation and ignoring tissue formation outcomes. In three-dimensional (3D) systems, fluid flow has been demonstrated to stimulate bone-like tissue growth including collagen formation (Akiva et al., 2021;Melke et al., 2018;Vetsch et al., 2016). One challenge in these 3D environments is that both increased mass transport and wall shear stress (WSS) as a result of fluid flow could have an effect (Hadida & Marchat, 2020;Wittkowske et al., 2016).
Cell and tissue anisotropy is then often observed in the tangential or circumferential direction of a pore (Bidan et al., 2012). Thus, to promote in vitro bone-like tissue growth and anisotropy, fluid flow, and curvature are likely two important factors. Therefore, in this study, bone-like tissue growth and anisotropy were evaluated in a concave channel in a 3D silk fibroin (SF) scaffold statically or under influence of directional fluid flow. It was assumed that when cells were oriented in the longitudinal direction of the channel, curvature was considered to be neglectable, while if cells aligned in the tangential direction of the channel, curvature was −0.67 mm −1 . Osteogenically stimulated human bone marrow-derived mesenchymal stromal cells (hBMSCs) were used because of their ability to proliferate and differentiate into osteoblasts and osteocytes (Akiva et al., 2021), making them an appropriate candidate for (personalized) in vitro bone models . As the cell response to fluid flow or curvature might change during osteogenic differentiation (Callens et al., 2020;Wittkowske et al., 2016), tissue growth and organization were studied over a period of 42 days with intermediate time points at Day 14 and Day 28 (Figure 1). In addition, before experiments were executed, computational simulations were performed. Fluid flow patterns and fluid shear stress magnitude at the channel wall were simulated with a computational fluid dynamics (CFD) model to (i) determine the optimal bioreactor settings for osteogenesis and bone-like tissue formation, and (ii) ensure only fluid flow at the channel wall in the longitudinal direction (i.e,. direction of the flow) to minimize an effect of mass transport and flow in the radial (i.e., into the scaffold) direction.
Lyophilized SF was dissolved in hexafluoro-2-propanol (003409, Fluorochem) at a concentration of 17% (w/v) and casted in cylindrical scaffold molds filled with NaCl granules with a size of <200 µm as templates for the pores. Molds were covered and after 3 h, covers were removed from molds, and hexafluoro-2-propanol was allowed to evaporate for 7 days whereafter β-sheets were induced by submerging SF-salt blocks in 90% MeOH for 30 min. SF-salt blocks were cut into discs of 3 mm height with an Accutom-5 (04946133, Struer). NaCl was dissolved from the scaffolds in ultra-pure water, resulting in porous sponges. These sponges were punched with a 9 mm diameter biopsy punch for the outer dimensions and a 3-mm diameter biopsy punch for the central channel with a fixed curvature of −0.67 mm −1 . The dimensions of the channel are based on previous research in which a 3-mm channel remained open over a period of 42 days (Vetsch et al., 2016), which is essential to enable studying the influence of directional fluid flow. Scaffolds were sterilized by autoclaving in phosphate-buffered saline (PBS) at 121°C for 20 min.

| Scaffold geometry
To obtain the detailed geometry of the scaffold, a microcomputed tomography scan (µCT) was acquired with a µCT100 imaging system (Scanco Medical). Scanning was performed for air-dried scaffolds with an isotropic voxel size of 3.9 µm, energy level of 45 kVp, intensity of 200 µA, integration time of 300 ms, and with twofold frame averaging.
To reduce part of the noise, a constrained Gaussian filter was applied with a filter support of 1 and a filter width sigma of 0.8 voxel. Filtered images were segmented at a global threshold of 55% of the maximum grayscale value. Unconnected objects smaller than 50 voxels were removed through component labeling. A distance transformation function was used to determine the pore size distribution at four regions of interest (location S1, S3, S5, and S7 of Figure 2a).

| Microscale: Scaffold permeability calculation
As the region of interest was the scaffold's central channel, the porous region around the channel was homogenized using the permeability that was determined based on the reconstructed geometry from the µCT scan as previously described (Zhao et al., 2019). This was done to reduce the excessively high computational cost caused in modeling the irregular micro-struts of the porous SF scaffold. The homogenized scaffold permeability was determined based on the geometry of 8 representative volumetric elements (RVEs) with a diameter of 500 µm (>4 times average pore size) and a height equal to the height of the scaffold (Figure 2a). The RVEs were selected from the total scaffold geometry that was reconstructed using Seg3D software (University of Utah). The fluid domain of RVEs meshed using the same strategy as in (Zhao et al., 2019) with global maximum and minimum element sizes if 20 µm and 0.2 µm, respectively.
The RVEs' permeability was determined from Darcy's law (Equation 1): where p ∆ is the pressure drop over the scaffold height H determined by solving the CFD model for each RVE; Q is the prescribed flow rate, A the cross-sectional area to the flow, μ the dynamic viscosity of the culture medium (μ = 1.09 mPa•s for cell culture medium (Maisonneuve et al., 2013)), and κ the permeability.
In the CFD model, the medium flow was defined as incompressible Newtonian and described by the Navier-Stokes equation that was solved using the finite volume method (FVM). As a convergence criterion it was required that the root-mean-square residual of the mass and momentum was smaller than a fixed threshold set at 10 −4 .

| Macro-scale: Wall shear stress calculation
The macro-scale model representing the full scaffold and perfusion bioreactor domain was used for the fluid shear stress calculations on the scaffold channel wall (Figure 2c). In this macro-structural model, the where q is the Darcy velocity and p the pressure.
The remaining central bioreactor channel without tissue formation was defined as free fluid (incompressible, Newtonian laminar flow), and described by the Navier-Stokes equations: where v is the fluid velocityvector and ρ the medium density (ρ = 1000 kg/m 3 ).
The top and bottom surfaces of the porous media domain were defined as boundaries with continuity of mass flux, and the scaffold internal channel wall was defined as a non-slip wall boundary as it is assumed that there will be only minor fluid following across the wall.
At the inlet of the CFD model, a constant flow rate of 1.5 mL/min was prescribed according to the experimental condition. At the outlet, a relative pressure of 0 Pa was applied. The macro-scale CFD model was meshed with 1,720,090 tetrahedral elements, solved by FVM using ANSYS CFX (ANSYS, Inc.) and the same convergence criteria as for micro-scale model.
To check the assumption that only little fluid flows through the interface between the channel and the scaffold, a CFD model in which a scaffold with idealized cubic pore shape and uniform pore size, was utilized. As this model only served to support the assumptions made for the main macro-model as described above, details on this model can be found in Supporitng Information (Section 2).

F I G U R E 2
Multi-scale CFD model parameters. (a) Scaffold geometry was obtained with µCT scanning whereafter the reconstructed scaffold was discretized into eight RVEs for permeability determination using the micro-model. (b) Calculated permeability of each scaffold RVE. (c) Geometry and boundary conditions of the CFD model. CFD, computational fluid dynamics; µCT, microcomputed tomography; RVE, representative volumetric element.

| Biochemical content analyses
To

| Statistical analyses
Statistical analyses were performed, and graphs were prepared in GraphPad Prism (version 9.3.0, GraphPad) and R (version 4.    (Figure 4f). At the channel wall, no differences between statically and dynamically cultured scaffolds were found. This suggests that the observed thickening in cell layer is the result of ECM production by cells at the channel wall, rather than cell proliferation.

| Organic matrix growth and mineralization
Collagen deposition was visualized in the horizontal and vertical plane ( Figure S4). Vertical plane images revealed collagen formation after picrosirius red staining through the entire scaffold for both statically  (Figure 6a,b). Although nonsignificant, more mineralization seemed present in statically cultured scaffolds.
As we were mostly interested in the scaffold channel wall, the scaffold channel volume was also analyzed for the presence of mineralization. Interestingly, at the channel wall differences between statically and dynamically cultured scaffolds could not be observed.

| Cell and tissue organization
Over the culture period progression, no clear trend in cell and tissue organization was observed for statically and dynamically cultured scaffolds (Figure 7a,c). From the actin fiber distributions, no consistent influence of directional fluid flow was observed (Figure 7b,d). On Day 28, cells tended to align more in the tangential direction of the channel for both statically and dynamically cultured scaffolds. This was, however, not consistent for all scaffolds in the dynamically cultured group.

| DISCUSSION
With the transition in the application of bone tissue engineering strategies from bone regeneration to 3D in vitro models, the challenge to create an organized bone ECM has been identified Although the mechanism by which bone ECM gains its dense and organized structure is not fully understood, mechanical loading like directional fluid flow and curvature (especially concavities) have been identified as potential contributors. In this study, we aimed at evaluating 3D cell and tissue growth and organization in a concave channel with and without directional fluid flow stimulation over a period of 42 days to include the contribution of cell differentiation.
As a result, directional fluid flow improved organic ECM growth but not organization. After 28 days of culture, when osteogenic differentiation of the cells was likely accomplished, they tended to have a small preference for orienting themselves in the tangential direction of the channel. Even when fluid flow was applied in the perpendicular direction, most samples showed cells with a preference for alignment in the tangential direction of the channel which might be attributed by its curvature (Bidan et al., 2012).
In this study, a CFD model was used to calculate the WSS conditions might be considered static (i.e., limited mass transport).
Although bone cells are highly sensitive to fluid flow (Wittkowske et al., 2016), differences in collagen formation between statically and dynamically cultured scaffolds were only observed at the channel wall, which might be explained by the relative static conditions within the scaffold. Another assumption for determining the WSS on the cells was their attachment. Only if cells have a flat attachment to the channel wall, WSS on cells is comparable to the calculated WSS.
When cells bridge pores, fluid flow not only induces shear but also strain (Mccoy & O'Brien, 2010). After seeding, cells indeed bridged the pores at the channel wall. As such, the calculated WSS magnitude might have been an underestimation of the by the cells experienced mechanical load. Cells also covered the channel wall already directly after seeding. It is expected that once they produce ECM, the irregular channel wall gets covered with a more homogeneous tissue layer in which cells experience less strain and are stimulated by mostly WSS (Hadida & Marchat, 2020). However, substantial tissue growth in the channel will likely also change the fluid flow-induced shear stress (Zhao et al., 2020). Thus, the fluid flow-induced mechanical load is expected to change over time which hinders the interpretation of the obtained results and therefore is a limitation of the present study. In this study, tissue growth and mineralization were already monitored over the entire culture period. Future studies would benefit from including these tissue growth and mineralization parameters in their models to get a more realistic estimation of the change in stress over time and potentially adapt the input flow accordingly (Giorgi et al., 2016;Hadida & Marchat, 2020).
Interestingly, while other studies have reported increased mineralization under the influence of WSS (Akiva et al., 2021;Melke et al., 2018), in our study, an opposite effect was observed. Studies with a similar set-up have also found more mineralization in statically cultured scaffolds than dynamically cultured scaffolds (Vetsch et al., 2016;Vetsch et al., 2017). In the used perfusion bioreactor set-up, only half of the medium volume can be replaced, whereas in statically cultured bioreactors all the medium can be replaced. To account for this, osteogenic supplements were added in a double concentration under the assumption that they are either consumed or degraded before the next medium change. However, this way of medium replacement might also induce a difference between the groups in protein concentration derived from FBS or in soluble factors produced by the cells (Vis et al., 2020). Recently, the impact of alkaline phosphatase in FBS on mineralization has been shown (Ansari et al., 2022). We, therefore, suggest that the difference in mineralization is attributed to the bioreactor system and its practical limitations, something that needs to be considered for future experiments using this bioreactor system. At the channel wall, differences in mineralization between statically and dynamically cultured constructs were absent. This might indicate that in dynamically cultured constructs, mechanically stimulated cells at the channel wall contributed more to mineralization than cells within the scaffold that likely sensed no to limited shear stress. To confirm this hypothesis, additional research is needed in which the medium replacement method is similar for both the static and dynamic condition.
In our effort to improve cell and tissue organization in 3D, directional fluid flow was applied in a concave channel. Fluid flow has been shown to stimulate cellular alignment in 2D (Xu et al., 2018), while (mainly concave) curvature has been shown to induce anisotropic collagen formation in 3D (Bidan et al., 2012;Bidan et al., 2013;Callens et al., 2020). By applying fluid flow in the longitudinal direction of a concave channel in a 3D scaffold, we attempted to identify the potential driver of 3D cell and tissue organization. When cells were oriented in the longitudinal direction of the channel, curvature was considered to be neglectable, while if cells aligned in the tangential or circumferential direction, curvature was −0.67 mm −1 . Over the entire culture period, no clear influences of directional fluid flow were observed. Only a small preference for the tangential channel direction was observed after 28 days for both statically and dynamically cultured scaffolds which might be attributed to the channel curvature.
To confirm that curvature was indeed the main cue for cell and tissue organization, additional research is needed in which curvature magnitude and channel shape are varied under controlled fluid shear stress. However, challenges to studying this in a 3D environment are (i) an irregular channel wall, (ii) the differentiation state of the cells, and (iii) the channel diameter and thus the curvature magnitude. First, while in this scaffold, the smallest possible pores were produced to maximize the channel-to-pore size ratio, scaffold pores might still have induced small and local changes in curvature, which could have locally influenced initial cell orientation. This might also explain why after 28 days cells tended to orient into the tangential direction of the channel, as once cells have formed a monolayer and produced their own ECM, curvature becomes a more dominant factor than local scaffold properties (Kommareddy et al., 2010). However, one would then also expect to see cell alignment in constructs cultured for 42 days which could not be detected in this study. Second, previous research has shown that undifferentiated hBMSCs prefer to avoid curvature and would therefore align in the longitudinal direction of the channel (Callens et al., 2020). Therefore, cells might have changed their orientation during their differentiation process. Third, bone-like tissue growth is mainly stimulated with higher concavities (Vetsch et al., 2016).
However, to avoid closing of the channel, which would have induced unpredictable flow patterns, a trade-off between curvature and channel diameter had to be made.
Recently, prostaglandin E2 signaling by osteocytes was identified as a potential inducer of osteoblast alignment (Matsuzaka et al., 2021).
Within this setup, studying the contribution of osteocytes to ECM anisotropy might be difficult, as osteocyte differentiation from hBMSCs is highly challenging and requires specific loading conditions and long-term culture (Akiva et al., 2021;Zhang et al., 2022).
Osteoblast-osteocyte co-cultures might be performed to overcome

| CONCLUSION
In the present study, we presented a computationally informed 3D model for bone-like tissue growth. In our attempt to improve ECM density and anisotropy, cell organization and tissue growth were evaluated under influence of directional fluid flow in a concave channel. Based on the results obtained within this study, we suggest that organic tissue production but not anisotropy can be stimulated through the application of directional fluid flow. As such, an attempt was made to improve the resemblance of in vitro produced bone-like ECM to the physiological bone ECM.

ACKNOWLEDGMENTS
We gratefully acknowledge the financial support of the research program TTW with project number TTW 016.Vidi.188.021, which is (partly) financed by the Dutch Research Counsil (NWO).