Front. Bioeng. Biotechnol. Frontiers in Bioengineering and Biotechnology Front. Bioeng. Biotechnol. 2296-4185 Frontiers Media S.A. 1256825 10.3389/fbioe.2023.1256825 Bioengineering and Biotechnology Original Research Fluid flow shear stress and tissue remodeling—an orthodontic perspective: evidence synthesis and differential gene expression network analysis Nile et al. 10.3389/fbioe.2023.1256825 Nile Mustafa 1 Folwaczny Matthias 2 Wichelhaus Andrea 1 Baumert Uwe 1 * † Janjic Rankovic Mila 1 † 1 Department of Orthodontics and Dentofacial Orthopedics, LMU University Hospital, LMU Munich, Munich, Germany 2 Department of Conservative Dentistry and Periodontology, LMU University Hospital, LMU Munich, Munich, Germany

Edited by: Mario Rothbauer, Medical University of Vienna, Austria

Reviewed by: Marcel Rodrigues Ferreira, São Paulo State University, Brazil

Stefaan Verbruggen, Queen Mary University of London, United Kingdom

 *Correspondence: Uwe Baumert, uwe.baumert@med.uni-muenchen.de

These authors share senior authorship

 18 09 2023 2023 11 1256825 17 07 2023 28 08 2023 Copyright © 2023 Nile, Folwaczny, Wichelhaus, Baumert and Janjic Rankovic. 2023 Nile, Folwaczny, Wichelhaus, Baumert and Janjic Rankovic

This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

Introduction: This study aimed to identify and analyze in vitro studies investigating the biological effect of fluid-flow shear stress (FSS) on cells found in the periodontal ligament and bone tissue.

Method: We followed the PRISMA guideline for systematic reviews. A PubMed search strategy was developed, studies were selected according to predefined eligibility criteria, and the risk of bias was assessed. Relevant data related to cell source, applied FSS, and locus-specific expression were extracted. Based on this evidence synthesis and, as an original part of this work, analysis of differential gene expression using over-representation and network-analysis was performed. Five relevant publicly available gene expression datasets were analyzed using gene set enrichment analysis (GSEA).

Result: A total of 6,974 articles were identified. Titles and abstracts were screened, and 218 articles were selected for full-text assessment. Finally, 120 articles were included in this study. Sample size determination and statistical analysis related to methodological quality and the ethical statement item in reporting quality were most frequently identified as high risk of bias. The analyzed studies mostly used custom-made fluid-flow apparatuses (61.7%). FSS was most frequently applied for 0.5 h, 1 h, or 2 h, whereas FSS magnitudes ranged from 6 to 20 dyn/cm2 depending on cell type and flow profile. Fluid-flow frequencies of 1 Hz in human cells and 1 and 5 Hz in mouse cells were mostly applied. FSS upregulated genes/metabolites responsible for tissue formation (AKT1, alkaline phosphatase, BGLAP, BMP2, Ca2+, COL1A1, CTNNB1, GJA1, MAPK1/MAPK3, PDPN, RUNX2, SPP1, TNFRSF11B, VEGFA, WNT3A) and inflammation (nitric oxide, PGE-2, PGI-2, PTGS1, PTGS2). Protein-protein interaction networks were constructed and analyzed using over-representation analysis and GSEA to identify shared signaling pathways.

Conclusion: To our knowledge, this is the first review giving a comprehensive overview and discussion of methodological technical details regarding fluid flow application in 2D cell culture in vitro experimental conditions. Therefore, it is not only providing valuable information about cellular molecular events and their quantitative and qualitative analysis, but also confirming the reproducibility of previously published results.

 periodontal ligament cells osteoblasts osteocytes mesenchymal stem cells tissue remodeling fluid flow shear stress mechanobiology In vitro models section-at-acceptance Biomechanics 1 Introduction

The term “mechanotransduction” indicates processes through which cells sense mechanical stimuli and convert them into biological signals. Fluid flow shear stress (FSS) is a mechanobiological stimulus that is exerted by flowing fluids on cells within microporous tissues. Such stimulation, by virtue of altering the local interstitial and blood fluid flow, affects the behavior of cells by inducing biochemical and biophysical changes such as gene expression, activating signaling pathways and protein synthesis, which in turn would influence cell behavior such as proliferation, migration, differentiation, and apoptosis.

Interstitial fluid is considered a major component of the body mass and is commonly found in the extracellular matrix (Wittkowske et al., 2016). The term “interstitial fluid” is so-called due to plasma leakage from blood capillaries into the interstitial space, which is then returned to blood circulation by lymphatic drainage. Research has suggested that the main function of this fluid is to facilitate cell nutrition and waste removal (Cowin and Cardoso, 2015). Osmotic and hydrostatic pressure differences between blood capillaries, interstitial space, and lymphatics are the main driving forces behind this slow but constant interchange. Tensile stresses generated within the tissue matrix during function induce interstitial fluid flow changes, which elect signaling cascades by the surrounding mechanosensing cells (Tang et al., 2014; Jiang et al., 2016; Wittkowske et al., 2016; Li et al., 2019a).

Bone is a mechanosensing tissue and has a remarkable ability in sensing mechanical cues to maintain its structural integrity and functions. Mechanical stimulation perceived by skeletal cells, in particular interstitial FSS, guide bone remodeling, which is a physiological lifelong process. Another tissue known for its mechanosensing properties and closely related to bone remodeling is the periodontal ligament (PDL), a specialized tissue sandwiched between the tooth root and the alveolar tooth socket (Naveh et al., 2018). Earlier research on PDL mechanical properties highlighted the importance of circulatory pressure in maintaining the mechanical properties of the PDL (Bien, 1966). Recent research has suggested that fluid-filled PDL porosities and spaces give it a hydromechanical coupling property (Ashrafi et al., 2020), which is not only an important part of healthy masticatory function but may also be related to the tissue response during orthodontic treatment (Krishnan and Davidovitch, 2009). Orthodontic treatment is widely used in dentistry to correct dental malocclusion, improve function and dentofacial aesthetics. Moving a tooth using an orthodontic appliance is force- and time-dependent. When a constant force is applied to a tooth, several mechanical and biological changes occur within the periodontium around the tooth root. One of these mechanical changes is the disruption of interstitial fluid flow balance within and between the PDL and the bone. The disruption of fluid flow movement happens due to the stretching and compression of fluid-filled spaces and capillaries in the periodontal region, thus leading to differences in pore pressure (Ashrafi et al., 2020). Mathematical and in silico studies have suggested that porosity and permeability, properties of poro-elasticity, are influenced by pore pressure generated in tissues by mechanical loading (Weinbaum et al., 1994; Zeng et al., 1994; Chen et al., 1998; Ashrafi et al., 2020) (Figure 1A). Consequently, FSS is generated within porous tissues of the periodontium due to the flow of interstitial fluid along the pore pressure gradient (Wittkowske et al., 2016) (Figure 1B). The cellular changes that occur during orthodontic tooth movement (OTM) appear as bone resorption in direction of tooth movement at the compression side and as bone formation at the tension side, thus leading to tooth movement (Feller et al., 2015).

Generation of fluid flow in orthodontic tooth movement, its calculation, and modalities. (A) The effect of compression and stretching of a poro-elastic tissue on porosity size and pore pressure. (B) Illustration of fluid flow shear stress (τ) within a rectangular chamber (Q, flow rate; µ, viscosity of medium; height, h, and width, w, of the parallel flow chamber). Based on this equation, a change in the dimension of porosities influences the magnitude of FSS generated on cells. (C) Three common fluid flow patterns: steady laminar (D), pulsatile laminar (E), and oscillatory laminar flow (F). The fluid flow profile of each is determined by its flow rate and flow direction.

The PDL is a complex connective tissue composed of different cell types including but not limited to fibroblasts, stem cells, endothelial cells and macrophages (Marchesan et al., 2011). The presence of bone osteoblasts, osteocytes, and osteoclasts near the PDL makes these cells not less important in initiating the process of bone remodeling during OTM. It has been suggested that the loss of mechanical function is associated with the atrophy of PDL and extensive resorption of alveolar bone (Cohn, 1965), whereas excessive force leads to a loss of alveolar bone support (Nogueira et al., 2014). Such interrelated cellular adaptation to mechanical forces requires a biological system that senses these forces. There is increasing evidence that biological fluid flow mechanisms are involved in the adaptation of both soft and hard tissue (Buga et al., 1991; Klein-Nulend et al., 1996; Li et al., 2019a). The interlink between the PDL and bone through collagen fibers and blood supply allows the transmission of physical strain across the extracellular matrix of bone, thus resulting in a fluid movement within the bone tissue (McGarry et al., 2005; Feller et al., 2015; Wittkowske et al., 2016; Ashrafi et al., 2020). FSS generated by fluid movement in response to mechanical load is sensed by cells of the PDL and bone (van der Pauw et al., 2000; Wittkowske et al., 2016; Li et al., 2019a). Cell deformation caused by different mechanical stresses excites different cell signaling pathways responsible for cell remodeling due to cell deformation (Wittkowske et al., 2016). Anatomical wise, different cell signaling mechanisms are found in cells, at the level of the cell membrane, cell attachment and cell process. Moreover, mechanoreceptors such as Ruffini nerve endings may also be involved in periodontal tissue remodeling. Research has found that stretching and FSS activate different cell signaling pathways (McGarry et al., 2005): FSS is sensed first by cell membrane, while stretching is applied directly through cell attachment (Mullender et al., 2004).

To understand the biological effect of FSS on cell behavior, different in vitro setups have been utilized to determine the effect of shear stress on cells in a well-controlled cell culture environment. Since the magnitude of FSS is still debatable and may vary at different parts of the human body, it would be useful to study the different FSS profiles (steady laminar, pulsatile laminar and oscillatory laminar; Figures 1C–F) and magnitudes applied on similar cell types.

To get an overview on tissue remodeling, including the one during OTM, this article focused on the cellular behavior in response to FSS including cells from both bone and periodontal ligament tissue. Therefore, the aim of this review was to i) to provide overview of findings from 2D in vitro experiments utilizing FSS on bone and PDL cells, summarizing experimental settings (e.g., FSS parameters, in vitro apparatuses and methods used in determining FSS-related cell behavior) and gene and protein/metabolite expression data. ii) To elucidate the biological behavior of PDL and bone cells subjected to FSS, collected information was further used for enrichment analysis and protein-protein interaction (PPI) network construction and analysis. Additionally, the methodological and reporting quality of the included studies were assessed using an adapted risk of bias tool (Sun et al., 2021).

 2 Materials and methods

A systematic identification of relevant original research was performed based on the guidelines of PRISMA 2020 (Moher et al., 2009; Page et al., 2021). Since only in vitro studies were regarded, a registration of the study protocol in the PROSPERO database was not possible.

2.1 Eligibility criteria

The inclusion criteria were formulated based on the P.I.C.O. model (Schardt et al., 2007).

1. P(atient)/cells used: foremost target were primary human cell types that are known to be involved in orthopedics/orthodontics-related bone remodeling including PDL cells, osteoblasts, osteocytes, and mesenchymal stem cells. Additionally, primary mouse osteoblasts, osteocytes and PDL cells, and the mouse osteocyte-like cell line MLO-Y4 were considered. This cell line derived from long bones of a transgenic mouse expressing the immortalizing T antigen and resembles a mature osteocyte phenotype (Kato et al., 1997). Though this cell line is known not to present the full primary osteocyte phenotype (Zhang et al., 2019), a huge amount of evidence comes from studies using this cell line (Robling and Bonewald, 2020).

2. I(ntervetion): in vitro fluid-flow shear stress (FSS) applied to cells in 2D culture; no limitations concerning flow characteristics.

3. C(ontrol): Cells not subjected to FSS.

4. O(utcome): FSS parameters (apparatus, type of flow, magnitude, duration, and frequency); reporting expression of genes, proteins and/or metabolites related to FSS using quantitative methods including but not limited to quantitative and semi-quantitative reverse-transcriptase polymerase chain reaction (RT-qPCR, sqPCR), ELISA, Western blotting (WB), radioimmunoassay (RIA) and immunofluorescence (IF).

The exclusion criteria were adopted as follows.

1. Review article/short communication/protocol study.

2. Primary cells independent of origin, not fulfilling the inclusion criteria.

3. In vivo, in silico, ex vivo, or in situ models.

4. 3D cell culture models.

5. No FSS-related reporting on gene or protein expression or experiment.

6. Unclear/unrelated in vitro fluid-flow model design.

7. Article not published in English.

8. Unhealthy/sick patient or controls.

9. Co-culture studies unless monoculture data is provided.

2.2 Search plan and article selection

The search strategy aimed to identify studies that used in vitro cell models simulating bone and periodontal ligament tissues remodeling under exposure against FSS. Specific keywords were used to establish a search strategy, taking into consideration the field of study, mechanical stimulation applied, cells responsible for remodeling and the setting of the study (Table 1). The search was completed on 19.07.2021. All results were imported to EndNote® X9.3.1 (Clarivate Analytics, Philadelphia, Pennsylvania, United States).

Search strategy design used for the systematic review. The search strategy was built according to four main variables based on the topic and the scope of the project: field, force, cell type, setting. Different key words were used under each main variable to increase the detection capability.

Field Mechanical stimulation Cells Setting
orthodont* OR tooth movement OR ortho* OR Dent* OR Oral OR Mechanotransduction OR bone remodeling OR Stomatology OR mechanosens* OR calcium signaling OR osteogenesis OR orthodontic tooth movement OR osteoblastogenesis OR Periodont* Extracellular Fluid flow OR interstitial fluid flow OR "Shear Strength" [MeSH] OR Stress, "Mechanical" [MeSH] OR fluid flow OR pulsatile flow OR flow chamber OR shear stress OR FSS OR fluid shear stress OR shear strength OR laminar shear stress OR pulsatile fluid flow PDL OR hPDLCs OR hPDLFs OR hPDLs OR bone cells OR osteo* OR bone OR periodontal ligament OR fibroblast* OR MLO-Y4 OR human osteoblast-like cells OR human fetal osteoblast-like cells OR mesenchymal stem cells OR primary human osteoblasts OR progenitor cell* OR stem cell* OR human PDL-cells OR human PDL-fibroblasts OR human PDLFs OR human PDLs cell* OR model OR in vitro

Search: ("bone remodeling"[tw] OR "calcium signaling"[tw] OR "orthodontic tooth movement"[tw] OR "tooth movement"[tw] OR dent* OR mechanosens*[tw] OR mechanotransduction[tw] OR oral OR ortho* OR orthodont* OR osteoblastogenesis OR osteogenesis[tw] OR periodont* OR stomatology) AND ("Extracellular Fluid flow"[tw] OR "flow chamber"[tw] OR "fluid flow"[tw] OR "fluid shear stress"[tw] OR "interstitial fluid flow"[tw] OR "laminar shear stress"[tw] OR "pulsatile flow"[tw] OR "pulsatile fluid flow"[tw] OR "Shear Strength"[MeSH] OR "shear stress"[tw] OR "Stress, Mechanical"[MeSH] OR FSS[tw]) AND ("bone cells"[tw] OR "human fetal osteoblast-like cells"[tw] OR "human osteoblast-like cells"[tw] OR "human PDL-cells"[tw] OR "human PDL-fibroblasts"[tw] OR "human PDLFs"[tw] OR "human PDLs"[tw] OR "mesenchymal stem cells"[tw] OR "MLO-Y4"[tw] OR "periodontal ligament"[tw] OR "primary human osteoblasts"[tw] OR "progenitor cell*"[tw] OR "stem cell*"[tw] OR bone[tw] OR fibroblast*[tw] OR hPDLCs[tw] OR hPDLFs[tw] OR hPDLs[tw] OR osteo*[tw] OR PDL[tw]) AND ("in vitro"[tw] OR cell*[tw] OR model[tw])

Filters: from 1000/1/1 ‐ 2021/7/19 Sort by: Publication Date

The studies were initially filtered by reading the abstract and the title according to the established eligibility criteria. Then, full-text reading was done for the remaining studies. After full-text reading, the articles not fulfilling the inclusion criteria were excluded (Supplementary Data Sheet S1), while the remaining were considered for data extraction (Supplementary Data Sheet S2). All steps were discussed between the authors (M.N., M.J.R., U.B.) till reaching a unified decision.

 2.3 Risk of bias assessment

The risk of bias (RoB) of the included in vitro studies was evaluated using two separate assessment sheets specific for methodological or reporting RoB as previously published (Sun et al., 2021) (Supplementary Data Sheet S3). A score was applied to each criterion depending on its evaluation as previously published (Sun et al., 2021): “+” (low risk of bias), “–” (high risk of bias), “?” (incomplete/unclear risk of bias), and “n.a.” (not applicable). The data entry was carried out as previously reported (Sun et al., 2021), and for each cell group the RoB assessment was tabulated (Supplementary Data Sheet S3). All steps were discussed between the authors (M.N., M.J.R., U.B.) till reaching a unified decision.

 2.4 Data extraction

For each cell-type, data extraction tables were designed to collect data related to the experimental design and FSS-related expression outcome (Supplementary Data Sheet S2). Therein, we extracted study data as originally reported related to (i) cells used (hPDLCs, hMSCs, hOst, hOcyt, mPDLCs, mOst, mOcyt) including cell-related information (origin of cells, age, number and sex of donor, health status, tooth type, isolation method, passages used, and cell density/confluency used); (ii) fluid flow type (steady laminar, pulsatile laminar, oscillatory laminar); (iii) apparatus used in each experiment; (iv) shear stress magnitudes/frequency/duration; (v) genes and analytes (official gene symbol if applicable); (vi) the methods used to measure their expression (RT-qPCR, sqPCR, ELISA, WB, RIA, EMSA, IF), and (vii) the measured FSS-related differential expression. All steps were discussed between the authors (M.N., M.J.R., U.B.) till reaching a unified decision. The pattern of gene expression was determined as minimum/maximum expression and reported as fold change, relative expression or as a calculated ratio. Gene, protein and/or metabolite regulation was calculated and reported as ratio in relation to the corresponding control when applicable.

2.4.1 Information related to FSS

The type of fluid flow was defined as “steady laminar”, “pulsatile laminar” and “oscillatory laminar” based on the information given by the author of each study. The fluid flow character was confirmed by taking “flow rate” and “flow direction” into consideration, that was delivered by the motorized apparatuses to apply FSS to cells. The whole design of the fluid flow circuit such as, e.g., the utilization of the gravity force and the usage of axillary pulse dampers were taken into consideration. Characteristics of the fluid flow apparatuses and fluid flow circuits were extracted from the publication or references mentioned therein including manufacturer websites (Supplementary Data Sheet S2). Information was regarded as “not given” (n.g.) only after carefully following the previously mentioned steps. The fluid flow frequencies were extracted in Hertz (Hz) and FSS as dyn/cm2 or Pascal (Pa). Herein, FSS magnitudes are reported as dyn/cm2, which can be converted to Pascal according to the following relationship: 10 dyn/cm2 ≈ 1 Pa. FSS apparatuses identified in each study were grouped and quantified in tables based on cell type and fluid flow profiles (Supplementary Data Sheet S4). A quantitative summary of the most used fluid flow apparatuses and fluid flow profiles in each cell group was then provided (Tables 3, 4). Finally, the most frequently used FSS durations and magnitudes were quantified and summarized in tables using the relationship between dyn/cm2 and Pascal as mentioned (Supplementary Data Sheet S4).

 2.4.2 Information related to genes and proteins

All genes and molecules were extracted from the included studies without changing. The official name/symbol of the identified genes was determined using the HUGO Gene Nomenclature Committee (HGNC) database 1 for human genes, the Mouse Genome Informatics (MGI) 2 website for mouse genes, and the Gene Database at the National Center for Biotechnology Information 3 . Primer-BLAST 4 was used to verify the specificity of PCR primers. Specificity of protein data, antibodies or ELISA was verified using information provided in the publication and/or suppliers mentioned therein. If possible, official gene symbols according to HGNC or MGI were used. If antibody specificity was not sufficient, their targets according to the given manufacturer were recorded (Sun et al., 2021).

Gene and protein expression-related data were either directly acquired from the text source or extracted from the graphs using Engauge Digitizer Software (version 12.1) (Mitchell et al., 2019). To provide a comprehensive comparative approach of gene/protein expression data between different studies, relevant data and ratios were calculated if they were not mentioned by the study-authors. “Fold changes” were extracted if they were calculated according to Livak and Schmittgen (2001) or their calculation was reported as being “delta delta Ct” (ΔΔCt). Additionally, the following terms were introduced: “relative gene expression” was entitled to percentages or gene expression ratios normalized to the control, and not calculated by ΔΔCt; “ratio-calc” indicated manual calculation by dividing intervention/control; “ratio” indicated ratios given by study-authors such as normalization to control in case of small molecules data or in case of gene expression ratios, e.g., the ratio of RANKL/OPG or Bcl-2/Bax.

 2.4.3 Gene and small molecules analysis

Using the data from Supplementary Data Sheet S2, genes/metabolites were considered for further analysis, which were extracted from at least three (≥3) studies (Supplementary Data Sheet S5). This threshold aimed to accommodate variations in study numbers per cell group and ensure representative data. Studies reporting similar genes/metabolites were categorized by fluid flow profile and regulation direction (up or down). Upregulation included “increase” or “increase with plateau” patterns; downregulation encompassed “decrease” or “decrease with plateau”. All other patterns were labeled as “other” regulation changes. Repeated studies on the same gene/metabolite were counted as one. The fraction/percentage of upregulation of the same gene/metabolite was calculated as part of the total (Supplementary Data Sheet S6).

 2.5 Evidence-synthesis-based bioinformatical analysis of differential gene expression data 2.5.1 Over-representation and network analysis

Based on the data extracted (Supplementary Data Sheet S2), two species-specific lists with differential expressed genes (DEGs) were compiled according to the previously defined criteria (Sun et al., 2021). The “human” list included DEGs from reports studying the effect of FSS on human cells (hMSCs, hOst, hPDLCs), and the “mouse” list from studies using mouse cells (mOst, mOcyt). Both DEG lists were used to generate protein-protein interaction (PPI) networks and for over-representation analysis (ORA) as described (Sun et al., 2021) (Figure 2B) with modifications described below.

Workflows applied in this systematic review. (A) PRISMA2020 flow diagram for the whole process of study selection according to Moher et al. (2009) and Page et al. (2021). (B) Over representation and network analysis of gene lists derived from the review process. Gene lists were compiled, listing examined genes, proteins or metabolites showing differential gene expression after FSS application in primary human or mouse cell lines (Table 3, Supplementary Data Sheets S2, 6). From these, differential expressed gene (DEG) lists were generated according to the specified criteria. For both human and mouse DEG lists protein-protein interaction (PPI) networks were constructed using STRING-DB (Szklarczyk et al., 2023) (Figure 5). Pathway analysis was done using “GeneOntology/Biological Process” (Ashburner et al., 2000), KEGG (Kanehisa and Goto, 2000), and WikiPathways (Pico et al., 2008) (Supplementary Data Sheet S7). Subnetworks were identified in both PPI networks using the “GLay community detection algorithm” implemented in clusterMaker2 (Morris et al., 2011) (Table 6, Figure 5). Essential nodes, so-called “hub genes”, were detected using cytoHubba (Chin et al., 2014) (Figure 5, Supplementary Data Sheet S7). Altogether 5 high-throughput gene expression studies applying microarray or RNA-seq methodology were identified studying the effect of FSS on either human MSCs or the mouse osteocyte cell line MLO-Y4. Raw data were downloaded from the respective repositories and individually processed using bioinformatics workflows. Pre-ranked gene lists were generated, and gene set enrichment analysis (GSEA) was applied to each using easyGSEA (Cheng et al., 2021). Cell-type specific comparisons were done using easyVizR (Cheng et al., 2021) (Figure 6, Supplementary Data Sheet S8).

Potential interactions between DEGs at protein level were analyzed using protein-protein interaction (PPI) networks, which were generated with the “Search Tool for the Retrieval of Interacting Genes/Proteins” database (STRING-DB) (v11.5)5 (Szklarczyk et al., 2019; Szklarczyk et al., 2023) using the stringApp plugin version 2.0.1 (Doncheva et al., 2019) with Cytoscape version 3.9.1 (Shannon et al., 2003; Su et al., 2014). Only high confidence interactions were included in the predicted networks, i.e., a minimum required combined score of 0.7 was applied. Cluster analysis and cluster visualization of both human and mouse PPI networks was done by applying the “GLay Community Detection” algorithm (Su et al., 2010) as implemented in the clusterMaker2 app version 2.3.4 (Morris et al., 2011) with default settings. Hub genes, i.e., essential genes in a network, were identified with cytoHubba version 0.1 (Chin et al., 2014) using default settings. For each node, a total score was calculated based on the eleven different local and global topological methods calculated by this plugin (Sun et al., 2021) (Supplementary Data Sheet S7). A node was considered a hub node if its total score was at least twofold higher than the mean total score of all nodes of that network (Sun et al., 2021). Networks were visualized with Cytoscape version 3.9.1.

Over-representation analysis was applied to the network with the stringApp plugin using the “GeneOntology/Biological Process” (Ashburner et al., 2000), WikiPathways (Pico et al., 2008; Martens et al., 2021), and the “Kyoto Encyclopedia of Genes and Genomes” (KEGG) (Kanehisa and Goto, 2000; Du et al., 2014) databases. As previously described, a cut-off ratio of ≥0.05 (i.e., at least 5% of the queried genes must be included in a specific database entry) was applied to increase the specificity of the enrichment terms (Sun et al., 2021). If available, the 15 most significant terms or pathways were reported. General pathways related to cancer and/or infectious diseases like “Chagas disease” (hsa05142) or “Pertussis” (hsa05133) were removed, and their accession numbers were reported. Gene set enrichment was applied to the complete networks, and to each identified cluster (Supplementary Data Sheet S7).

 2.5.2 High-throughput gene expression studies and pre-ranked gene set enrichment analysis (GSEA)

Though data extraction identified several high-throughput gene expression studies, comparable gene expression raw data (microarray or RNA-seq) from three of the identified studies only (Govey et al., 2014; Govey et al., 2015; Li et al., 2019b) was publicly available (Table 2). Data descriptions and the related publications were evaluated concerning experimental conditions (cell type characteristics, FSS parameters, controls) and the number of samples for each (FSS or corresponding control). Only studies and experimental conditions represented by at least two replicates each were considered. To identify additional data specific to our preselected cell types and in relation to FSS, hand-searching of Gene Expression Omnibus (GEO) 6 , ArrayExpress 7 , and the Sequence Read Archive 8 was done. Two additional studies were identified: GSE18114 from GEO and E-MEXP-1923 from ArrayExpress (Glossop and Cartmell, 2009). Since no publication was identified related to the GSE18114 data, and the ArrayExpress dataset passed selection criteria, our reanalysis was finally based on five high-throughput gene expression studies (Table 2). In two microarray studies (GSE82269, E-MEXP-1923) FSS was applied to hMSCs (Glossop and Cartmell, 2009; Diaz et al., 2017), and in three studies FSS was applied to the mouse osteocytic-cell line MLO-Y4. The latter consists of one microarray study (GSE42874) (Govey et al., 2014) and two RNA-seq studies (SRP060567/GSE70667 and SRP212008) (Govey et al., 2015; Li et al., 2019b) (Figure 2B). To ease reanalysis, uniformly processed RNA-seq gene data from both studies were downloaded from the recount3 repository 9 (Wilks et al., 2021). The detailed workflow is outlined in Supplementary Data Sheet S8. Shortly, the raw data were downloaded from the specific repositories and analyzed separately with workflows specific to the platform used (microarray, RNA-seq) (Figure 2B). Pre-ranked gene lists were exported, and gene set enrichment analysis (GSEA) was applied using easyGSEA 10 (Cheng et al., 2021) using the provided gene set definitions from KEGG, GO/Biological Process, and WikiPathways. For each cell type (hMSCs, mouse MLO-Y4), shared features were identified separately using easyVizR using default filtering settings (p < 0.05) (Figure 2B) (Cheng et al., 2021). All comparison results with FDR <0.1 in either of the easyGSEA results were retained. No further filtering was applied (Supplementary Data Sheet S8). As mentioned above, general pathways related to cancer and/or infectious diseases like “C. disease” (hsa05142) or “Pertussis” (hsa05133) were removed, and their accession numbers reported.

Summary of the high-throughput gene expression studies analyzed with pre-ranked GSEA.

Study Id References Platform Type Cell type FSS parameters Samples
E-MEXP-1923 (AE) Glossop and Cartmell (2009) Affymetrix GeneChip Human Genome U133 Plus 2.0 [HG-U133_Plus_2] Microarray (probe) Human bone marrow MSCs from 1 female donor (23 years old; commercial source) FSS: 1 dyn/cm2 for 1 h and 2 h post-FSS; Controls: sham treated 2 samples FSS, 2 samples control
GSE82269 (GEO) Diaz et al. (2017) Illumina HumanHT-12 V4.0 Expression BeadChip (GPL10558) Microarray (BeadArray) Human bone marrow MSCs from 3 donors (1 female, 20years; 2 male, 22 and 25 years; commercial source) FSS: 15 dyn/cm2 for 6 h; Controls: sham treated 3 samples FSS, 3 samples control (batch correction: donor)
GSE42874 (GEO) Govey et al. (2014) Affymetrix Mouse Genome 430A 2.0 Array (GPL8321) Microarray (probe) Mouse MLO-Y4 osteocyte-like cell line FSS: sinusoidally oscillating FF w/peak shear stress of 10 dyn/cm2 @ 1 Hz for 2 h and post-FSS incubation for 0/2/8/24 h; Controls: paired sham treated 3 samples FSS (2 h FSS and 2 h post-FSS), 3 samples control
GSE70667 (recount3: SRP060567) Govey et al. (2015) Illumina HiSeq 2500 (Mus musculus) (GPL17021) RNA-seq Mouse MLO-Y4 osteocyte-like cell line FSS: sinusoidally oscillating FF w/peak shear stress of 10 dyn/cm2 @ 1 Hz for 2 h and post-FSS incubation for 2 h; Controls: paired sham treated 3 samples FSS, 3 samples control
SRP212008 (recount3) Li et al. (2019b) Illumina NextSeq 550 RNA-seq Mouse MLO-Y4 osteocyte-like cell line FSS: 15 dyn/cm2 oscillatory FSS @ 1 Hz for 2 h; Controls: sham treated 3 samples FSS, 3 samples control

GEO, Gene Expression Omnibus (https://www.ncbi.nlm.nih.gov/geo/; AE, ArrayExpress (https://www.ebi.ac.uk/biostudies/arrayexpress/); recount3 study explorer (https://jhubiostatistics.shinyapps.io/recount3-study-explorer/).

 3 Results 3.1 Study selection

The study selection process was summarized in the PRISMA 2020 flow chart (Page et al., 2021) (Figure 2A). The finalized search strategy identified 6,974 unique studies. No articles were added using manual search of specific journals or reference chaining. After duplicate removal, title/abstract of all remaining records were screened, and 6,756 articles were excluded according to the predefined criteria. The remaining 218 reports were assessed for eligibility by full-text reading, and 98 articles were excluded in obedience to the pre-defined exclusion criteria (Figure 2A; Supplementary Data Sheet S1). Finally, 120 studies were included in this systematic review for data extraction and subsequent differential gene expression network analysis. No studies were identified related to mouse PDL cells.

 3.2 Risk of bias assessment

The assessment of methodological quality was based on 15 predefined criteria (Supplementary Data Sheet S3). Overall, a large variability was observed between cell groups (Figure 3). Generally, a high risk of bias score was mostly assigned for “Statistical analysis” and “Sample size determination”. For both, primary mouse osteoblasts (mOst) and mouse osteocytes (mOcyt) a high risk of bias ratio in “Reporting bias: selective outcome data” was identified. Publications reporting on human primary osteoblasts (hOst) cell group showed an increased risk of bias ratio for the criterion “Statement conflict of interest/funding” compared to other cell groups. “Incomplete/unclear risk of bias” figure was most frequently identified for “Test organism/system”, “Test substance/treatment details”, “confounding bias” and “selection bias: allocation concealment”.

Summaries of the Risk of Bias (RoB) assessment of the included studies. In the first row, the RoB of the studies on human primary cells (MSCs, osteoblasts, and PDLCs), and in the second row the studies on primary mouse cells (osteoblasts, osteocytes) are shown. For each cell type, both RoB categories (reporting and methodological RoB) are given. RoB scoring was applied as follows: “+”, low risk of bias (LoB); “-”“, high risk of bias (HoB); “?”, not applicable or incomplete/unclear risk of bias.

Concerning the reporting quality assessment (Figure 3), which was based on 10 predefined criteria (Supplementary Data Sheet S3), a high risk of bias figure was mostly assigned for the criterion “Ethical statement”, while incomplete/unclear risk of bias was assigned for “Cell maintenance conditions” in most of the studies. Due to the differences in the number of studies included in each cell group, a groupwise comparison was considered not conclusive.

 3.3 Fluid flow characteristics

Fluid flow characteristics were defined based on fluid flow direction and rate as well as the fluid flow systems used to create the fluid movement in each study (Figure 1C, Table 3). Oscillatory laminar fluid flow was the most investigated fluid flow type in mOcyt- (21/58; 32.6%) and hMSCs-related studies (10/21; 47.6%). Both, steady laminar and pulsatile laminar fluid flow characters were equally investigated in mOst-related studies (9/22; 40.9%). Nevertheless, pulsatile laminar fluid flow was mainly used in hOst cell group (10/14; 71.4%). Steady laminar fluid flow was predominantly used in human periodontal ligament (hPDLCs) cell group (5/7; 71.4%), while pulsatile laminar fluid flow was used by the only identified study of human osteocytes (hOcyt).

A numerical summary of fluid flow pattern in relation to each cell group type (Figures 1C–F; Supplementary Data Sheet S4).

Cell type (abbreviation) No. of studies Flow profiles No. of flow profiles used
Steady laminar Pulsatile laminar Oscillatory laminar Unknown
Human MSCs (hMSCs) 21 8 2 10 1 21
Human osteoblasts (hOst) 14 2 10 2 0 14
Human osteocytes (hOcyt) 1 0 1 0 0 1
Human PDL (hPDLCs) 7 5 2 0 0 7
Mouse osteoblasts (mOst) 22 9 9 4 0 22
Mouse osteocytes (mOcyt) 57 15 19 21 3 58
Total 122 39 43 37 4 123
3.4 Fluid flow apparatus for shear stress application

Fluid flow systems were analyzed for each cell group by taking into consideration the whole fluid-flow experimental setting including but not limited to the type of motor used, accessory fluid flow pulsation damper and the utilization of gravity forces. Afterwards, the FSS systems were compared with the fluid flow description by each author (Table 4).

A numerical summary of fluid flow chamber type in relation to cell group type (for full details refer to Supplementary Data Sheet S4).

Chamber type Mouse Human Total Fluid flow type
Osteocytes Osteoblasts MSCs Osteoblasts PDLCs Osteocytes
Custom-made 35 11 12 10 7 1 76 Steady laminar, oscillatory laminar, pulsatile laminar
Rocking culture system 1 1 1 1 0 0 4 Oscillatory laminar
Streamer, Streamer Gold, FlexFlow (Flexcell Inc.) 8 0 4 0 0 0 12 Steady laminar, oscillatory laminar, pulsatile laminar
IBIDI 2 2 2 0 0 0 6 Steady laminar, oscillatory laminar, pulsatile laminar
Rotational orbital shaker 1 0 1 1 0 0 3 Oscillatory laminar
Cytodyne a 2 4 0 2 0 0 8 Steady laminar, pulsatile laminar
PeCon parallel plate (PeCon GmbH) a 1 0 0 0 0 0 1 Oscillatory laminar
Focht Chamber System 2 (Bioptechs Inc.) a 1 0 0 0 0 0 1 Pulsatile
Bioflux system (Fluxion Biosciences) 0 1 0 0 0 0 1 Steady laminar
Not given 7 3 1 0 0 0 11 Laminar, oscillatory laminar, pulsatile laminar
Total 58 22 21 14 7 1 123

No recent information available.

Different fluid flow systems were used to deliver steady laminar, pulsatile laminar and oscillatory laminar FSS on cells. Most of the identified fluid flow systems were custom-made (76/123; 61.7%) and were used to deliver either of the three fluid flow types. To deliver all previously mentioned fluid flow types, commercial systems from two companies were identified: Flexcell Inc (12/123; 9.7%) and ibidi (6/123; 4.8%). Flexcell Inc (Flexcell Corporation, Hillsborough, NC, United States) provides two fluid flow systems called FlexFlow™ and FlexCell Streamer®. FlexFlow™ consists of a single chamber and can be used either for tension or FSS applications, whereas the FlexCell Streamer® is used to apply FSS to up to 6 culture slides in parallel. Both fluid flow variants provide FSS from 0–35 dyn/cm2. Ibidi GmbH (Gräfelfing, DE) provides a fluid flow system consisting of a specific pump system, a fluidic unit and channel slides called “µ-Slide”. In some identified studies, an ibidi channel slide was integrated into a custom-made fluid flow circuit driven by pumps from various manufacturers (e.g., Fu et al., 2008; Seref-Ferlengez et al., 2016; Lee et al., 2017). A “Cytodyne flow chamber” was used in 8/123 (6.5%) of the identified studies. We were not able to find any current information related to the manufacturer or this product.

Rocking culture systems (4/123; 3.2%) and rotational orbital shakers (3/123; 2.4%) were used to create an oscillatory fluid flow culturing environment. The PeCon parallel plate system (PeCon GmbH, Erbach, DE), the Focht Chamber System 2 (Bioptechs Inc., Butler, PA, United States) and the Bioflux system (Fluxion Biosciences, Oakland, CA, United States) are fluid flow chambers used for live imaging and were identified once each.

 3.5 Fluid flow shear stress

A summary of the most frequently used FSS magnitudes and durations in hMSCs, hOst, hPDLCs, mOst and mOcyt is provided in Supplementary Data Sheet S4. Herein, we provide a short overview of the main findings. If several FSS durations were listed, the most frequently reported one is mentioned first and the following in descending order.

3.5.1 Human cells

In the hMSCs cell group, the most frequently used FSS durations with a steady laminar fluid flow profile were 0.5, 6 and 24 h and the most frequently used FSS magnitude was 2 dyn/cm2 for 2 h. Nevertheless, the mostly used FSS duration in the hPDLCs cell group was 2 h and the mostly used FSS magnitude was 6 dyn/cm2 for 2 h or 4 h. For oscillatory fluid flow, the mostly applied FSS durations in the hMSCs cell group, were 2, 1 and 0.5 h and the mostly used FSS magnitudes were 20 dyn/cm2 for 0.05 h, 5 dyn/cm2 for 0.05 h, 10 dyn/cm2 for 0.05 and 2 h. For pulsatile laminar fluid flow, the most frequently used FSS duration in the hMSCs cell group was 1 h, and the most frequently used FSS magnitude was 6 dyn/cm2 for 1 h. Similarly, the mostly used FSS duration in the hOst cell group was 1 h, while the most frequently applied FSS magnitude/duration combination was 7 dyn/cm2 for 1 h.

 3.5.2 Mouse cells

In the mOcyt cell group, the most frequently used FSS durations using a steady laminar fluid flow profile were 2 and 0.5 h. Magnitudes of 16 dyn/cm2 (for 2, 0.5, 4, 1 and 24 h) and 10 dyn/cm2 (for 0.16 h and 2 h) were the ones most commonly applied. In contrast, the mostly used FSS durations for the same fluid flow profile in the mOst cell group were 0.5, 1, 0.16 h and the mostly used FSS magnitude was 12 dyn/cm2 for 1 h. For oscillatory fluid flow, the mostly used FSS durations in the mOcyt cell group were 1 and 2 h and the mostly used FSS magnitudes were 10 dyn/cm2 for 2 and 1 h and 20 dyn/cm2 for 1 and 0.15 h. For the same fluid flow profile, the mostly used FSS durations in the mOst cell group were 0.5 and 1 h. For pulsatile laminar fluid flow, the most frequently used FSS durations in the mOcyt cell group, were 1 and 0.5 h and the mostly used FSS magnitudes were 7 dyn/cm2 for 1 h, 10 dyn/cm2 for 0.5, 0.16 and 1 h and 4 dyn/cm2 for 0.5, 1 and 2 h. Similarly, the most frequently used FSS durations listed, were 1, 0.5 h in the mOst cell group. However, the mostly used FSS magnitude was 6 dyn/cm2 for 0.083, 0.5 and 1 h.

 3.6 Fluid flow frequencies

The most frequently used fluid flow frequencies were 1 Hz in the hMSCs cell group and 5 Hz in hOst, hOcyt and hPDLCs cell groups. A frequency of 5 Hz was reportedly used two times more than 1 Hz in the mOst cell group, while equally used in the mOcyt cell group (Supplementary Data Sheet S4).

 3.7 Analysis of gene, protein, or metabolite expression

Though species-specific data was extracted and analyzed, in the following sections gene symbols were written according to human gene nomenclature independent of the species originally reported in the relevant study (i.e., all uppercase). Information on the experimental conditions related to cell culture and FSS application and the differential expression of genes, proteins, or metabolites after FSS application was summarized in Supplementary Data Sheet S2. Due to the huge body of information, the focus in this section was placed on the most frequently studied genes, proteins, or metabolites (Supplementary Data Sheets S5, 6; Figure 4).

Proportion of studies reporting an upregulation of genes or metabolites in human and mouse cell groups. Altogether 22 genes or metabolites were identified, which were investigated by at least three studies in the mentioned cell groups (hMSCs, human MSCs; hOst, human osteoblasts; hPDLCs, human PDL cells; mOst, mouse osteoblasts; mOcyt, mouse osteocytes). Bubble size and its color represents the number of studies (● 3, ● 4–7, ● 8–11 or ● ≥ 12). The fraction of studies reporting an upregulation is given for each gene or metabolite. It was calculated as total number of studies reporting an upregulation divided by the total number of studies reporting the specific gene or metabolite (Supplementary Data Sheet S6).

Expression of genes, proteins, or metabolites related to FSS application was extracted from the 120 included studies, of which 43 studies reported on human and 79 studies reported on mouse cell types, including two studies that investigated both species (Figure 2B). The human list contained 108 and the mouse list contained 111 (Table 5). Genes or proteins that were not clearly assigned to a specific gene symbol due to ambiguities in the reported PCR primers or antibodies used in WB or ELISA assays were included (e.g., “MAPK3/MAPK1”, “GSK3A/GSK3B”) (Table 5). Additionally, information on posttranslational modifications like phosphorylation or protein translocation between cytoplasm and nucleus were extracted (Supplementary Data Sheet S2), but not included in the following considerations.

Genes, proteins, and metabolites analyzed in the included studies in primary human or mouse cell lines (Supplementary Data Sheet S2). Differential expressed genes (DEGs) used for the subsequent over-representation analysis (ORA) are given in bold (Figures 2B, 5; Supplementary Data Sheet S7).

Human Mouse
ADIPOQ, AKT1, ALPL, BAX, BCL2, BGLAP, BMP2, BMP7, Calcium, CCL2 (?), CCL3, CCL5, CD44, COL1A1, CTNNB1, CXCL8, CYP24A1, CYP27B1, DUSP6, EGF, EGR1, ENG, FGF2, FOS, FOSB, FOSL1, FOSL2, GADD45B, GLI1, GSK3B, HGF, HIF1A, HMOX1, IBSP, IER3, IGF1, IGF2, IGFBP1, IL1B, IL1RN, IL6, ITGB1, JUN, KDR, MAP3K8, MAPK14/MAPK11/MAPK12, MAPK3/MAPK1, MAPK8, MMP1, MMP2, MRTFA, MTOR, MYH11, MYH2, NANOG, NFATC1, NFKB1, Nitric oxide, NOS1, NOS2, NOS3, PDGFA, PDGFB, PDGFRA, PDGFRB, PECAM1, PGE-2, PGF2α, PGI-2, PIK3CB, PON1, POU5F1, PPARA, PPP1R12A, PTCH1, PTGES, PTGS1, PTGS2, PTK2, ratio (BCL2/BAX), ratio (RANKL/OPG), RCAN1, RHOA, RPS6KB1, RUNX2, S100A4, S100A7, S100A8, SOX2, SP7, SPP1, SRC, TAGLN, TGFB1, TGFB3, THY1, TIMP1, TIMP2, TNAP, TNFAIP6, TNFRSF11B, TNFSF11, TRPV1, TRPV4, VDR, VEGFA, VWF, WIF1 6-keto-PGF1α, ACTB, AKT1, ALPL, APC, AXIN1, BAD, BAX, BCL2, BGLAP, BMPR1A, Calcium, cAMP, CASP3, CASP3/CASP7, CCL2, CCL3, CCL4, CCL5, CCL7, CCN1, CCND1, CD44, CSF2, CTNNB1, CXCL1, CXCL2, cyclic AMP, DKK1, DKK2, DLX1, DLX5, DMP1, EFNA2, EFNB2, EPHA2, EPHB4, ERAL1, ESR1, FABP4, FLT1, FOS, FZD6, GJA1, GJC1, GSK3A, GSK3A/GSK3B, GSK3B, HGF, IGF1, IGF1R, IL11, IL6, ITGAV, ITGB1, ITGB3, JAK2, JUN, KDR, LEF1, LEPR, LIF, LIMK2, LPL, LRP5, MAPK3/MAPK1, MEPE, MKI67, MMP14, NCOA1, Nitric oxide, Nos (unspecific), NOS2, NPY, NRP1, PDPN, PGE-2, PGF2α, PHEX, PIEZO1, PIEZO2, PKD1, PKD2, PPARG, PTGER2, PTGS1, PTGS2, PTHLH, ratio (BCL2/BAX), ratio (RANKL/OPG), RUNX2, SFRP1, SFRP2, SFRP4, SOST, SP7, SPP1, STAT3, TGFB1, TJP1, TNF, TNFRSF11B, TNFSF11, TUBA1A, VCL, VEGFA, WNT1, WNT3A, WNT4, WNT5A, ZIC1

A total of 22 genes or proteins and metabolites were investigated by at least three studies in relation to FSS in the five included cell groups (Figure 4, Supplementary Data Sheet S6). These include (gene/metabolite; cell type): alkaline phosphatase (ALP/ALPL/ALPP; hMSCs, hPDLCs); BGLAP (osteocalcin; hMSCs, mOcyt); BMP2 (bone morphogenetic protein 2; hMSCs); calcium (Ca2+; hMSCs, mOcyt); COL1A1 (collagen 1A1; hMSCs); MAPK1/MAPK3 (mitogen-activated protein kinase 1/3, aka ERK1/2; hMSCs, mOst, mOcyt); NO (nitric oxide; hMSCs, hOst, mOst, mOcyt); PGE-2 (prostaglandin E2; hOst, mOst, mOcyt); PGI-2 (prostacyclin/prostaglandin I2; hOst); PTGS1 (prostaglandin-endoperoxide synthase 1; hOst); PTGS2 (prostaglandin-endoperoxide synthase 2; hMSCs, hOst, mOst, mOcyt); RUNX2 (Runt-related transcription factor 2; hMSCs, mOcyt); SPP1 (secreted phosphoprotein 1; hMSCs, hOst, mOst, mOcyt); VEGFA (vascular endothelial growth factor A; hMSCs, mOcyt). Additionally, in mouse osteocytes the following loci were analyzed: AKT1 (AKT serine/threonine kinase 1), CTNNB1 (β1-catenin), GJA1 (gap junction protein α1), PDPN (podoplanin), SOST (sclerostin), TNFRSF11B (TNF receptor superfamily member 11b), TNFSF11 (tumor necrosis factor superfamily member 11), and WNT3A (Wnt family member 3A).

3.7.1 Human cells

Four to seven studies reported an FSS-related increase of Ca2+, VEGFA and MAPK1/MAPK3 in the hMSCs cell group, and nitric oxide and PGE-2 in the hOst cell group. A complete upregulation has been reported of nitric oxide and COL1A1 in the hMSCs cell group, PGI-2, PTGS1 and PTGS2 in the hOst cell group and ALP/ALPL/ALPP in the hPDLCs cell group by three studies each. While the reported upregulation percentages were 83% for ALP, 86% for PTGS2, 67% for BMP2 and 67% for SPP1 in the hMSCs cell group, BGLAP and RUNX2 were the least upregulated genes in response to FSS at 57% and 60%, respectively.

 3.7.2 Mouse cells

At least twelve studies reported prevalent upregulation of Ca2+, PGE-2 and PTGS2 in the mOcyt cell group in response to FSS. More than three studies reported dominant upregulation in nitric oxide, PTGS2 and PGE-2 in the mOst cell group. FSS induced upregulation of nitric oxide, BGLAP, VEGFA, RUNX2, SPP1, MAPK1/MAPK3, AKT1, PDPN, CTNNB1, GJA1, TNFRSF11B and WNT3A by at least 65% in the mOcyt cell group. Nevertheless, TNFSF11 and SOST were the least expressed genes in the mOcyt cell group at 64% and 25%, respectively.

 3.8 Over-representation and network analysis

Over-representation analysis (ORA) and consecutive network analysis were done as described (Figure 2B). ORA was applied to gene lists compiled from studies that applied FSS to human (hMSCs, hOst, hPDLCs) and mouse (mOst, mOcyt) cells. Altogether, 80 of 108 genes (∼74.1%) from the human and 76 of 111 genes (∼68.5%) from the mouse list (Table 5) were identified as differentially expressed genes (DEGs) as defined above.

3.8.1 PPI networks

From both gene lists protein-protein interaction (PPI) networks were generated. The network statistics for both networks were calculated (Figure 5). The human network consisted of 73/80 DEGs; 7 DEGs were placed outside the PPI network (Table 6; Figure 5). The mouse PPI network included 74/76 DEGs; 2 DEGs were placed outside (Table 6; Figure 5). Most of the human DEGs were reported to be upregulated after FSS application in most of the cell types (e.g., VEGFA, PTGS2, ALPL, RUNX2, FGF2, TIMP1). Only five genes (PUSF1, PDGFRA, PDGFRB, PTCH1, GLI1) were reported to be downregulated. SPP1, BGLAP and IL1B were the only three extracted DEGs, that showed a cell-type specific, contrasting regulation pattern of either up- or downregulation. Overall, a similar pattern was observed in the network related to mouse DEGs. Most of the genetic loci were upregulated in either of both cell types; PTGS2, WNT1, WNT3A, and TNFSF11 were upregulated in both cell types. Contrasting regulation was found for LRP5 (mOst: upregulation; mOcyt: downregulation), whereas only DKK1 was downregulated in both mouse cell types after FSS application. Other genes (DKK2, FZD6, SFRP1, SFRP2, DLX5, HGF, TNF, TJP1, LPL, FABP4, BAD, BAX) were downregulated in either mOst or mOcyt.

Human (left column) and mouse (right column) protein-protein interaction (PPI) networks were created using lists of DEGs (Table 5). (A, B) For each species PPI networks and their basic network statistics were created using the stringApp and Cytoscape. (C, D) Identified sub-networks identified from the complete networks shown above (A, B) (Table 6, Supplementary Data Sheet S7). The legend applies to all four networks: the node size corresponds to the number of reports identified for the corresponding underlying gene. The edges’ line style depicts the number of STRING sources for the given connection. Additionally, hub genes identified with cytoHubba (Chin et al., 2014) were colored red (Table 6). Pathway analysis was applied to each of the complete networks (A, B) and to each of the subnetworks (C, D) (Supplementary Data Sheet S7).

Differential expressed genes (DEGs) from the human and mouse gene lists, their contributions to the respective protein-protein interaction networks and the “GLay community” clusters. Identified hub genes were underlined. The corresponding networks and clusters were depicted in Figure 5.

Gene list Genes (n) in STRING network Cluster number Number of genes (n) and genes in cluster Cluster score
H.s Genes in the H.s. PPI network (73): ADIPOQ, ALPL, BAX, BCL2, BGLAP, BMP2, BMP7, CCL3, CCL5, CD44, COL1A1, CXCL8, CYP24A1, CYP27B1, EGR1, ENG, FGF2, FOS, FOSB, FOSL1, FOSL2, GLI1, HIF1A, HMOX1, IBSP, IGF1, IGFBP1, IL1B, IL1RN, IL6, ITGB1, JUN, KDR, MAP3K8, MMP1, MMP2, MRTFA, MYH11, MYH2, NANOG, NFKB1, NOS1, NOS2, NOS3, PDGFA, PDGFRA, PDGFRB, PECAM1, POU5F1, PPARA, PTCH1, PTGES, PTGS1, PTGS2, RHOA, RUNX2, S100A4, S100A7, S100A8, SOX2, SP7, SPP1, TAGLN, TGFB1, THY1, TIMP1, TIMP2, TNFAIP6, TNFRSF11B, TNFSF11, VDR, VEGFA, VWF Genes outside the H.s. PPI network (7): DUSP6, GADD45B, IER3, PON1, TRPV1, TRPV4, WIF1 Hub genes (cut-off: SUMscore ≥4596.6; n = 11): VEGFA, IL6, FGF2, IL1B, MMP2, TGFB1, MMP1, IGF1, TIMP1, RUNX2, JUN 1 (20) ALPL, BAX, BCL2, BGLAP, BMP2, BMP7, COL1A1, CYP24A1, CYP27B1, GLI1, IBSP, IGF1, IGFBP1, PTCH1, RUNX2, SP7, SPP1, TNFRSF11B, TNFSF11, VDR 0.191
2 (26) ADIPOQ, CCL3, CCL5, CXCL8, HMOX1, IL1B, IL1RN, IL6, JUN, MAP3K8, MMP1, MMP2, NFKB1, NOS1, NOS2, NOS3, PPARA, PTGES, PTGS1, PTGS2, S100A7, TGFB1, TIMP1, TIMP2, TNFAIP6, VEGFA 0.239
3 (22) CD44, ENG, FGF2, HIF1A, ITGB1, KDR, MRTFA, MYH11, MYH2, NANOG, PDGFA, PDGFRA, PDGFRB, PECAM1, POU5F1, RHOA, S100A4, S100A8, SOX2, TAGLN, THY1, VWF 0.211
4 (5) EGR1, FOS, FOSB, FOSL1, FOSL2 0.044
M.m Genes in the M.m. PPI network (74): ACTB, ALPL, APC, AXIN1, BAD, BAX, BGLAP, BMPR1A, CASP3, CCL7, CCNA2, CCND1, CD44, CTNNB1, CXCL1, CXCL2, DKK1, DKK2, DLX1, DLX5, DMP1, ESR1, FABP4, FOS, FZD6, GJA1, GJC1, HGF, IGF1, IGF1R, IL11, IL6, ITGAV, ITGB3, JUN, LEF1, LEPR, LIF, LPL, LRP5, MEPE, MKI67, MMP14, NCOA1, NOS2, NPY, PDPN, PHEX, PIEZO1, PKD1, PKD2, PPARG, PTGER2, PTGS1, PTGS2, RUNX2, SFRP1, SFRP2, SFRP4, SOST, SP7, SPP1, TGFB1, TJP1, TNF, TNFRSF11B, TNFSF11, TUBA1A, VCL, VEGFA, WNT1, WNT3A, WNT4, WNT5A Genes outside the M.m. PPI network (2): ERAL1, PIEZO2 Hub genes (cut-off: SUMscore ≥12,622.9; n = 12): CTNNB1, WNT1, WNT3A, WNT5A, DKK1, WNT4, LRP5, SFRP1, IL6, TNF, ACTB, JUN 1 (21) APC, Axin1, BAD, BAX, CASP3, CCN1, CCND1, DKK1, DKK2, ESR1, FZD6, LEF1, LRP5, MKI67, SFRP1, SFRP2, SFRP4, WNT1, WNT3A, WNT4, WNT5A 0.223
2 (33) ACTB, CCL7, CD44, CXCL1, CXCL2, FABP4, FOS, HGF, IGF1, IGF1R, IL11, IL6, ITGAV, ITGB3, JUN, LEPR, LIF, LPL, MMP14, NCOA1, NOS2, NPY, PDPN, PPARG, PTGER2, PTGS1, PTGS2, SPP1, TGFB1, TNF, TNFRSF11B, TNFSF11, VEGFA 0.25
3 (10) BMPR1A, CTNNB1, GJA1, GJC1, PIEZO1, PKD1, PKD2, TJP1, TUBA1A, VCL 0.06
4 (10) ALPL, BGLAP, DLX1, DLX5, DMP1, MEPE, PHEX, RUNX2, SOST, SP7 0.088

Eleven hub genes were identified in the human network, and twelve in the mouse network (Figure 5; Table 6, Supplementary Data Sheet S7). In each network, four highly connected gene clusters were identified (Figure 5; Table 6, Supplementary Data Sheet S7). The smallest one (cluster #4 from the human network) contained 5 genes from the human DEG list. Cluster #2 identified in the mouse network was the largest one, compromising 33 genes from the mouse DEG list.

 3.8.2 Over-representation analysis (ORA)

ORA was applied to identify biological processes and signaling pathways involved in regulation of gene expression after FSS application as described. For each database, the top 15 enriched terms from both DEG lists were ranked according to false discovery rate (FDR) and reported (Supplementary Data Sheet S7). ORA was also applied to all clusters separately using the same cut-offs (Supplementary Data Sheet S7).

ORA of the human DEG list (Supplementary Data Sheet S7, section 7.2). Most enriched terms were related to cytoskeleton/extracellular matrix reorganization (GO:0030335, GO0030198, hsa04510, WP306, WP3932) and ossification (GO:0001503, GO:0001649, GO:0030278, WP322, WP474, WP4808). Signaling pathways related to ERK1/ERK2 (GO:0070372, GO:0070374), MAPK (hsa04010), EGFR (hsa01521), Toll-like receptor-related signaling (hsa04620, WP75), and HIF-1 signaling (hsa04066) were also identified. Community detection revealed 4 clusters compromising altogether 73 DEGs. Cluster-specific ORA revealed that cluster #1 contained genes participating in ossification, or more generally in development/differentiation of mineralized tissue including bone and tooth. Being the largest cluster, cluster #2 contained genes relevant to the proliferation and migration of smooth muscle cells and other cell types, but also included pathways related to inflammation (i.e., rheumatoid arthritis, IL18, IL17 and IL1R signaling pathways), and genes from the TNF- and HIF1-signaling pathways. Though cluster #3 contained 22 DEGs, only 4 to 7 of these were assigned to enrichment terms. The terms with ≥5 DEGs assigned were “cardiac progenitor differentiation” (7 DEGs, WP2406), “endothelial cell migration” (5 DEGs; GO:0043542), and “EGFR tyrosine kinase inhibitor resistance” (5 DEGs; hsa01521 and WP4806). Since the smallest cluster (cluster #4) only contained 5 DEGs, a clear attribution to one of the three databases (GO/Biological Process, KEGG, WikiPathways) was not achieved. The 11 identified hub genes were mostly related to clusters #1 (IGF1, RUNX2) and #2 (IL1B, IL6, JUN, MMP1, MMP2, TGFB1, TIMP1, VEGFA). FGF2 was the only hub gene in cluster #3, and cluster #4 contained none.

ORA of the mouse DEG list (Supplementary Data Sheet S7, section 7.3). Enriched terms from GO/Biological Process were mostly related to organ or tissue development (GO:0048732, GO:2000027, GO:0050678, GO:0061138, GO:0060562, etc.), ossification (GO:0030278, GO:0001503) and WNT-related signaling (GO:0016055, GO:0060070). Ossification and WNT-related signaling were also identified with KEGG and WikiPathways. Additional pathways were identified related to HIPPO- (mmu04390) and TNF-signaling (mmu04668), and focal adhesion (WP85, mmu04510). Community detection revealed 4 clusters compromising altogether 74 DEGs. The first cluster contained 21 DEGs, that were mainly assigned to WNT and BMP signaling, but also to general developmental processes (i.e., embryonic or dermatome development) and HIPPO signaling. In all enrichment terms, at least one of the genes of the WNT family (WNT1, WNT3A, WNT4, WNT5A) was listed. Cluster #2, consisting of 33 DEGs, contained enrichment terms describing cell proliferation and/or migration of different cell types (i.e., bone, smooth muscle cells, leukocytes). Other enrichment terms were related to bone remodeling processes (bone remodeling, bone resorption, endochondral ossification, osteoclast signaling). Clusters #3 and #4 each contained 10 DEGs. With the DEGs contained in cluster #3 none of the enrichment terms satisfied the applied criteria (min of 4 DEGs from the list). In cluster #4, only the enrichment term “biomineral tissue development” (GO:0031214) satisfied the applied criteria. The 12 identified hub genes were mostly related to clusters #1 (DKK1, LRP5, SFRP1, WNT1, WNT3A, WNT4, WNT5A) and #2 (ACTB, IL6, JUN, TNF), whereas one hub gene was found in cluster #3 (CTNNB1).

 3.9 High-throughput gene expression studies and gene set enrichment analysis (GSEA)

Five different human and mouse high-throughput gene expression studies were identified with raw data available in public repositories as described (Figure 2B). In two studies, FSS was applied to human bone marrow-derived mesenchymal stem cells (hMSCs). The effect of FSS on the mouse osteocyte cell line MLO-Y4 was analyzed in the other three studies. In this section, we will only focus on those clusters that were shared by either both (hMSC group) or at least two (MLO-Y4 group) GSEA results.

Human mesenchymal stem cells (hMSCs) (Figure 6; Supplementary Data Sheet S8, section 8.4). Both studies shared numerous pathways that were significantly enriched. These were related to inflammation (e.g., WP2865, WP4493, WP4481, HSA04657, WP4754, WP2637), but also to TNF- (HSA04668, WP 2036, HSA04668), MAPK- (WP382, HSA04010), and Toll-like receptor signaling (WP75, HSA04620). Several pathways related to metabolism (e.g., vitamin B12, WP1533; folate, WP176; arachidonic acid, HSA00590) were also significantly enriched, as well as the cytosolic DNA-sensing pathway (WP4655, HSA04623).

Pre-ranked gene set enrichment analysis (GSEA) of high-throughput gene expression studies investigating the effect of FSS on human mesenchymal stem cells (human MSCs, left column) or mouse MLO-Y4 osteocytes (right column). Differential gene expression was evaluated for each study separately. Ranked gene lists were generated and used for GSEA using easyGSEA; for each cell type the GSEA results were compared with easyVizR (Cheng et al., 2021) and the shared gene sets were further analyzed (yellow column in each of the upset-plots). The shared gene sets with p < 0.05, and FDR<0.1 in any of the individual studies were retained. The top 50 gene sets from hMSCs and all gene sets from mouse MLO-Y4 were visualized using dot plots. NES, normalized enrichment score; FDR, false discovery rate.

Mouse osteocyte cell line MLO-Y4 (Figure 6; Supplementary Data Sheet S8, section 8.5). The easyVizR comparison of the results obtained from easyGSEA revealed, that at least two of the three studies shared pathway terms that were significantly enriched. Terms that were found to be shared between two studies were related to ribosomal RNA processing (GO:0006364) or more generally to ribosome biogenesis (GO:0042254, GO:0042274), but also to inflammatory response (GO:0006954, GO:0050729), or chemotaxis (GO:0006935, GO:0030593). Terms enriched in all three studies were related to the NF-KappaB- (MMU04064), TNF- (MMU04668), IL17-signaling (MMU04657) and cytokine-cytokine receptor interaction (MMU04060).

 4 Discussion

The purpose of this systematic review was to identify studies that investigated the effect of fluid flow shear stress (FSS) on different cell types involved in bone remodeling, to analyze their molecular findings and to summarize commonly used fluid flow apparatuses, FSS magnitudes, durations and mechanobiological responses such as gene expression and metabolites release using protein-protein interaction (PPI) networks and enrichment analysis. Defined criteria were used to evaluate the risk of bias of each study using tailored clear definitions.

4.1 Commonly used apparatuses to apply FSS

Generally, a fluid flow device is considered a system that consists of a fluid-flow generating part (e.g., cone and plate, orbital shaker, rocking plate, peristaltic/gear/syringe/osmotic pump with or without gravity support), a fluid flow chamber, a pulse dampener, and other auxiliary parts such as bubble trap, digital bubble, or flow sensor. For in vitro fluid flow experiments, FSS apparatus selection may be influenced by the pattern of fluid flow, the fluid generating part, or the downstream analytics. For instance, different flow patterns may require different flow driving forces that necessitates a specific kind of pump or use a modified fluid flow system to convert one flow profile into another flow profile. The use of a pulse dampener will convert a pulsating fluid flow pattern into a steady laminar fluid flow pattern. Methods applied for downstream analysis may also play a role in selecting an appropriately sized fluid flow chamber for a given experiment. Chambers used in microscopy tend to be smaller and more delicate compared to larger parallel flow chambers or rocking and orbital shaking platforms. These provide a larger surface for cell monolayer culturing yet yield higher output of mechanically stimulated samples for further analysis. However, smaller parallel flow chambers are capable of exerting uniform and higher FSS magnitudes compared to rocking and orbital shaking platforms, which generate lower but non-uniform (variable) FSS magnitudes (Zhou et al., 2010).

Due to the complexity and unicity of each of the apparatus systems used in each experiment, our identification strategy depended firstly on the kind of flow chamber used. The other parts of the fluid flow apparatus controlling the fluid flow pattern such as the type of driving force (pump or/and gravity), the use of pulse dampener and the general design of the system as defined above were considered secondary.

The most frequently used fluid flow systems from the identified studies consisted of either custom-made or commercially available fluid flow chambers (Flexcell, PeCon, Fluxion Bioscience, Focht or Cytodyne). In both cases, a commercially available driving part (e.g., peristaltic/syringe pump) was used. According to our review, self-designed fluid flow chambers were the most frequently used experimental setup to deliver either of the three flow types (steady laminar, pulsatile laminar and oscillatory laminar fluid flow) independent of the cell type. The possible reason behind this might be that self-designed fluid flow chambers can be customized to serve the required experimental purpose such as geometry, cell culturing area, and/or fluid flow type (Sonam et al., 2016). For primary mouse and human osteoblasts, the Cytodyne fluid flow chamber was the second most frequently used to deliver steady laminar and pulsatile laminar FSS. Chambers from Flexcell Inc. Were the second most preferred for mOcyt and hMSCs to deliver all three fluid flow patterns. Moreover, the three fluid flow systems used for microscopy (Bioflux system, Focht Chamber System 2, and PeCon parallel plate system) were used to deliver either steady laminar, pulsatile laminar or oscillatory laminar fluid flow, respectively. Nevertheless, additional information regarding the Cytodyne fluid flow system and the PeCon parallel plate chamber were not identified.

Finally, the fluid flow profile applied should be defined more precisely. Using the term “laminar” to describe the fluid flow profile only can be confusing since three different fluid flow patterns (Figure 1; “steady laminar”, “pulsatile laminar” and “oscillatory laminar”) are considered in biomechanics. The term “laminar” describes that the flow is smooth with layers or lamina, and it can be adjusted to either steady or unsteady. The impact of the fluid flow patterns on cellular expression will be discussed in the next section.

 4.2 FSS magnitudes and durations

Independent of the cell type, in any included study either of the three fluid flow types (Figure 1) was applied at different FSS magnitudes. This selection was based on the study’s aim and the required in vitro parameters like cell type, FSS magnitude and type of fluid flow to mimic an actual clinical or physiological situation, or to investigate the effect of FSS on gene expression. The application of these parameters in vitro simulations facilitates a more physiologically relevant in vitro setup and thus a robust analysis. In this systematic review, we found that the FSS parameters were chosen according to.

1) Theoretical estimation of the physiological FSS in the canaliculi-lacunae network in bone tissue and blood vessels. These studies have calculated physiological levels of fluid flow rate depending on the porous geometry of tissue, poro-elasticity and fluid properties using mathematical models (Table 7). Different levels of shear stress acting on the various cell types have been identified and were attributed to the differences in their micro-geometry (blood vessels, bones, muscles, ligament) (Table 7) and differences in the circulating fluids (e.g., blood versus serum) (Weinbaum et al., 1994; Zeng et al., 1994; Davies, 1995; Chen et al., 1998; Wittkowske et al., 2016). The average wall shear stress in large arteries with uniform geometry ranges between 20 and 40 dyn/cm2 (Davies, 1995). Other studies have estimated the interstitial FSS in the tunica media of an artery to ∼1 dyn/cm2 (Tada and Tarbell, 2000; Wang et al., 2007). On the bone level, a wall shear stress of 8–30 dyn/cm2 in lacunar-canalicular canals of a trabecular bone model was calculated (Weinbaum et al., 1994), whereas in the canaliculi of an osteon it ranged between 6 and 30 dyn/cm2 (Zeng et al., 1994). However, while these studies did not take into consideration other bone cell types such as osteoblasts and osteoclasts, giving that they occur at sites with larger porosities (surface of newly formed or removed bone), it was suggested that these cells may encounter shear stresses less than 8 dyn/cm2 (Wittkowske et al., 2016). Though fluid flow in tissues of defined geometry (e.g., blood vessels, bone lacunar canaliculi) are easier to predict mathematically, complex tissue geometry (e.g., ligaments, tendons) may require computational finite element analysis to estimate the true mechanical effect of fluid flow within these tissues.

2) Finite element analysis (FEA) represents an important method to simulate mechanobiological interaction using 3D-models or pre-existing measurements to build accurate biological geometry. Recently, a model of the lacunar-canalicular space from rat tibiae was obtained using confocal microscopy and 3D reconstruction (Verbruggen et al., 2014). An interstitial FSS of ∼110 dyn/cm2 in the osteocytes’ mechanical microenvironment after simulated vigorous physical activity was calculated using FEA. This approach has also been used to develop a model of fibrillar arrays in ligaments and tendons and a maximum FSS of 12.1 dyn/cm2 within these tissues was estimated (Chen et al., 1998) (Table 7).

3) Ex vivo studies. Based on an ex vivo model of murine tibia and mathematical modelling, it was shown that osteocytes elastically deform due to the application of a physiological shear stress of ∼50 dyn/cm2 (Price et al., 2011). Advances in microscopic life imaging enabled to determine the effect of FSS on osteocytes ex vivo (Verbruggen et al., 2015) or in vivo (Lewis et al., 2017; Lewis et al., 2021) in murine models. The bone was mechanically loaded using a 3-point-bending setup. Either strain distribution and deformation of osteocytes was determined (Verbruggen et al., 2015) or the increases in calcium signaling in osteocytes directly resulting from FSS (Lewis et al., 2017; Lewis et al., 2021) have been measured giving direct evidence of FSS effects on bone cell activity.

4) In vitro studies have been applied to study the sensitivity of different cell types including endothelial cells and osteocytes. Similar to endothelial cells in vitro, osteocytes respond to magnitudes of shear stress as low as 6 dyn/cm2 (Klein-Nulend et al., 1995b). This confirms the theoretical estimation of FSS in the canaliculi of osteons (Zeng et al., 1994). In addition, a positive correlation between FSS strength and bone cellular response has been suggested (Rath et al., 2010; Li et al., 2012).

5) Fluid flow profiles in human tissues have been an ongoing research topic due to their relationship to normal and abnormal physiological activity. Using a simple representation of the walking cyle, a reversed direction of fluid flow in bone canaliculi during relaxation and compression was predicted (Piekarski and Munro, 1977). These authors also highlighted the change in hydrostatic pressure within the lacunar-canalicular system during loading and unloading. This change of pattern within human tissues governs the fluid flow profile of fluid movement, e.g., a oscillatory fluid flow profile correlates with the normal cyclic activity of the skeletal system such as walking and running (Lu et al., 2012). In contrast, the change in human posture from sitting to standing would result in a physiological but unusual unidirectional steady fluid flow (Jacobs et al., 1998). In periodontal ligament (PDL), compressive force application leads to an increase of pore pressure and an outflow of interstitial fluid from the PDL porosities, whereas the inflow of interstitial fluid in tension areas resulted from a decrease in pore pressure (Ashrafi et al., 2020). Therefore, a change in interstitial flow direction in these tissues during cyclic loading and unloading is closely resembling the oscillatory fluid flow pattern.

A summary of estimated fluid flow shear stress in different fluid contained porous tissues.

Porous tissue
Blood vessels Bone Bone Bone Ligaments and tendon Periodontal ligament (PDFL)
Cell type Endothelial cells Osteocytes (trabecular bone model lacunar canaliculi) Osteocytes (osteon canaliculi model) Osteoblasts Fibroblasts PDL cells
References Davies (1995) Weinbaum et al. (1994) Zeng et al. (1994) Wittkowske et al. (2016) Chen et al. (1998) Unknown
Wall shear stress 20–40 dyn/cm2 8–30 dyn/cm2 6–30 dyn/cm2 <8 dyn/cm2 12.1 dyn/cm2 (max. Shear stress) Unknown
Reasoning Blood is directly pumped by heart The osteocytic processes are extended in tiny channels surrounded by calcified bone matrix (canalicular walls) The osteocytic processes are extended in tiny channels surrounded by calcified bone matrix (canalicular walls) These cells are located at sites of big porosities on the surface of newly formed bone n.g Unknown
Medium Blood Interstitial fluid Interstitial fluid Interstitial fluid Interstitial fluid Interstitial fluid
Type of flow Big blood vessels: pulsatile laminar Small blood vessels: steady laminar Constant loading: steady laminar Intermittent loading: oscillatory laminar Constant loading: steady laminar Intermittent loading: oscillatory laminar Constant loading: steady laminar Intermittent loading: oscillatory laminar Constant loading: steady laminar Intermittent loading: oscillatory laminar Constant loading: steady laminar Intermittent loading: oscillatory laminar
Note Was estimated Was estimated ** physiological frequency for loading and unloading in long bone that correspond with oscillatory fluid flow is 1–20 Hz Salvi et al. (2010) Was estimated Difficult to estimate the FSS values: Constant remodeling Unknown mechanical properties of soft osteoid Was estimated

Evidence on physiological fluid flow profiles within human tissue is important to understand the cellular biophysical environment. Oscillatory fluid flow patterns have been investigated by in vitro studies as being the most likely fluid flow profile experienced by bone cells in vivo (Jacobs et al., 1998; Ponik et al., 2007). Though several studies emphasized the importance of oscillatory flow in maintaining bone integrity (Genetos et al., 2007; Ponik et al., 2007; You et al., 2008; Li et al., 2012), there is also evidence, that the oscillatory fluid flow profile is less stimulatory to bone cells than steady laminar or pulsatile fluid flow profiles (Jacobs et al., 1998). Moreover, bone cell adaptation to common mechanical stress makes them more sensitive to excessive mechanical stress patterns (Turner et al., 2002; Lu et al., 2012). Thus, new FSS magnitudes or profiles are more stimulating to osteocytes than regular FSS magnitudes/profiles.

Herein, oscillatory, and pulsatile fluid flow were the most frequently applied fluid flow profiles in hOst, mOst, mOcyt, and hMSCs, whereas steady laminar fluid flow was mostly applied to hPDLCs. As pointed out, the selection of a fluid flow profile is dependent on the experimental setting, especially the presence of dynamic loading and unloading. In this context it should be mentioned, that the application of static loading to a tissue (e.g., the periodontal ligament), may result in a fluid flow profile that differs from applying cyclic loading and unloading commonly seen in daily human activity such as walking or food chewing (Ashrafi et al., 2020). Based on the in-/outflow of interstitial fluid during compression- or tension-related tissue deformation (Figure 1), it is suggested that the periodontium is subjected to a prominent one-directional dynamic fluid flow during the initial phase of orthodontic tooth movement (OTM) (Ashrafi et al., 2020), which becomes less dominant during the lag-phase of OTM (Reitan, 1960). The duration of the lag-phase is determined by the time required for periodontal remodeling and consequently determines the loading/unloading frequency on tissues, thus affecting fluid flow character (Reitan, 1960).

Independent of the studied cell type, the most frequently selected durations of FSS exposure were 0.5, 1 and 2 h. The main reasoning for this selection was based on the following criteria: (i) optimal time frame to determine cellular biological responses; (ii) technical issues. For instance, from a study on oscillatory fluid flow profile it was concluded that FSS durations of 1, 2, and 4 h resulted in an upregulation of PTGS2, TNFSF11, and TNFRSF11B gene expression favoring bone formation (Li et al., 2012). These authors also pointed out, that longer FSS applications than 4 h might require a different pH buffering system of the cell culture medium without negatively influencing the mechanical properties of the chamber. This and other drawbacks of parallel flow chambers (Huesa et al., 2010), like the formation of air bubbles that may impair the cellular in vitro environment (Anderson et al., 2006), or the negative impact of continuous cellular growth and matrix deposition on flow chamber geometry (Nauman et al., 2001) should be taken into consideration.

Intriguingly, the most used shear stress was 10 dyn/cm2 regardless of FSS duration and fluid flow profile, which corresponds with the physiological shear stresses (Weinbaum et al., 1994; Zeng et al., 1994). In addition, frequencies of 1 Hz and 5 Hz, which were commonly used in most cell groups, lie within the physiological range of bone loading and unloading caused by locomotion or maintaining posture (Weinbaum et al., 1994).

Finally, it is important to point out that the presence of endothelial cells, bone cells, mesenchymal stem cells and fibroblasts in proximity to each other may facilitate intercellular communication through body fluids including interstitial fluid (Cenaj et al., 2021).

 4.3 Genes or metabolites reported by at least three studies

Herein, we identified genes and metabolites (AKT1, alkaline phosphatase, BGLAP, BMP2, Ca2+, COL1A1, CTNNB1, GJA1, MAPK1/MAPK3, nitric oxide, PDPN, PGE-2, PGI-2, PTGS1, PTGS2, RUNX2, SOST, SPP1, TNFRSF11B, TNFSF11, VEGFA, WNT3A), which were investigated by at least three studies on the effect of FSS on cells of human and mouse origin such as MSCs, PDLCs, and bone cells. These genes/metabolites were grouped based on their functional role during periodontal ligament remodeling as being responsible for tissue formation (AKT1, alkaline phosphatase, BGLAP, BMP2, Ca2+, COL1A1, CTNNB1, GJA1, MAPK1/MAPK3, PDPN, RUNX2, SPP1, TNFRSF11B, VEGFA, WNT3A), tissue degradation (SOST, TNFSF11) or inflammation (nitric oxide, PGE-2, PGI-2, PTGS1, PTGS2).

4.3.1 Human cells 4.3.1.1 Tissue formation

The periodontium has a complex ultrastructure including periodontal ligament and bone that is created by various cell types (Marchesan et al., 2011). Periodontal cells are involved in the perpetual remodeling of the periodontium thereby interacting with each other. This process is mediated and regulated by numerous genes/metabolites involved in cell differentiation, proliferation and migration, extracellular matrix formation and tissue maturation.

In addition to the regulatory rule of Ca2+ during inflammation and its multiple functions in various biological systems, several studies have highlighted the importance of intracellular Ca2+ release in hMSCs proliferation and osteogenic differentiation (Riddle et al., 2006; Hu et al., 2017). Our systematic review has identified six studies in which all of them reported an increase of intracellular calcium release by hMSCs after FSS application. This is supported by recent studies, which show direct evidence of FSS effects on bone cell activity, i.e., increased calcium signaling in osteocytes directly resulting from FSS (Lewis et al., 2017; Lewis et al., 2021).

The activity of alkaline phosphatase (ALPL) in mesenchymal cells has been suggested as an indicator of osteogenesis (Lim et al., 2013). Herein six studies were identified that determined the effect of FSS on ALPL using sqPCR, RT-qPCR or ELISA in hMSCs. A large proportion of studies (>80%) reporting upregulation or release of ALPL in hPDLCs and hMSCs suggest that FSS induces osteogenic differentiation in these cells (Lim et al., 2013; Lim et al., 2014; Qi and Zhang, 2014; Tang et al., 2014) (Figure 4; Supplementary Data Sheets S5, 6). In addition, FSS has been found to upregulate COL1A1 in hMSCs, which provides structural support and strength to various tissues and plays a crucial role in bone formation and integrity (Lim et al., 2013).

Several included studies reported FSS-dependent upregulation of VEGFA, MAPK1/MAPK3 and RUNX2. However, while both MAPK1/MAPK3 and RUNX2 are important regulators of osteoblastic proliferation and regulation (Riddle et al., 2006; Kuo et al., 2015), the number of studies reporting MAPK1/MAPK3 upregulation was higher than reports on RUNX2 upregulation in FSS-stimulated hMSCs. Nevertheless, the relevant studies investigating VEGFA have shown 100% upregulation in hMSCs in relation to FSS (Figure 4). This result was expected since VEGFA is not only involved in vascular development and angiogenesis but also in bone healing (Street et al., 2002). Noteworthy, hMSCs were sensitive to very low oscillatory FSS of 0.01 dyn/cm2 lasting for up to 3 weeks, which led to an upregulation of genes involved in tissue formation such as ALPL, COL1A1, BGLAP, SPP1, VEGFA, and RUNX2 (Lim et al., 2013).

 4.3.1.2 Inflammation

Inflammation has a significant effect on bone remodeling through bone resorption (Epsley et al., 2020). It is also considered as one of the earliest reactions to orthodontic mechanical force on the cellular level (Li et al., 2018). Several molecules and metabolites such as PGE-2 (a derivative of arachidonic acid and product of cyclooxygenase enzymes), nitric oxide (NO; product of nitric oxide synthases, NOS1-3) and Ca2+ (intracellular calcium ions) are produced within tissues and play a crucial role in modulating the inflammatory responses. FSS has been suggested to increase the level of NO and PGE-2 in various biological systems such as osteoblasts, osteocytes, and periodontal ligament cells (van der Pauw et al., 2000; Bakker et al., 2013a; Bakker et al., 2013b; van der Meijden et al., 2016; Deepak et al., 2017). The regulation of NO and PGE-2 was suggested to be dependent on Ca2+ signaling transduction (Ajubi et al., 1999). Independent on the human cell type, at least three relevant studies reported exclusive upregulation of NO, Ca2+, and PGE-2 in response to FSS. The mobilization of these metabolites by bone cells in response to FSS has been found to regulate tissue remodeling through an inflammatory pathway thereby acting as a second messenger (NO, Ca2+) or locally as a hormone-like substance (PGE-2) (Klein-Nulend et al., 1995a; Hung et al., 1995; Reich et al., 1997; McAllister and Frangos, 1999; You et al., 2001).

PTGS1 and PTGS2 are two isoforms of the prostaglandin-endoperoxide synthase (PTGS; aka cyclooxygenase, COX) that are upregulated during inflammation. Both genes are upregulated after FSS application in hOst (Joldersma et al., 2001; Bakker et al., 2003). However, the upregulation of PTGS2 has been found to be stronger and persisted up to 24 h compared PTGS1. Though we identified 3/3 (100%) studies reporting an FSS-dependent upregulation of PTGS1 and PTGS2, it has been suggested that PTGS2 is the mechanosensitive form of two PTGS isoforms that regulates the inflammatory response of bone to mechanical loading (Bakker et al., 2003).

 4.3.2 Mouse cells 4.3.2.1 Tissue formation

Similarly, several studies on mouse bone cells have shown the FSS-dependent upregulation of genes and metabolites commonly expressed in hMSCs, hOst and hPDLCs. These include the upregulation of NO and BGLAP (100%) by mOst and mOcyt and Ca2+ (93%) by mOcyt in response to FSS.

Interestingly, the immortalized mouse osteocyte cell line MLO-Y4 was investigated in a larger number of included studies than primary human osteocytes. This could be related not only to its availability and reproducibility but also to the multiple roles of osteocytes in bone formation, maintenance, and remodeling. At least three of the studies have shown an upregulation (>80%) of AKT1, PDPN, CTNNB1, GJA1, TNRFSF11B, VEGFA and WNT3A expression after FSS (Figure 4). The expression of AKT1 has been considered essential for the functional opening of CX43 hemichannels and the release of central anabolic factors in response to FSS (Riquelme et al., 2021). It was also reported that PGE-2 is essential in promoting gap junction-mediated communication in osteocytes by controlling the up and downregulation of GJA1 (Xia et al., 2010).

TNRFS11B/TNRFS11 and the WNT/β-catenin signaling pathways have been considered to play a key role in regulating bone architecture and bone homeostasis (Bonewald and Johnson, 2008; Macsai et al., 2008; Baron and Kneissel, 2013; Liao et al., 2017; Yan et al., 2018). In addition, the upregulation of PDPN has been regarded as a reflection of the conversion of osteoblasts to osteoid cells or osteocytes (Zhang et al., 2006). Interestingly, mouse MLO-Y4 cells are sensitive even to very small pulsatile FSS of 2 dyn/cm2 (Huang et al., 2017). This FSS magnitude upregulated WNT3A expression and increased PGE-2 release (Huang et al., 2017).

 4.3.2.2 Tissue degradation

While an increasing number of reports have shown upregulation of the above-mentioned genes in mOcyt, only few studies have reported upregulation of genes that are responsible for bone resorption such as SOST1, TNFSF11 or TNFSF11B (Supplementary Data Sheets S5, 6). Only 25% (1/4) of studies reported a FSS-dependent upregulation of SOST1 in mOcyt leading to a downregulation of the Wnt/β-catenin pathway and thus inhibiting osteogenesis (Li et al., 2013). Variable expression of TNFSF11 or TNFSF11B was observed after application of different FSS profiles. High magnitudes of steady laminar FSS (16 dyn/cm2; 30 dyn/cm2) were shown to support the upregulation of TNFSF11 or TNFSF11B (Li et al., 2013) in comparison to low magnitudes (10 dyn/cm2) (Haugh et al., 2015; Yan et al., 2018). However, pulsatile fluid flow profiles at lower magnitudes (≈10 dyn/cm2) upregulated the expression of TNFSF11 or TNFSF11B (Kulkarni et al., 2010; Kulkarni et al., 2012; Bakker et al., 2013b; Liao et al., 2017). Interestingly, only 12.5% (1/8) of the studies using mouse osteocytes have shown an increase in the RANKL/OPG ratio and thus favoring bone resorption (Supplementary Data Sheets S5, 6).

 4.3.2.3 Inflammation

In addition to the inflammatory mediators mentioned above for human cells upregulation of PGE-2, PTGS2, and nitric oxide were found in the included studies. All inflammation-related genes and metabolites were solely upregulated by FSS independent of the specific mouse cell type.

 4.4 General observations

The categorization of the genes/proteins/metabolites used herein (inflammation, tissue formation, tissue degradation) was used as a structural framework to improve data description and analysis thus enhancing the comprehensibility of the complex interactions between the different molecular pathways involved on various levels. Still, the results derived from studies on both human and mouse cells suggest dominant tissue formation character of FSS. However, direct comparison of their results was not feasible. This was primarily due to the heterogeneity which exist in different aspects within and between the studies (Sun et al., 2021), such as: 1) differences related to the cells used (e.g., cell density, cell isolation methods, passaging number, donor-related differences (sex, age, etc.), cell culturing materials and cell maintenance conditions); 2) differences related to the applied analytical techniques (quantification methods, reference gene selection, differences in expression reporting, and post-transcriptional modification); 3) differences in the experimental set-up like timings of sample collection, fluid flow profile used, duration of FSS application/post-FSS incubation, FSS magnitude and/or frequency, and data collection and analysis.

The effect of donor age and donor sex on cell physiology is well established. Due to the high heterogeneity of studies along with partially missing reports on donor age and/or sex, these variables were not considered for stratification in data analysis. Nevertheless, age-related changes relevant to FSS application have been identified in the literature and been summarized herein to emphasize their significance for further research. Aging is a physiological process which is inevitably bound to various cellular, mechanical, and spatial changes within tissues, including bone (Brodt et al., 1999). In osteocytes, aging is associated with a marked reduction of trabecular (Qiu et al., 2002) and cortical osteocyte density (Carter et al., 2013a; Carter et al., 2013b). Lacunar canalicular permeability and organization was also found to be negatively affected by aging, especially in females, leading to a disruption in the bone microarchitecture, reduced fluid flow and thus impairs cellular communication between osteocytes (Currey, 1964; Okada et al., 2002; Rodriguez-Florez et al., 2014; Ashique et al., 2017; Schurman et al., 2021).

Also, the periodontal ligament (PDL) undergoes characteristic age-related cellular, structural, and functional changes. As mentioned above, the PDL is located between the tooth and the alveolar bone and consists of cells of different tissue types embedded in a matrix, which is composed mainly of collagen type 1 and 2 (Karimbux et al., 1992; Marchesan et al., 2011). With increasing age a reduction in PDL thickness is commonly observed (Barczyk et al., 2013). Additionally, in aged populations collagen fibers are poorer organized (Cho and Garant, 1984; Kawada and Komatsu, 2000) and cellular collagen expression is reduced (Oehmke et al., 2004). Therefore, aging also affects the mechanical properties of the periodontal ligament, which dependent on collagen. Nevertheless, the periodontal ligament cells may have impaired function with aging. These altered functions are closely related to periodontal ligament fibroblast, which is suggested to produce more inflammatory mediators and yet is more prone to damage (Abiko et al., 1998). Therefore, it is suggested that age-related cellular changes in bone and periodontal ligament might impair the response to treatment at higher age.

The relationship between fluid flow and age were addressed in two in vitro studies (Joldersma et al., 2000; Klein-Nulend et al., 2002). It was found that apart from reduced proliferation capacity of aged cells, the mechanosensitivity of these cells to fluid flow remained unaffected. However, these in vitro studies did not take into consideration age-related micromechanical (structural) bone changes, which may influence the magnitude of fluid flow within bone and resulting shear stress acting on these cells.

Despite several limitations, we suggest that during orthodontic treatment, interstitial fluid flow occurs within the periodontium at stretching and compression areas around the tooth (Ashrafi et al., 2020). However, stronger bone deposition occurs in tissues subjected to stretching compared to compression (Li et al., 2018). While the vascular supply is reduced in the compression areas, we assume that the effect of interstitial fluid flow, which depends on the blood and lymphatic circulation, remains undisrupted on the tension side (Li et al., 2018). The fluid flow environment and the metabolic similarities between FSS-induced and stretched tissues during OTM remain to be further elucidated in the future.

 4.5 Over-representation, gene set enrichment and network analysis

Due to the heterogeneity of the included studies, a direct comparison of FSS-dependent gene expression was not possible. In contrast to a recently published systematic review (Sun et al., 2021), five studies with publicly available high-throughput gene expression data (microarray, RNA-seq) were identified reporting FSS-dependent gene expression in hMSCs and mOcyt (MLO-Y4). Therefore, the previously published analytical workflow (Sun et al., 2021) was extended to gain additional insights into the molecular pathways regulated by FSS in OTM-related human and mouse cells. The bipartite workflow consisted of the following steps: (1) generation of protein-protein interaction (PPI) networks and application of over-representation analysis (ORA). (2) Re-analysis of five high-throughput gene expression studies. (3) To increase comparability within and between the workflow arms, identical databases were engaged for enrichment.

4.5.1 PPI networks and over-representation analysis (ORA)

Species-specific lists of FSS-dependent differential expressed genes (DEGs) were extracted from all studies identified herein. Reports from high-throughput expression studies (microarray, RNA-seq, proteomics) were excluded, since only the most significantly up- and downregulated genes were reported, and different cut-offs were applied. Instead, gene expression data derived from RT-qPCR reanalysis was extracted and included. To increase specificity, only those genes showing clear FSS-dependent expression were included (Sun et al., 2021). Due to its methodological heterogeneity protein expression data were excluded (ELISA vs. WB quantification vs. enzyme activity; antibody specificity).

From the human gene list, seven genetic loci were not integrated into the network (IER3, TRPV1, TRPV4, PON1, DVSP6, GADD45B, WIF1), whereas only two genetic loci (ERAL1, PIEZO2) were not integrated into the mouse network. Lowering the STRING-DB cut-off score to 0.4 integrated the first 6 singletons (IER3, TRPV1, TRPV4, PON1, DVSP6, GADD45B) into the human network. All newly established edges originated from text-mining results, and only a few interactions were additionally supported by co-expression or other experimental data. Lowering the cut-off to 0.3 also integrated WIF1 into the network, also due to text-mining results. In the mouse network, decreasing the cut-off score to 0.4 established a PPI-interaction between PIEZO2 and PKD2 of the existing network (text-mining score: 0.45; co-expression score: 0.061). Lowering the cut-off to 0.2 did not integrate ERAL1 into the network but instead established an additional PPI-interaction between PIEZO2 and PKD1 of the existing network (text-mining score: 0.366, co-expression core: 0.062). Generally, a lower cut-off score might integrate singletons into an existing network, but at the same time it might increase the number of PPI interactions between other nodes of the existing network. Though not affecting enrichment directly, the reduction of the cut-off score to 0.3 will affect several network statistics by increasing the total number of edges and the average number of neighbors and decreasing the network diameter and the characteristic path length (not shown). This, in turn, might will affect the number of clusters and hub genes and thus downstream enrichment analysis. With a cut-off score of 0.7, 4 clusters were identified in each network. After application of a cut-off score of 0.3 to both PPI networks, the same clustering algorithms identified 3 clusters in the human and 2 clusters and 1 singleton (Eral1) in the mouse PPI networks (results not shown).

Generally, the findings from the two species-specific networks were consistent with the research topics covered in the included studies (Sun et al., 2021). Enrichment terms related to ossification and cytoskeleton/extracellular matrix reorganization were identified in both species. Since constraints were applied to increase specificity in reporting, the list of enrichment terms might suggest, that terms found in one PPI network were not identified in other networks or its clusters, e.g., BMP-related pathways and biological processes in the mouse PPI network, or ERK1/ERK2-related signaling in the human PPI network. Indeed, these terms were also identified in the other PPI network, but due to the applied selection criteria, these terms did not pass this cut-off. For example, ERK1/ERK2 signaling (GO:0070372) was represented by 18 genes from the human gene list covering 6.2% of the genes belonging to this biological process (FDR: 4.16E-14). ORA identified ERK1/ERK2-related terms also in the mouse list, but none was reported due to either criteria: “Regulation of erk1 and erk2 cascade” (GO:0070372; 11 genes, 3.5%, FDR: 3.88E-7), “Positive regulation of erk1 and erk2 cascade” (GO:0070374; 8 genes, 4.1%, FDR: 2.19E-6), and “ERK1 and ERK2 cascade” (GO:0070371; 3 genes, 9.7%, FDR 2.50E-3). The alike also applied to the “Hippo signaling pathway” (mmu04390) identified in the mouse PPI network (13 genes, 8.3%, FDR: 6.32E-13). Pathways related to Hippo signaling were also identified in cluster#3 of the human PPI network but failed to be reported due to the applied constraints: “Hippo-Merlin signaling dysregulation” (WP4541; 5 genes, 4.2%, FDR:1.92E-5) and “H.ippo signaling regulation pathways” (WP4540; 4 genes, 4.1%, FDR: 2.90E-4).

 4.5.2 Selection of high-throughput gene expression studies

The search strategy used herein identified several studies that applied high-throughput expression platforms to study FSS-dependent expression in one of the selected cell types: hMSCs were analyzed with antibody arrays (Lim et al., 2013) or microarrays (Kuo et al., 2015). Microarray platforms were also applied to study hPDLCs (Qi and Zhang, 2014) and mOsts (Lau et al., 2006). Mouse osteocytes were analyzed with microarrays (Govey et al., 2014; Kitase et al., 2014), RNA-seq (Govey et al., 2015; Li X. et al., 2019; Du et al., 2020) and proteomics (Govey et al., 2014). In line with ORA, we focused on gene expression platforms (microarray, RNA-seq) and identified three studies with publicly available raw data (Govey et al., 2014; Govey et al., 2015; Li et al., 2019b). By hand-searching of public repositories two additional studies with available raw data were identified (Glossop and Cartmell, 2009; Diaz et al., 2017).

 4.5.3 Pre-ranked GSEA

Several enriched pathways and biological processes were identified as shared in each species, indicating that despite differences in experimental conditions FSS-related gene expression might be similar. Nevertheless, differences were observed concerning enrichment statistics (enrichment score, associated p values) and genes belonging to the leading edges.

During the reanalysis of the data sets used for GSEA it became apparent, that each step in DEG analysis will influence the pre-ranked gene lists from each study and thus the subsequent GSEA. Originally, the workflow based on the “limma” Bioconductor package (Ritchie et al., 2015) was considered for both microarray and RNA-seq data. Finally, “limma” was chosen for the microarray and the “DESeq2” Bioconductor package (Love et al., 2014) for RNA-seq data analysis. Additionally, the inclusion of different microarray platforms required adaptations to the originally considered workflows (Love et al., 2015; Klaus and Reisenauer, 2016; Law et al., 2016). Additional steps during the workflow related to normalization, unsupervised filtering of non- or low-expressed probes/genes, model-selection (design-matrix, covariates, blocking factor) and the method to assess differential gene expression will influence the DEG results. Ranking of the gene lists will also influence the results of GSEA (Reimand et al., 2019). Usually, a score is used, that describes the level of differential expression for each gene, e.g., according to the t-statistics from the limma-workflow, according to log2FoldChange, or combinations of statistics and log2FoldChange (sign (log2FC)×-log10(p-value) or log2FC×-log10(p-value)).

The source of enrichment terms (pathways, disease state, transcription factor binding sites, miRNA or lncRNA interaction data, etc.) will also affect the outcome of an enrichment analysis, independent if ORA or GSEA is applied. Additionally, though sometimes named (nearly) identical, terms derived from different sources contain different gene sets. As such the KEGG “Toll-like receptor signaling pathway” (hsa04620) and WikiPathways’ pathway with the same name (WP75) share 99 genes out of 104 (KEGG) respectively 103 (WikiPathways). Therefore, the application of different enrichment sets and implementations for the actual enrichment is frequently advised in the literature (Reimand et al., 2019). Based on the experience gained during reanalysis, one bottleneck in the described workflow is the availability of resources to compare the data from different sources. For example, easyGSEA results were compared with the online companion program easyVizR (Cheng et al., 2021). For this comparison, information on the leading edge’ genes, the enrichment score, the p-value and the FDR or adjusted p-value were needed for each enrichment term from each of the studies. As such, a group-wise comparison of the pre-ranked GSEA results generated by STRING-DB (Szklarczyk et al., 2019; Szklarczyk et al., 2023) was not possible, since the unadjusted p-value was not given with the results from STRING-DB.

 4.6 Strengths and limitations

Herein discussed evidence implicates that fluid flow is an important player in bone and periodontal tissue remodeling and especially in formation. To our knowledge, this is the first review giving a comprehensive overview and discussion of methodological technical details regarding fluid flow application in 2D cell culture in vitro experimental conditions. It is also providing valuable information about cellular molecular events and their quantitative and qualitative analysis.

Even though, the current analysis was done under an orthodontic point of view, specifically in terms of orthodontically induced bone remodeling, it must be clearly emphasized that most of the data considered herein have been reported from studies not exclusively addressing the orthodontic-related research questions issues. While this might be considered as limitation it is clearly indicating the need for more evidence on cells from the orofacial skeleton. Moreover, it is also providing a profound basis for these studies, concerning both molecular background and technical design.

Moreover, the considerable methodological and biological-related heterogeneity among studies should be considered when interpreting the current results. To increase the number of standardized and reproducible studies in the future there is a clear necessity for the development of an optimized fluid flow model, carefully adapted to the specific research requirements.

Even though this review tried to categorize all studies considering different fluid flow modes and cellular types, not all variables that affect bone remodeling were encompassed. This includes some donor-related characteristics that might affect bone remodeling, such as age and sex. Due to its complexity and importance for tissue remodeling, this topic has yet not sufficiently addressed and might therefore comprise an important subject for future work.

 Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

 Author contributions

MN: Conceptualization, Data curation, Formal Analysis, Methodology, Software, Validation, Visualization, Writing–original draft, Writing–review and editing, Investigation. MF: Writing–review and editing, Supervision. AW: Supervision, Writing–review and editing, Resources. UB: Resources, Writing–review and editing, Conceptualization, Data curation, Formal Analysis, Funding acquisition, Methodology, Project administration, Software, Validation, Visualization, Writing–original draft. MJ: Conceptualization, Investigation, Methodology, Project administration, Validation, Writing–original draft, Writing–review and editing.

 Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. MN was supported by a scholarship from the “Research Grants—Doctoral Programmes in Germany” scholarship program of the German Academic Exchange Service (Deutscher Akademischer Austauschdienst, DAAD; personal reference number: 91796042).

 Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

 Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

 Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fbioe.2023.1256825/full#supplementary-material

PRISMA 2020 exclusion criteria.

Data extraction from studies on human and mouse cells; separately including references.

Risk of Bias assessment; updated extraction sheets; RoB methodology and methodology results for human and mouse cells separately including references.

Summary tables on FSS magnitude/duration and apparatuses.

Analysis of commonly investigated genes and metabolite mentioned by at least three studies.

Commonly investigated genes and metabolite mentioned by at least three studies.

Detailed results from gene list analysis.

Detailed results from gene set enrichment analysis of the high-throughput studies.

https://www.genenames.org/

https://www.informatics.jax.org/genes.shtml

https://www.ncbi.nlm.nih.gov/gene/

https://www.ncbi.nlm.nih.gov/tools/primer-blast/

https://string-db.org

https://www.ncbi.nlm.nih.gov/geo/

https://www.ebi.ac.uk/biostudies/arrayexpress

https://www.ncbi.nlm.nih.gov/sra/

https://jhubiostatistics.shinyapps.io/recount3-study-explorer/

https://tau.cmmt.ubc.ca/eVITTA/(easyGSEA, easyVizR)

 References Abiko Y. Shimizu N. Yamaguchi M. Suzuki H. Takiguchi H. (1998). Effect of aging on functional changes of periodontal tissue cells. Ann. Periodontol. 3 (1), 350–369. 10.1902/annals.1998.3.1.350 Ajubi N. E. Klein-Nulend J. Alblas M. J. Burger E. H. Nijweide P. J. (1999). Signal transduction pathways involved in fluid flow-induced PGE2 production by cultured osteocytes. Am. J. Physiol. 276 (1), E171–E178. 10.1152/ajpendo.1999.276.1.E171 Anderson E. J. Falls T. D. Sorkin A. M. Knothe Tate M. L. (2006). The imperative for controlled mechanical stresses in unraveling cellular mechanisms of mechanotransduction. Biomed. Eng. Online 5, 27. 10.1186/1475-925X-5-27 Ashburner M. Ball C. A. Blake J. A. Botstein D. Butler H. Cherry J. M. (2000). Gene ontology: tool for the unification of biology. Nat. Genet. 25 (1), 25–29. 10.1038/75556 Ashique A. M. Hart L. S. Thomas C. D. L. Clement J. G. Pivonka P. Carter Y. (2017). Lacunar-canalicular network in femoral cortical bone is reduced in aged women and is predominantly due to a loss of canalicular porosity. Bone Rep. 7, 9–16. 10.1016/j.bonr.2017.06.002 Ashrafi M. Ghalichi F. Mirzakouchaki B. Zoljanahi Oskui I. (2020). Numerical simulation of hydro-mechanical coupling of periodontal ligament. Proc. Inst. Mech. Eng. H. 234 (2), 171–178. 10.1177/0954411919887071 Bakker A. D. Huesa C. Hughes A. Aspden R. M. van't Hof R. J. Klein-Nulend J. (2013a). Endothelial nitric oxide synthase is not essential for nitric oxide production by osteoblasts subjected to fluid shear stress in vitro . Calcif. Tissue Int. 92 (3), 228–239. 10.1007/s00223-012-9670-x Bakker A. D. Klein-Nulend J. Burger E. H. (2003). Mechanotransduction in bone cells proceeds via activation of COX-2, but not COX-1. Biochem. Biophys. Res. Commun. 305 (3), 677–683. 10.1016/s0006-291x(03)00831-3 Bakker A. D. Zandieh-Doulabi B. Klein-Nulend J. (2013b). Strontium ranelate affects signaling from mechanically-stimulated osteocytes towards osteoclasts and osteoblasts. Bone 53 (1), 112–119. 10.1016/j.bone.2012.11.044 Barczyk M. Bolstad A. I. Gullberg D. (2013). Role of integrins in the periodontal ligament: organizers and facilitators. Periodontol 63 (1), 29–47. 10.1111/prd.12027 Baron R. Kneissel M. (2013). WNT signaling in bone homeostasis and disease: from human mutations to treatments. Nat. Methods 19 (2), 179–192. 10.1038/nm.3074 Bien S. M. (1966). Fluid dynamic mechanisms which regulate tooth movement. Adv. Oral Biol. 2, 173–201. 10.1016/b978-1-4832-3118-1.50014-1 Bonewald L. F. Johnson M. L. (2008). Osteocytes, mechanosensing and Wnt signaling. Bone 42 (4), 606–615. 10.1016/j.bone.2007.12.224 Brodt M. D. Ellis C. B. Silva M. J. (1999). Growing C57Bl/6 mice increase whole bone mechanical properties by increasing geometric and material properties. J. Bone Min. Res. 14 (12), 2159–2166. 10.1359/jbmr.1999.14.12.2159 Buga G. M. Gold M. E. Fukuto J. M. Ignarro L. J. (1991). Shear stress-induced release of nitric oxide from endothelial cells grown on beads. Hypertension 17 (2), 187–193. 10.1161/01.hyp.17.2.187 Carter Y. Thomas C. D. Clement J. G. Peele A. G. Hannah K. Cooper D. M. (2013a). Variation in osteocyte lacunar morphology and density in the human femur--a synchrotron radiation micro-CT study. Bone 52 (1), 126–132. 10.1016/j.bone.2012.09.010 Carter Y. Thomas C. D. L. Clement J. G. Cooper D. M. L. (2013b). Femoral osteocyte lacunar density, volume and morphology in women across the lifespan. J. Struct. Biol. 183 (3), 519–526. 10.1016/j.jsb.2013.07.004 Cenaj O. Allison D. H. R. Imam R. Zeck B. Drohan L. M. Chiriboga L. (2021). Evidence for continuity of interstitial spaces across tissue and organ boundaries in humans. Commun. Biol. 4 (1), 436. 10.1038/s42003-021-01962-0 Chen C. T. Malkus D. S. Vanderby R. Jr. (1998). A fiber matrix model for interstitial fluid flow and permeability in ligaments and tendons. Biorheology 35 (2), 103–118. 10.1016/S0006-355X(99)80001-8 Cheng X. Yan J. Liu Y. Wang J. Taubert S. (2021). eVITTA: a web-based visualization and inference toolbox for transcriptome analysis. Nucleic Acids Res. 49 (W1), W207–W215. 10.1093/nar/gkab366 Chin C. H. Chen S. H. Wu H. H. Ho C. W. Ko M. T. Lin C. Y. (2014). cytoHubba: identifying hub objects and sub-networks from complex interactome. BMC Syst. Biol. 8 (4), S11. 10.1186/1752-0509-8-S4-S11 Cho M. I. Garant P. R. (1984). Formation of multinucleated fibroblasts in the periodontal ligaments of old mice. Anat. Rec. 208 (2), 185–196. 10.1002/ar.1092080205 Cohn S. A. (1965). Disuse atrophy of the periodontium in mice. Arch. Oral Biol. 10 (6), 909–919. 10.1016/0003-9969(65)90084-1 Cowin S. C. Cardoso L. (2015). Blood and interstitial flow in the hierarchical pore space architecture of bone tissue. J. Biomech. 48 (5), 842–854. 10.1016/j.jbiomech.2014.12.013 Currey J. D. (1964). Metabolic starvation as a factor in bone reconstruction. Acta Anat. (Basel). 59, 77–83. 10.1159/000142603 Davies P. F. (1995). Flow-mediated endothelial mechanotransduction. Physiol. Rev. 75 (3), 519–560. 10.1152/physrev.1995.75.3.519 Deepak V. Kayastha P. McNamara L. M. (2017). Estrogen deficiency attenuates fluid flow-induced [Ca(2+)](i) oscillations and mechanoresponsiveness of MLO-Y4 osteocytes. FASEB J. 31 (7), 3027–3039. 10.1096/fj.201601280R Diaz M. F. Vaidya A. B. Evans S. M. Lee H. J. Aertker B. M. Alexander A. J. (2017). Biomechanical forces promote immune regulatory function of bone marrow mesenchymal stromal cells. Stem Cells 35 (5), 1259–1272. 10.1002/stem.2587 Doncheva N. T. Morris J. H. Gorodkin J. Jensen L. J. (2019). Cytoscape StringApp: network analysis and visualization of proteomics data. J. Proteome Res. 18 (2), 623–632. 10.1021/acs.jproteome.8b00702 Du J. Yang J. He Z. Cui J. Yang Y. Xu M. (2020). Osteoblast and osteoclast activity affect bone remodeling upon regulation by mechanical loading-induced leukemia inhibitory factor expression in osteocytes. Front. Mol. Biosci. 7, 585056. 10.3389/fmolb.2020.585056 Du J. Yuan Z. Ma Z. Song J. Xie X. Chen Y. (2014). KEGG-PATH: kyoto encyclopedia of genes and genomes-based pathway analysis using a path analysis model. Mol. Biosyst. 10 (9), 2441–2447. 10.1039/c4mb00287c Epsley S. Tadros S. Farid A. Kargilis D. Mehta S. Rajapakse C. S. (2020). The effect of inflammation on bone. Front. Physiol. 11, 511799. 10.3389/fphys.2020.511799 Feller L. Khammissa R. A. Schechter I. Thomadakis G. Fourie J. Lemmer J. (2015). Biological events in periodontal ligament and alveolar bone associated with application of orthodontic forces. ScientificWorldJournal 2015, 1–7. 10.1155/2015/876509 Fu Q. Wu C. Shen Y. Zheng S. Chen R. (2008). Effect of LIMK2 RNAi on reorganization of the actin cytoskeleton in osteoblasts induced by fluid shear stress. J. Biomech. 41 (15), 3225–3228. 10.1016/j.jbiomech.2008.08.010 Genetos D. C. Kephart C. J. Zhang Y. Yellowley C. E. Donahue H. J. (2007). Oscillating fluid flow activation of gap junction hemichannels induces ATP release from MLO-Y4 osteocytes. J. Cell. Physiol. 212 (1), 207–214. 10.1002/jcp.21021 Glossop J. R. Cartmell S. H. (2009). Effect of fluid flow-induced shear stress on human mesenchymal stem cells: differential gene expression of IL1B and MAP3K8 in MAPK signaling. Gene Expr. Patterns. 9 (5), 381–388. 10.1016/j.gep.2009.01.001 Govey P. M. Jacobs J. M. Tilton S. C. Loiselle A. E. Zhang Y. Freeman W. M. (2014). Integrative transcriptomic and proteomic analysis of osteocytic cells exposed to fluid flow reveals novel mechano-sensitive signaling pathways. J. Biomech. 47 (8), 1838–1845. 10.1016/j.jbiomech.2014.03.022 Govey P. M. Kawasawa Y. I. Donahue H. J. (2015). Mapping the osteocytic cell response to fluid flow using RNA-Seq. J. Biomech. 48 (16), 4327–4332. 10.1016/j.jbiomech.2015.10.045 Haugh M. G. Vaughan T. J. McNamara L. M. (2015). The role of integrin αVβ3 in osteocyte mechanotransduction. J. Mech. Behav. Biomed. Mat. 42, 67–75. 10.1016/j.jmbbm.2014.11.001 Hu K. Sun H. Gui B. Sui C. (2017). TRPV4 functions in flow shear stress induced early osteogenic differentiation of human bone marrow mesenchymal stem cells. Biomed. Pharmacother. 91, 841–848. 10.1016/j.biopha.2017.04.094 Huang J. Romero-Suarez S. Lara N. Mo C. Kaja S. Brotto L. (2017). Crosstalk between MLO-Y4 osteocytes and C2C12 muscle cells is mediated by the Wnt/β-catenin pathway. JBMR Plus 1 (2), 86–100. 10.1002/jbm4.10015 Huesa C. Helfrich M. H. Aspden R. M. (2010). Parallel-plate fluid flow systems for bone cell stimulation. J. Biomech. 43 (6), 1182–1189. 10.1016/j.jbiomech.2009.11.029 Hung C. T. Pollack S. R. Reilly T. M. Brighton C. T. (1995). Real-time calcium response of cultured bone cells to fluid flow. Clin. Orthop. Relat. Res. 313, 256–269. Jacobs C. R. Yellowley C. E. Davis B. R. Zhou Z. Cimbala J. M. Donahue H. J. (1998). Differential effect of steady versus oscillating flow on bone cells. J. Biomech. 31 (11), 969–976. 10.1016/s0021-9290(98)00114-6 Jiang N. Guo W. Chen M. Zheng Y. Zhou J. Kim S. G. (2016). Periodontal ligament and alveolar bone in health and adaptation: tooth movement. Front. Oral Biol. 18, 1–8. 10.1159/000351894 Joldersma M. Burger E. H. Semeins C. M. Klein-Nulend J. (2000). Mechanical stress induces COX-2 mRNA expression in bone cells from elderly women. J. Biomech. 33 (1), 53–61. 10.1016/s0021-9290(99)00172-4 Joldersma M. Klein-Nulend J. Oleksik A. M. Heyligers I. C. Burger E. H. (2001). Estrogen enhances mechanical stress-induced prostaglandin production by bone cells from elderly women. Am. J. Physiol. Endocrinol. Metab. 280 (3), E436–E442. 10.1152/ajpendo.2001.280.3.E436 Kuo Y. C. Chang T. H. Hsu W. T. Zhou J. Lee H. H. Hui-Chun Ho J. (2015). Oscillatory shear stress mediates directional reorganization of actin cytoskeleton and alters differentiation propensity of mesenchymal stem cells. Stem Cells 33 (2), 429–442. 10.1002/stem.1860 Kanehisa M. Goto S. (2000). Kegg: kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 28 (1), 27–30. 10.1093/nar/28.1.27 Karimbux N. Y. Rosenblum N. D. Nishimura I. (1992). Site-specific expression of collagen I and XII mRNAs in the rat periodontal ligament at two developmental stages. J. Dent. Res. 71 (7), 1355–1362. 10.1177/00220345920710070201 Kato Y. Windle J. J. Koop B. A. Mundy G. R. Bonewald L. F. (1997). Establishment of an osteocyte-like cell line, MLO-Y4. J. Bone Min. Res. 12 (12), 2014–2023. 10.1359/jbmr.1997.12.12.2014 Kawada J. Komatsu K. (2000). In vitro effects of collagenase on biomechanical properties and morphological features of the rat molar periodontal ligament. Jpn. J. Oral Biol. 42 (3), 193–205. 10.2330/joralbiosci1965.42.193 Kitase Y. Lee S. Gluhak-Heinrich J. Johnson M. L. Harris S. E. Bonewald L. F. (2014). CCL7 is a protective factor secreted by mechanically loaded osteocytes. J. Dent. Res. 93 (11), 1108–1115. 10.1177/0022034514553008 Klaus B. Reisenauer S. (2016). An end to end workflow for differential gene expression using Affymetrix microarrays. F1000Res 5, 1384. 10.12688/f1000research.8967.2 Klein-Nulend J. Semeins C. M. Ajubi N. E. Nijweide P. J. Burger E. H. (1995a). Pulsating fluid flow increases nitric oxide (NO) synthesis by osteocytes but not periosteal fibroblasts--correlation with prostaglandin upregulation. Biochem. Biophys. Res. Commun. 217 (2), 640–648. 10.1006/bbrc.1995.2822 Klein-Nulend J. Semeins C. M. Burger E. H. (1996). Prostaglandin mediated modulation of transforming growth factor-beta metabolism in primary mouse osteoblastic cells in vitro . J. Cell. Physiol. 168 (1), 1–7. 10.1002/(sici)1097-4652(199607)168:1<1::Aid-jcp1>3.0.Co;2-t Klein-Nulend J. Sterck J. G. Semeins C. M. Lips P. Joldersma M. Baart J. A. (2002). Donor age and mechanosensitivity of human bone cells. Osteoporos. Int. 13 (2), 137–146. 10.1007/s001980200005 Klein-Nulend J. van der Plas A. Semeins C. M. Ajubi N. E. Frangos J. A. Nijweide P. J. (1995b). Sensitivity of osteocytes to biomechanical stress in vitro . FASEB J. 9 (5), 441–445. 10.1096/fasebj.9.5.7896017 Krishnan V. Davidovitch Z. (2009). On a path to unfolding the biological mechanisms of orthodontic tooth movement. J. Dent. Res. 88 (7), 597–608. 10.1177/0022034509338914 Kulkarni R. N. Bakker A. D. Everts V. Klein-Nulend J. (2010). Inhibition of osteoclastogenesis by mechanically loaded osteocytes: involvement of MEPE. Calcif. Tissue Int. 87 (5), 461–468. 10.1007/s00223-010-9407-7 Kulkarni R. N. Bakker A. D. Everts V. Klein-Nulend J. (2012). Mechanical loading prevents the stimulating effect of IL-1β on osteocyte-modulated osteoclastogenesis. Biochem. Biophys. Res. Commun. 420 (1), 11–16. 10.1016/j.bbrc.2012.02.099 Lau K. H. Kapur S. Kesavan C. Baylink D. J. (2006). Up-regulation of the Wnt, estrogen receptor, insulin-like growth factor-I, and bone morphogenetic protein pathways in C57BL/6J osteoblasts as opposed to C3H/HeJ osteoblasts in part contributes to the differential anabolic response to fluid shear. J. Biol. Chem. 281 (14), 9576–9588. 10.1074/jbc.M509205200 Law C. W. Alhamdoosh M. Su S. Dong X. Tian L. Smyth G. K. (2016). RNA-seq analysis is easy as 1-2-3 with limma, Glimma and edgeR. Glimma edgeR 5, 1408. 10.12688/f1000research.9005.3 Lee H. J. Diaz M. F. Ewere A. Olson S. D. Cox C. S. Jr. Wenzel P. L. (2017). Focal adhesion kinase signaling regulates anti-inflammatory function of bone marrow mesenchymal stromal cells induced by biomechanical force. Cell. Signal. 38, 1–9. 10.1016/j.cellsig.2017.06.012 Lewis K. J. Cabahug-Zuckerman P. Boorman-Padgett J. F. Basta-Pljakic J. Louie J. Stephen S. (2021). Estrogen depletion on in vivo osteocyte calcium signaling responses to mechanical loading. Bone 152, 116072. 10.1016/j.bone.2021.116072 Lewis K. J. Frikha-Benayed D. Louie J. Stephen S. Spray D. C. Thi M. M. (2017). Osteocyte calcium signals encode strain magnitude and loading frequency in vivo . Proc. Natl. Acad. Sci. U. S. A. 114 (44), 11775–11780. 10.1073/pnas.1707863114 Li J. Rose E. Frances D. Sun Y. You L. (2012). Effect of oscillating fluid flow stimulation on osteocyte mRNA expression. J. Biomech. 45 (2), 247–251. 10.1016/j.jbiomech.2011.10.037 Li M. Zhang C. Yang Y. (2019a). Effects of mechanical forces on osteogenesis and osteoclastogenesis in human periodontal ligament fibroblasts: A systematic review of in vitro studies. Bone Jt. Res. 8 (1), 19–31. 10.1302/2046-3758.81.BJR-2018-0060.R1 Li X. Han L. Nookaew I. Mannen E. Silva M. J. Almeida M. (2019b). Stimulation of Piezo1 by mechanical signals promotes bone anabolism. Elife 8, e49631. 10.7554/eLife.49631 Li X. Liu C. Li P. Li S. Zhao Z. Chen Y. (2013). Connexin 43 is a potential regulator in fluid shear stress-induced signal transduction in osteocytes. J. Orthop. Res. 31 (12), 1959–1965. 10.1002/jor.22448 Li Y. Jacox L. A. Little S. H. Ko C. C. (2018). Orthodontic tooth movement: the biology and clinical implications. Kaohsiung J. Med. Sci. 34 (4), 207–214. 10.1016/j.kjms.2018.01.007 Liao C. Cheng T. Wang S. Zhang C. Jin L. Yang Y. (2017). Shear stress inhibits IL-17A-mediated induction of osteoclastogenesis via osteocyte pathways. Bone 101, 10–20. 10.1016/j.bone.2017.04.003 Lim K. T. Hexiu J. Kim J. Seonwoo H. Choung P. H. Chung J. H. (2014). Synergistic effects of orbital shear stress on in vitro growth and osteogenic differentiation of human alveolar bone-derived mesenchymal stem cells. Biomed. Res. Int. 2014, 1–18. 10.1155/2014/316803 Lim K. T. Kim J. Seonwoo H. Chang J. U. Choi H. Hexiu J. (2013). Enhanced osteogenesis of human alveolar bone-derived mesenchymal stem cells for tooth tissue engineering using fluid shear stress in a rocking culture method. Tissue Eng. Part C Methods 19 (2), 128–145. 10.1089/ten.tec.2012.0017 Livak K. J. Schmittgen T. D. (2001). Analysis of relative gene expression data using real-time quantitative PCR and the 2−ΔΔCT method. Methods 25 (4), 402–408. 10.1006/meth.2001.1262 Love M. I. Anders S. Kim V. Huber W. (2015). RNA-seq workflow: gene-level exploratory analysis and differential expression. F1000Res 4, 1070. 10.12688/f1000research.7035.1 Love M. I. Huber W. Anders S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15 (12), 550. 10.1186/s13059-014-0550-8 Lu X. L. Huo B. Park M. Guo X. E. (2012). Calcium response in osteocytic networks under steady and oscillatory fluid flow. Bone 51 (3), 466–473. 10.1016/j.bone.2012.05.021 Macsai C. E. Foster B. K. Xian C. J. (2008). Roles of Wnt signalling in bone growth, remodelling, skeletal disorders and fracture repair. J. Cell. Physiol. 215 (3), 578–587. 10.1002/jcp.21342 Marchesan J. T. Scanlon C. S. Soehren S. Matsuo M. Kapila Y. L. (2011). Implications of cultured periodontal ligament cells for the clinical and experimental setting: A review. Arch. Oral Biol. 56 (10), 933–943. 10.1016/j.archoralbio.2011.03.003 Martens M. Ammar A. Riutta A. Waagmeester A. Slenter D. N. Hanspers K. (2021). WikiPathways: connecting communities. Nucleic Acids Res. 49 (D1), D613–D621. 10.1093/nar/gkaa1024 McAllister T. N. Frangos J. A. (1999). Steady and transient fluid shear stress stimulate NO release in osteoblasts through distinct biochemical pathways. J. Bone Min. Res. 14 (6), 930–936. 10.1359/jbmr.1999.14.6.930 McGarry J. G. Klein-Nulend J. Mullender M. G. Prendergast P. J. (2005). A comparison of strain and fluid shear stress in stimulating bone cell responses--a computational and experimental study. FASEB J. 19 (3), 1–22. 10.1096/fj.04-2210fje Mitchell M. Muftakhidinov B. Winchen T. Wilms A. van Schaik B. badshah400 (2019). Engauge digitizer software. version 12.1 ed.). Moher D. Liberati A. Tetzlaff J. Altman D. G. Group P. (2009). Preferred reporting items for systematic reviews and meta-analyses: the PRISMA statement. BMJ 339, b2535. 10.1136/bmj.b2535 Morris J. H. Apeltsin L. Newman A. M. Baumbach J. Wittkop T. Su G. (2011). clusterMaker: a multi-algorithm clustering plugin for Cytoscape. BMC Bioinforma. 12, 436. 10.1186/1471-2105-12-436 Mullender M. El Haj A. J. Yang Y. van Duin M. A. Burger E. H. Klein-Nulend J. (2004). Mechanotransduction of bone cells in vitro: mechanobiology of bone tissue. Med. Biol. Eng. Comput. 42 (1), 14–21. 10.1007/bf02351006 Nauman E. A. Satcher R. L. Keaveny T. M. Halloran B. P. Bikle D. D. (2001). Osteoblasts respond to pulsatile fluid flow with short-term increases in PGE(2) but no change in mineralization. J. Appl. Physiol. (1985) 90 (5), 1849–1854. 10.1152/jappl.2001.90.5.1849 Naveh G. R. S. Foster J. E. Silva Santisteban T. M. Yang X. Olsen B. R. (2018). Nonuniformity in ligaments is a structural strategy for optimizing functionality. Proc. Natl. Acad. Sci. U. S. A. 115 (36), 9008–9013. 10.1073/pnas.1807324115 Nogueira A. V. Nokhbehsaim M. Eick S. Bourauel C. Jager A. Jepsen S. (2014). Regulation of visfatin by microbial and biomechanical signals in PDL cells. Clin. Oral Investig. 18 (1), 171–178. 10.1007/s00784-013-0935-1 Oehmke M. J. Schramm C. R. Knolle E. Frickey N. Bernhart T. Oehmke H. J. (2004). Age-dependent changes of the periodontal ligament in rats. Microsc. Res. Tech. 63 (4), 198–202. 10.1002/jemt.20027 Okada S. Yoshida S. Ashrafi S. H. Schraufnagel D. E. (2002). The canalicular structure of compact bone in the rat at different ages. Microsc. Microanal. 8 (2), 104–115. 10.1017/s1431927601020037 Page M. J. McKenzie J. E. Bossuyt P. M. Boutron I. Hoffmann T. C. Mulrow C. D. (2021). The PRISMA 2020 statement: an updated guideline for reporting systematic reviews. BMJ 372, n71. 10.1136/bmj.n71 Pico A. R. Kelder T. van Iersel M. P. Hanspers K. Conklin B. R. Evelo C. (2008). WikiPathways: pathway editing for the people. PLoS Biol. 6 (7), e184. 10.1371/journal.pbio.0060184 Piekarski K. Munro M. (1977). Transport mechanism operating between blood supply and osteocytes in long bones. Nature 269 (5623), 80–82. 10.1038/269080a0 Ponik S. M. Triplett J. W. Pavalko F. M. (2007). Osteoblasts and osteocytes respond differently to oscillatory and unidirectional fluid flow profiles. J. Cell. Biochem. 100 (3), 794–807. 10.1002/jcb.21089 Price C. Zhou X. Li W. Wang L. (2011). Real-time measurement of solute transport within the lacunar-canalicular system of mechanically loaded bone: direct evidence for load-induced fluid flow. J. Bone Min. Res. 26 (2), 277–285. 10.1002/jbmr.211 Qi L. Zhang Y. (2014). The microRNA 132 regulates fluid shear stress-induced differentiation in periodontal ligament cells through mTOR signaling pathway. Cell. Physiol. biochem. 33 (2), 433–445. 10.1159/000358624 Qiu S. Rao D. S. Palnitkar S. Parfitt A. M. (2002). Age and distance from the surface but not menopause reduce osteocyte density in human cancellous bone. Bone 31 (2), 313–318. 10.1016/s8756-3282(02)00819-0 Rath A. L. Bonewald L. F. Ling J. Jiang J. X. Van Dyke M. E. Nicolella D. P. (2010). Correlation of cell strain in single osteocytes with intracellular calcium, but not intracellular nitric oxide, in response to fluid flow. J. Biomech. 43 (8), 1560–1564. 10.1016/j.jbiomech.2010.01.030 Reich K. M. McAllister T. N. Gudi S. Frangos J. A. (1997). Activation of G proteins mediates flow-induced prostaglandin E2 production in osteoblasts. Endocrinology 138 (3), 1014–1018. 10.1210/endo.138.3.4999 Reimand J. Isserlin R. Voisin V. Kucera M. Tannus-Lopes C. Rostamianfar A. (2019). Pathway enrichment analysis and visualization of omics data using g:Profiler, GSEA, Cytoscape and EnrichmentMap. Nat. Protoc. 14 (2), 482–517. 10.1038/s41596-018-0103-9 Reitan K. (1960). Tissue behavior during orthodontic tooth movement. Am. J. Orthod. 46 (12), 881–900. 10.1016/0002-9416(60)90091-9 Riddle R. C. Taylor A. F. Genetos D. C. Donahue H. J. (2006). MAP kinase and calcium signaling mediate fluid flow-induced human mesenchymal stem cell proliferation. Am. J. Physiol. Cell Physiol. 290 (3), C776–C784. 10.1152/ajpcell.00082.2005 Riquelme M. A. Gu S. Hua R. Jiang J. X. (2021). Mechanotransduction via the coordinated actions of integrins, PI3K signaling and Connexin hemichannels. Bone Res. 9 (1), 8. 10.1038/s41413-020-00126-w Ritchie M. E. Phipson B. Wu D. Hu Y. Law C. W. Shi W. (2015). Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43 (7), e47. 10.1093/nar/gkv007 Robling A. G. Bonewald L. F. (2020). The osteocyte: new insights. Annu. Rev. Physiol. 82, 485–506. 10.1146/annurev-physiol-021119-034332 Rodriguez-Florez N. Oyen M. L. Shefelbine S. J. (2014). Age-related changes in mouse bone permeability. J. Biomech. 47 (5), 1110–1116. 10.1016/j.jbiomech.2013.12.020 Salvi J. D. Lim J. Y. Donahue H. J. (2010). Increased mechanosensitivity of cells cultured on nanotopographies. J. Biomech. 43 (15), 3058–3062. 10.1016/j.jbiomech.2010.07.015 Schardt C. Adams M. B. Owens T. Keitz S. Fontelo P. (2007). Utilization of the PICO framework to improve searching PubMed for clinical questions. BMC Med. Inf. Decis. Mak. 7, 16. 10.1186/1472-6947-7-16 Schurman C. A. Verbruggen S. W. Alliston T. (2021). Disrupted osteocyte connectivity and pericellular fluid flow in bone with aging and defective TGF-β signaling. Proc. Natl. Acad. Sci. U. S. A. 118 (25), e2023999118. 10.1073/pnas.2023999118 Seref-Ferlengez Z. Maung S. Schaffler M. B. Spray D. C. Suadicani S. O. Thi M. M. (2016). P2X7R-Panx1 complex impairs bone mechanosignaling under high glucose levels associated with type-1 diabetes. PLoS One 11 (5), e0155107. 10.1371/journal.pone.0155107 Shannon P. Markiel A. Ozier O. Baliga N. S. Wang J. T. Ramage D. (2003). Cytoscape: A software environment for integrated models of biomolecular interaction networks. Genome Res. 13 (11), 2498–2504. 10.1101/gr.1239303 Sonam S. Sathe S. R. Yim E. K. Sheetz M. P. Lim C. T. (2016). Cell contractility arising from topography and shear flow determines human mesenchymal stem cell fate. Sci. Rep. 6, 20415. 10.1038/srep20415 Street J. Bao M. deGuzman L. Bunting S. Peale F. V. Jr. Ferrara N. (2002). Vascular endothelial growth factor stimulates bone repair by promoting angiogenesis and bone turnover. Proc. Natl. Acad. Sci. U. S. A. 99 (15), 9656–9661. 10.1073/pnas.152324099 Su G. Kuchinsky A. Morris J. H. States D. J. Meng F. (2010). GLay: community structure analysis of biological networks. Bioinformatics 26 (24), 3135–3137. 10.1093/bioinformatics/btq596 Su G. Morris J. H. Demchak B. Bader G. D. (2014). Biological network exploration with Cytoscape 3. Curr. Protoc. Bioinforma. 47, 8 13 11–24. 10.1002/0471250953.bi0813s47 Sun C. Janjic Rankovic M. Folwaczny M. Otto S. Wichelhaus A. Baumert U. (2021). Effect of tension on human periodontal ligament cells: systematic review and network analysis. Front. Bioeng. Biotechnol. 9, 695053. 10.3389/fbioe.2021.695053 Szklarczyk D. Gable A. L. Lyon D. Junge A. Wyder S. Huerta-Cepas J. (2019). STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res. 47 (D1), D607–D613. 10.1093/nar/gky1131 Szklarczyk D. Kirsch R. Koutrouli M. Nastou K. Mehryary F. Hachilif R. (2023). The STRING database in 2023: protein-protein association networks and functional enrichment analyses for any sequenced genome of interest. Nucleic Acids Res. 51 (D1), D638–D646. 10.1093/nar/gkac1000 Tada S. Tarbell J. M. (2000). Interstitial flow through the internal elastic lamina affects shear stress on arterial smooth muscle cells. Am. J. Physiol. Heart Circ. Physiol. 278 (5), H1589–H1597. 10.1152/ajpheart.2000.278.5.H1589 Tang M. Peng Z. Mai Z. Chen L. Mao Q. Chen Z. (2014). Fluid shear stress stimulates osteogenic differentiation of human periodontal ligament cells via the extracellular signal-regulated kinase 1/2 and p38 mitogen-activated protein kinase signaling pathways. J. Periodontol. 85 (12), 1806–1813. 10.1902/jop.2014.140244 Turner C. H. Robling A. G. Duncan R. L. Burr D. B. (2002). Do bone cells behave like a neuronal network? Calcif. Tissue Int. 70 (6), 435–442. 10.1007/s00223-001-1024-z van der Meijden K. Bakker A. D. van Essen H. W. Heijboer A. C. Schulten E. A. Lips P. (2016). Mechanical loading and the synthesis of 1,25(OH)2D in primary human osteoblasts. J. Steroid Biochem. Mol. Biol. 156, 32–39. 10.1016/j.jsbmb.2015.11.014 van der Pauw M. T. Klein-Nulend J. van den Bos T. Burger E. H. Everts V. Beertsen W. (2000). Response of periodontal ligament fibroblasts and gingival fibroblasts to pulsating fluid flow: nitric oxide and prostaglandin E2 release and expression of tissue non-specific alkaline phosphatase activity. J. Periodontal Res. 35 (6), 335–343. 10.1034/j.1600-0765.2000.035006335.x Verbruggen S. W. Mc Garrigle M. J. Haugh M. G. Voisin M. C. McNamara L. M. (2015). Altered mechanical environment of bone cells in an animal model of short- and long-term osteoporosis. Biophys. J. 108 (7), 1587–1598. 10.1016/j.bpj.2015.02.031 Verbruggen S. W. Vaughan T. J. McNamara L. M. (2014). Fluid flow in the osteocyte mechanical environment: A fluid-structure interaction approach. Biomech. Model. Mechanobiol. 13 (1), 85–97. 10.1007/s10237-013-0487-y Wang Y. McNamara L. M. Schaffler M. B. Weinbaum S. (2007). A model for the role of integrins in flow induced mechanotransduction in osteocytes. Proc. Natl. Acad. Sci. U. S. A. 104 (40), 15941–15946. 10.1073/pnas.0707246104 Weinbaum S. Cowin S. C. Zeng Y. (1994). A model for the excitation of osteocytes by mechanical loading-induced bone fluid shear stresses. J. Biomech. 27 (3), 339–360. 10.1016/0021-9290(94)90010-8 Wilks C. Zheng S. C. Chen F. Y. Charles R. Solomon B. Ling J. P. (2021). recount3: summaries and queries for large-scale RNA-seq expression and splicing. Genome Biol. 22 (1), 323. 10.1186/s13059-021-02533-6 Wittkowske C. Reilly G. C. Lacroix D. Perrault C. M. (2016). In vitro bone cell models: impact of fluid shear stress on bone formation. Front. Bioeng. Biotechnol. 4, 87. 10.3389/fbioe.2016.00087 Xia X. Batra N. Shi Q. Bonewald L. F. Sprague E. Jiang J. X. (2010). Prostaglandin promotion of osteocyte gap junction function through transcriptional regulation of connexin 43 by glycogen synthase kinase 3/β-Catenin signaling. Mol. Cell. Biol. 30 (1), 206–219. 10.1128/mcb.01844-08 Yan Z. Wang P. Wu J. Feng X. Cai J. Zhai M. (2018). Fluid shear stress improves morphology, cytoskeleton architecture, viability, and regulates cytokine expression in a time-dependent manner in MLO-Y4 cells. Cell Biol. Int. 42 (10), 1410–1422. 10.1002/cbin.11032 You J. Reilly G. C. Zhen X. Yellowley C. E. Chen Q. Donahue H. J. (2001). Osteopontin gene regulation by oscillatory fluid flow via intracellular calcium mobilization and activation of mitogen-activated protein kinase in MC3T3-E1 osteoblasts. J. Biol. Chem. 276 (16), 13365–13371. 10.1074/jbc.M009846200 You L. Temiyasathit S. Lee P. Kim C. H. Tummala P. Yao W. (2008). Osteocytes as mechanosensors in the inhibition of bone resorption due to mechanical loading. Bone 42 (1), 172–179. 10.1016/j.bone.2007.09.047 Zeng Y. Cowin S. C. Weinbaum S. (1994). A fiber matrix model for fluid flow and streaming potentials in the canaliculi of an osteon. Ann. Biomed. Eng. 22 (3), 280–292. 10.1007/BF02368235 Zhang C. Bakker A. D. Klein-Nulend J. Bravenboer N. (2019). Studies on osteocytes in their 3D native matrix versus 2D in vitro models. Curr. Osteoporos. Rep. 17 (4), 207–216. 10.1007/s11914-019-00521-1 Zhang K. Barragan-Adjemian C. Ye L. Kotha S. Dallas M. Lu Y. (2006). E11/gp38 selective expression in osteocytes: regulation by mechanical strain and role in dendrite elongation. Mol. Cell. Biol. 26 (12), 4539–4552. 10.1128/mcb.02120-05 Zhou X. Liu D. You L. Wang L. (2010). Quantifying fluid shear stress in a rocking culture dish. J. Biomech. 43 (8), 1598–1602. 10.1016/j.jbiomech.2009.12.028 Glossary
AKT1 AKT serine/threonine kinase 1
ALP/ALPL/ALPP Alkaline phosphatase
BGLAP Bone gamma-carboxyglutamate protein
BMP2 Bone morphogenetic protein 2
Ca2+ Calcium ion
COL1A1 Collagen type 1 alpha 1
CTNNB1 Catenin beta 1
DEG/DEGs Differential expressed gene/genes
ELISA Enzyme-linked Immunosorbent Assay
FEA Finite element analysis
FSS Fluid flow shear stress
GEO Gene Expression Omnibus
GJA1 Gap junction protein alpha 1
GO GeneOntology
GSEA Gene Set Enrichment Analysis
hMSCs Human mesenchymal stem cells
hOcyt Human osteocytes
hOst Human osteoblasts
hPDLCs Human periodontal ligament cells
IF mmunofluorescence
KEGG Kyoto Encyclopedia of Genes and Genomes
MAPK1 Mitogen-Activated Protein Kinase 1
MAPK3 Mitogen-Activated Protein Kinase 3
mOcyt Mouse osteocytes
mOst Mouse osteoblasts
NO Nitric oxide
ORA Over Representation Analysis
OTM Orthodontic tooth movement
PDL Periodontal ligament
PDLCs Periodontal ligament cells
PDPN Podoplanin
PGE-2 Prostaglandin E2
PTGS1 Prostaglandin-endoperoxide synthase 1 (aka cyclooxygenase-1, Cox-1)
PTGS2 Prostaglandin-endoperoxide synthase 2 (aka cyclooxygenase-2, Cox-2)
RoB Risk of bias
RT-qPCR Reverse-transcriptase quantitative real-time polymerase chain reaction
RUNX2 Runt-related transcription factor 2
SPP1 Secreted phosphoprotein 1 (aka osteopontin, OPN)
sqPCR Semiquantitative reverse-transcriptase polymerase chain reaction
STRING Search Tool for Retrieval of Interacting Genes/Proteins
TNFRSF11B TNF receptor superfamily member 11b (aka osteoprotegerin, OPG)
TNFSF11 TNF superfamily member 11 (aka Receptor Activator of NF-κB Ligand; RANKL)
VEGFA Vascular endothelial growth factor A
WB Western blotting
WNT3A Wnt family member 3A