The small world of osteocytes: connectomics of the lacuno-canalicular network in bone

Network structures are ubiquitous in Nature and often fulfill important functions in transport and signal processing. In humans, examples include the airways (Metzger et al 2008), the vasculature (Blinder et al 2013, Pries and Secomb 2014) and the nervous system (Eguiluz et al 2005, Bullmore and Sporns 2009). The organization into networks can already occur on the level of the individual cells, with neurons using their extended dendrites and synapses to connect with other cells as the prototypical example (Cajal and May 1928, Helmstaedter et al 2013, Takemura et al 2013). It has also been known for a long time that the mineralized bone tissue of most vertebrates is densely populated by cells called osteocytes. These cells are embedded in the extracellular bone matrix during bone deposition, and are linked with each other and with blood vessels through a highly interconnected network of cell processes (Bonewald 2011) in appearance and size similar to the neuronal system (Buenzli and Sims 2015). Osteocyte bodies and their processes reside in hollow lacunae and narrow canals termed canaliculi, respectively, that together comprise the osteocyte lacuno-canalicular network (OLCN) (Franz-Odendaal et al 2006). The functional relevance of this network and the details of its architecture are still under debate, despite a large number of recent studies in this field (Asada et al 2013, Thi et al 2013, Hesse et al 2015, Milovanovic et al 2015, Sano et al 2015, Nango et al 2016).

One assumed function of osteocytes is to orchestrate the process of structural adaptation and material renewal in bone. These processes are thought to be mechanically controlled, with the result that new bone is added where mechanically needed and resorbed where local loading is low. The role of the osteocytes is to locally initiate bone (re)modeling by generating biochemical signals that activate precursor cells in response to mechanical stimuli (Nakashima et al 2011). In this case, the role of the OLCN would be that of a signal amplifier: small deformations or cracks within the bone tissue are sensed by osteocyte processes and transduced into a biochemical response that is then spatially integrated and triggers the recruitment of osteoblast and osteoclast precursor cells (Wang et al 2007). Besides this putative mechano-sensory role, osteocytes are suspected to be able to locally remodel the bone matrix that surrounds them (osteocytic osteolysis) (Belanger 1969, Teti and Zallone 2009). This would enable them to participate in the calcium and phosphate metabolism of the body by making the mineral reservoir deep in the bone tissue accessible through the network of canaliculi. The ability of osteocytes to secrete endocrine factors involved in mineral homeostasis strongly supports this possibility (Dallas et al 2013). For both, mechanosensing and mineral exchange, the efficiency of the network in distributing molecules and signals between cells and throughout the matrix plays a central role, as suggested by in vivo tracer studies and theoretical modeling (Tate et al 1998). Finally, the function of the canalicular network for the osteocytes is to supply nutrients and remove waste by maintaining contact to blood vessels (Bonewald 2011).

The efficiency of real-world networks such as roadmaps, metabolic networks, or the internet can be assessed using methods and measures derived from graph theory (Boccaletti et al 2006). By representing these diverse systems as a collection of edges connected by nodes, their structural and topological properties can be analyzed and compared. Common features of real-world networks include exponential or scale-free degree distributions, high clustering coefficients (CC), and hierarchical organization. In contrast to random graphs, empirical networks often exhibit small-world properties such that most points in the network can be reached from anywhere by a small number of steps, which is very efficient for transport and communication. The application of complex network theory to biological networks resulted in many important insights into the structure and function of these networks (Bullmore and Sporns 2009, Pries and Secomb 2014).

To date, such a characterization could not be carried out for the osteocyte network in bone, in part due to a lack of structural data that provide sufficient resolution and at the same time span large enough volumes. Direct imaging of osteocytes and their processes in living bone is difficult due to the limited accessibility and transparency of the bone matrix; therefore, most techniques for imaging the osteocyte network rely on imaging the OLCN in bone samples after cells and soft tissue have been removed. Mainly, three classes of imaging modalities have been reported, based on x-ray micro computed tomography (μCT), scanning electron microscopy (SEM), or confocal laser scanning microscopy (CLSM). TEM and TEM-based electron tomography are not well suited for imaging larger parts of the OLCN due to limitations in sample thickness and preparation. While conventional μCT scanners can image individual lacunae including their shape and orientation, phase-resolved synchrotron x-ray tomography has recently been applied to derive 3D images of the network of canaliculi around a single osteocyte at a few nm resolution (Dierolf et al 2010, Hesse et al 2015). SEM can only image surfaces due to the low penetration depth of electrons in dense tissue, but two approaches have been reported to derive 3D data of the OLCN using SEM. In one approach, the hollow spaces in bone are filled with a resin, and the surrounding bone tissue is then etched away, exposing the resin cast and hence the structure of the OLCN as a projection (Pazzaglia and Congiu 2013). The other approach is FIB-SEM slice-and-view, where the surface is repeatedly imaged with SEM and milled using a focused ion beam. This variant produced the most detailed and accurate 3D images of the OLCN (Schneider et al 2011) and its embedment into the collagen matrix (Reznikov et al 2014) to date, yet both SEM-based techniques inevitably result in the destruction of the sample. The third method, CLSM, requires fluorescent staining of the network while the sample itself is preserved (Anderson et al 2008, Kerschnitzki et al 2011). The resolution of CLSM is limited by the optical diffraction limit and is intrinsically anisotropic and depth-dependent, and therefore, in contrast to XRPT and FIB-SEM, cannot accurately resolve the diameter of individual canaliculi. Scaled-up models from CLSM stacks have e.g. been used to study transport properties of the OLCN, without detailed quantification of the canalicular network (Anderson et al 2008).

We previously showed that by image processing, the network can be reconstructed from CLSM image stacks as long as the resolution is sufficient to separate neighboring canaliculi (Kerschnitzki et al 2013). This study revealed a relationship between the proximity of the bone matrix to canaliculi and the nanoscale properties of the mineral crystals in the matrix, suggesting active participation of osteocytes in bone remodeling and mineral homeostasis. Here, we present an improvement of this imaging and analysis technique and use it to quantify the architecture of the osteocyte network in different bone types from mouse and sheep. We use fibrolamellar bone from sheep, which grows more slowly but in a more controlled way on top of a template surface (Ferretti et al 2002, Liu et al 2010, Kerschnitzki et al 2011), and woven bone from mice, which grows faster but in the absence of an organizing surface. Based on these data, we report the first detailed characterization of its topological and statistical properties on the subcellular level using measures from complex network theory. The OLCN reconstructed from CLSM image stacks is described as a network graph, i.e. by the locations of nodes and edges, and the adjacency matrix, which contains information on how nodes are connected through edges. This makes it possible to calculate generic statistical and topological properties such as edge length and degree distributions, and relationships between these properties. With this proof-of-principle approach, we aim to address the question whether the known differences in tissue organization and appearance of the OLCN in different types of bone from different species are reflected in the organization of the network of cell processes, or if the network at the cellular and subcellular level is independent of bone type and level of organization. By developing a quantitative characterization of the topological and statistical properties of the OLCN in different bone types, we hope to learn more about what purpose this network may have evolved for, how its formation is controlled, and how these network measures relate to bone material quality during physiological development or pathological conditions of age, disease and pharmaceutical intervention.

2.1. Sample preparation and imaging

The workflow from sample preparation to image analysis and quantification is summarized in figure 1(a). Bone samples, sample preparation, and staining with Rhodamine-6G were the same as previously described (Kerschnitzki et al 2011, 2013). We used fibrolamellar ovine bone from the mid-diaphysis of the femur of a 5 year old sheep as well as murine bone from the mid-diaphysis of the femur from a 12 month-old mouse, without initial aldehyde fixation. Bone samples were cut either in cross sections or in longitudinal sections with an initial thickness of 200 μm and polished from both sides in an automatic polisher (Logitech PM5, Logitech Ltd, Glasglow, UK) until a final thickness of 80 μm. Samples were kept wet throughout the whole preparation process. Until further usage, samples were wrapped in cotton gauze soaked in phosphate buffered saline solution (Sigma Aldrich GmbH, Germany) and stored at 4 °C. For each bone type (murine woven bone and ovine fibrolamellar bone) and orientation (longitudinal and transversal), five different non-overlapping fields of view far from the sample border were randomly selected and imaged, amounting to 20 image stacks with 51 images per stack in total. Image acquisition was done on a Leica DM IRBE confocal microscope (Leica Microsystems, Wetzlar) using a 100x oil immersion objective (NA = 1.4) using a voxel size of ∼0.2 μm × 0.2 μm × 0.2 μm. Total volume size was 512 × 512 × 51 voxels or 100 μm × 100 μm × 10 μm.

Zoom In Zoom Out Reset image size

2.2. Image processing

2.3. Cell lacunae detection and skeletonization

2.4. Spatial network analysis

The total void fraction is defined as the number of all foreground (lacunar or canalicular) voxels in the binary image divided by the total number of voxels. Canalicular density (Ca.Dn) is a measure for the density of cell processes within the matrix and is defined as the number of network voxels without lacunae divided by the total number of voxels without lacunae (Repp et al 2017). The proximity of a matrix voxel to cell lacunae, d_L, and to the entire lacuno-canalicular network d_LC, was measured by calculating the 3D Euclidean distance transform of the background to the foreground voxels in the binarized images of lacunae only, or of the entire network including lacunae. The distance d_net of any point within the network of canaliculi to reach the closest cell through the network was determined by iteratively tagging all canalicular voxels starting from cell lacunae. The (dimensionless) gain factor G_net of transport efficiency between cells and matrix due to presence of the canalicular network was calculated as the actual distance to the next lacuna through the network divided by the effective distance if the transport through canaliculi is k = v_network/v_matrix times faster compared to passive diffusion through the bone matrix,

$$\begin{eqnarray*}&&{G}_{{\rm{n}}{\rm{e}}{\rm{t}}}=({d}_{{\rm{L}}{\rm{C}}}+{d}_{{\rm{n}}{\rm{e}}{\rm{t}}})/({d}_{{\rm{L}}{\rm{C}}}+1/k\,\times \,{d}_{{\rm{n}}{\rm{e}}{\rm{t}}}),\end{eqnarray*}$$

averaged over all matrix voxels. Here, we define transport efficiency in the sense of how long signals and minerals travel between cells and matrix, assuming that travel time is proportional to distance. Since d_LC is the proximity of matrix to the lacuna-canalicular network and d_net the distance to the closest lacuna through the network, the sum of both is proportional to the travel time of a molecule from any point in the matrix to the next lacuna (figure 2(d)). G_net is a lower boundary for the actual time gain, since the scaling of transport time with distance is expected to be faster for directed transport through the network compared to passive diffusion through the matrix.

Zoom In Zoom Out Reset image size

2.5. Topological network analysis

2.6. Statistics

3.1. Structural quantification

3.2. Statistical properties

3.2.1. Network size and edge density

The representation of the OLCN as a graph enabled us to apply complex network analysis to determine statistical properties of the network, and to compare them between different bone types and with other biological networks. The most general properties of a graph are its size, defined as the number of nodes N, and its edge density, which is the fraction of existing edges over all possible edges. The size of the networks was larger in the more ordered fibrolamellar bone compared to woven bone (figure 3(a)). There was a clear difference between the two sample groups within both fibrolamellar and woven bone. This correlates with the apparent network density in the respective images (figure 1(c)). Since all samples had identical volume, differences in graph size N represent differences in the number of nodes per volume. The average node density was 4 × 10⁷ per mm³ in woven bone and 7 × 10⁷ per mm³ in fibrolamellar bone.

Zoom In Zoom Out Reset image size

Remarkably, the network graphs for the different samples do not exhibit a size-independent edge density, as would be the case for a random graph (figure 3(a)). Rather, the universal invariant in the osteocyte network is the number of edges per node. When number of nodes (excluding cells and endpoints) and edges are plotted against each other, all examined samples regardless of species and bone type fall onto a single straight line, as shown in the inset to figure 3(a). The number of edges per node averaged over all samples is 3.28 ± 0.01 if only nodes of degree ≥ 3 are included, and there is no statistically significant difference between the four groups. This predicts an edge density that scales with the number of nodes as 3.28 × N⁻¹. This power law precisely fits the data in figure 3(a). The equivalent Euler characteristic for these networks is 3.28 × N. If also endpoints, i.e. nodes with degree 1, are included the average degree drops to 2.8 and becomes significantly different between sample groups due to the different percentage of endpoints (see below).

3.2.2. Universal edge length and degree distributions

We next determined the cumulative statistical distributions of edge length and node degree, which are two common measures to classify and compare complex networks. Cumulative edge lengths were exponentially distributed in all examined networks (figure 3(b)) with a typical decay parameter of −2/3. The decay parameter (or slope of the distribution in a semi-logarithmic survival plot, shown as a gray line in figure 3(b)) was not significantly different between any of the four sample groups.

Cumulative degree distributions for all four sample groups are shown in figure 3(c). All four groups followed an exponential distribution with a decay (or slope in a semi-logarithmic plot, indicated by the gray line) around −4/3. Again, there was no significant difference in the decay parameter between the four groups. Taken together, all examined networks, despite their differences in appearance and density, universally follow the same exponential edge length and degree distributions.

In a spatial network, the 'importance' of a node depends not only on the number of edges that connect to it, but also on the length, or 'weight', of the edges. To quantify this, the 'weighted degree' WD is calculated as the sum of the length of all edges connected to a node. The cumulative weighted degree distributions are shown in figure 3(d). In this evaluation, cells were included and plotted together with all other nodes in one graph. The weighted degree distribution follows an exponential decay similar to the degree distribution, but with a long tail above about 200 μm due to the presence of cells. For comparison, the gray lines indicate an exponential decay of −1/2 and a power law decay with an exponent of 3/2, respectively.

3.3. Topology and connectivity

3.3.1. Tree-like local topology

We next looked at the local topological organization of the network of interconnected canaliculi. Here, it is not only interesting to look at the number of neighboring nodes to which a node is connected (its degree), but also at the probability of connections between these neighboring nodes, expressed in the CC. A CC of one means that all neighbors of a node are connected to each other, while a CC of zero means that there are no connections between neighbors of a node. The first case corresponds to a fully connected dense region of the network, whereas the latter case indicates a tree-like organization of this region of the network (figure 4(a)).

Zoom In Zoom Out Reset image size

To characterize the local organization of the canalicular network, we classified nodes according to their CC and their degree: nodes with a CC larger than 0.5 (independent of their degree) were identified as 'cluster nodes', whereas nodes with degree = 3 and CC = 0 were defined as 'tree nodes', indicating a dendritic-like organization of the network where a cell process branches into two sub-processes (figures 4(a) and (b)). We found in all cases that the majority of nodes belonged to the category of tree nodes, suggesting that the network of cell processes is largely organized as a dendritic branching network. All investigated networks exhibited only a small number of cluster nodes. Woven bone samples had a higher percentage of cluster nodes (1.48 ± 0.5%) and a smaller percentage of tree nodes (43.0 ± 3.9%) compared to fibrolamellar bone (0.76 ± 0.2% and 53.5 ± 1.7%, respectively).

A third relevant class of vertices are the end points of branches, defined as nodes with degree = 1 that are not at the boundary of the volume (i.e. no cut-off canaliculi). The percentage of such nodes in the network is a measure for how often canaliculi end as 'one-way roads' rather than being connected to other parts of the network. We found that on average, 29.8 ± 4.2% of all nodes in woven bone are such end points, whereas only 17.5 ± 1.9% of all nodes in fibrolamellar bone are end points (figure 4(c)).

3.3.2. ASP and betweenness centrality

Our investigation of the local topology showed that the network of osteocyte cell processes can be described as a dendritic branching network interspersed with a few high-CC nodes with many connections between neighbors. Woven bone had a higher percentage of highly clustered nodes compared to fibrolamellar bone. We next asked if these local differences in topology coincide with differences in global measures for topology and transport efficiency (e.g. for nutrients or chemical and mechanical signals). One possibility to express the efficiency of transport across a network is the ASP from one node to any other in the network. The ASP averaged over all samples was 48.2 ± 3.7 μm, with no significant differences between sample groups despite the larger overall size of the network in fibrolamellar bone. For any given node, the betweenness centrality, defined as the fraction of shortest paths in the network running through that node, is a measure for the importance of that node. Figure 4(d) shows two color-coded maps of the betweenness centrality of nodes for fibrolamellar and woven bone, respectively. It is immediately clear that in both bone types, high-centrality nodes align along preferential 'highways' across the network through which many shortest paths are running.

3.3.3. Small-world properties

We next investigated if the osteocyte lacuna-canalicular system is a small world network. Small world networks are a class of networks where the average path length between any two nodes is much smaller than in a random network of same size and density, resulting in a more efficient communication (Amaral et al 2000, Boccaletti et al 2006). Many biological, technical and social networks have been found to belong to this class (Amaral et al 2000, Eguiluz et al 2005). Although there is no standard definition of a small world network, a number of practical measures for small-worldness have been proposed. The classical measure is that the ASP is much smaller than in a random network of same size, but the CC is similar (Amaral et al 2000, Boccaletti et al 2006). While the ASP was not significantly different between all subgroups, CC was slightly higher in woven bone than in fibrolamellar bone (p < 0.05, figure 4(g)). A more rigid definition proposed by Humphrey et al (Humphries and Gurney 2008) is small-worldness S, defined as the ratio of ASP and CC divided by the same ratio for a ER random graph of the same size and edge density. We calculated this ratio for each of the 20 networks and found S to be 31.9 ± 4.6 in woven bone and 44.8 ± 5.0 for fibrolamellar bone (figure 4(e)). We observe a higher S for the larger fibrolamellar bone networks and a roughly linear relationship between N and S (figure 4(e)). We compared the network size and small-worldness S to those of >30 real small-world networks from (Humphries and Gurney 2008), including the C.Elegans neuronal connectome, electric circuit maps, and train networks (figure 4(f)). We find that the osteocyte network, with its intermediate size compared to the other networks, fits well into the linear relationship between S and N.

In this study, we report an extensive quantification of the spatial and topological properties of the OLCN in bones of different structural organization, i.e. woven bone from mouse and fibrolamellar bone from sheep. Woven bone has a higher overall void fraction but a smaller density of canaliculi compared to fibrolamellar bone (figure 1(c)). Concurrently, its matrix is on average further from the next canalicular surface (figure 2(a)), whereas the distance to the next lacuna both through the matrix as well as through the network is shorter (figure 2(c)). In woven bone, these mean distances exhibit more pronounced variation between samples compared to fibrolamellar bone (figure 2(b)). A smaller variation between samples suggests that this parameter is more controlled during development. Hence, the spacing and organization of the canaliculi appears to be more tightly regulated in fibrolamellar bone compared to woven bone. The higher degree of organization of the network in fibrolamellar bone is furthermore supported by the 2.5-fold higher gain in transport efficiency due to the presence of the network (figure 2(d)). Taken together, the spatial organization of the network in fibrolamellar bone is more efficient compared to woven bone in terms of controlling the bone matrix and accessing the stored mineral.

While the quantification of the spatial architecture of the OLCN revealed pronounced differences between all four sample groups, they exhibit surprisingly universal statistical properties on the level of network topology. The number of edges per node is 3.3 in all examined networks independent of size. Both edge length and node degree are exponentially distributed, with no statistically significant differences between the four sample groups. The presence of power-law tails of the weighted degree distribution including cell lacunae emphasizes the role of cells as 'hubs' within the network.

Investigation of the local topology reveals that the network is mainly organized in a tree-like fashion, with only a small percentage of nodes with high CC. This tree-like organization is the expected outcome of a dendritic outgrowth process with successively branching cell processes (Bonewald 2005, Buenzli and Sims 2015). The exponential degree distribution and presence of high-degree nodes in the network does not contradict such a branching process as the main organizing mechanism during growth, since high-degree nodes could also represent dense clusters of branching points (tree nodes) within a small volume below the optical imaging resolution (Wittig et al 2016). This would not affect small-world properties of the network since no links are added or removed. All examined networks seem to exhibit similar local organization despite pronounced differences in the total network size.

The topological architecture of the network allows quantifying the efficiency of the network for information transfer. An important measure here is the ASP between any two nodes in the network, and the number of such paths that run through a given node (betweenness centrality). We found that the shortest paths between nodes in the network run along a few 'highways' on which most high-centrality nodes are clustered together (figure 4(c)). In a completely homogeneous spatial network, the betweenness centrality would be highest in the center and decrease towards the edge (Barthelemy 2011). Our results show no such simple relationship, but instead the clustering along important paths linking distant regions of the network reveals the small-world character of the OLCN. By comparing the ratio of ASP and CC to that of a random graph of the same size and density, we derive the small-worldness S and find it to be significantly higher in fibrolamellar compared to woven bone networks, and proportional to network size N. Proportionality between S and N has been reported for a large number of real-world networks as diverse as neuronal connectomes, transport and signaling networks with S ∼ 0.023*N^0.96 (Humphries and Gurney 2008). Despite a different proportionality observed for our data, this probably explains the observed difference in S for the different bone types. When we compare the average S and N for all samples to the data in Humphries and Gurney (2008), we find that the osteocyte network exhibits a similar S as other real-world complex networks of similar size (figure 4(f)). These include both spatial as well as non-spatial networks, despite differences such as the cost associated with link length in spatially embedded networks such as the OLCN (Barthelemy 2011). Since small-worldness is a measure for how efficiently the network is organized to transmit signals between distant points in the network, we conclude that the small-world like topology of the osteocyte network, with interconnected 'trees' emanating from lacunae, is a result of adaptation towards efficient organization. A completely random or homogeneous spatial network of same size and edge density would be missing these short connections between distant regions.

How are these findings linked to the functions of the osteocyte network? As outlined in the introduction, the OLCN is believed to be important for bone mechanosensing, signal transduction, and nutrient supply of osteocytes. In all three of these functional scenarios, the efficiency of the network in distributing molecules and signals between cells and throughout the matrix plays a central role. The optimal strategy would then be to access as much bone matrix volume as possible with as few cells as possible, while keeping the effective distance between cells through the network as short as possible. Our results show that fibrolamellar bone, despite containing fewer cells per volume, exhibits a more optimized organization of the canalicular network with respect to these functional criteria as compared to woven bone.

This advantage might be a consequence of the more controlled but slower growth of highly organized bone against a pre-existing endogenous scaffold, whereas woven bone is rapidly formed in a less controlled manner and in the absence of an organizing surface (Ferretti et al 2002, Kerschnitzki et al 2011). The guidance by the scaffold surface results in long-ranged alignment of cells and extracellular matrix fibers during bone deposition. This oriented organization does not only improve the mechanical performance of bone, but serves as a blueprint for the osteocyte network. Cell processes align with the topographical cues of the surrounding matrix (Dunn and Heath 1976, Curtis and Wilkinson 1997) and direct the orientation of subsequently deposited extracellular matrix fibers (Wang et al 2003, Lamers et al 2010). This continuous feedback between cells and their environment gradually produces more coordinated and integrated higher-level structures. Together with the longer time available to the fibrolamellar network to grow and reorganize before the tissue structure is fixed due to mineralization, this allows for a more dense, well-organized and therefore efficient network. The higher degree of organization of the fibrolamellar OLCN together with the superior mechanical performance justifies why woven bone is gradually replaced by fibrolamellar bone, despite the high cost of resorbing and depositing bone tissue. Once the network is in place, the newly remodeled bone can more efficiently contribute to the existing mechano-sensory network and mineral depot of the skeleton.

Our results show that the organization of the osteocyte network can serve as a measure for bone quality (Seeman and Delmas 2006) due to its direct functional relevance. The quantification developed here may be useful in assessing bone quality during physiological development or pathological conditions of age, disease and pharmaceutical intervention, complementary to existing parameters such as bone mineral density. Although we did not apply our analysis to compare healthy to diseased bone, our choice of different bone types reflecting different degrees of organization demonstrates the potential of our method to quantify differences in efficiency. Pathological deviations in matrix protein structure, matrix deposition, or bone homeostasis likely result in a less efficiently organized network. Conversely, deficits in network organization lead to less efficient mineral exchange, impaired mechanotransduction and overall reduced osteocyte viability due to limited nutrient supply. By quantifying the parameters of the OLCN as introduced here, such two-way interactions between network organization and bone quality can be easily assessed by comparing healthy and diseased bone samples, and the results can e.g. be used to add a more accurate description of the canalicular network to theoretical flow models (Mishra and Tate 2003).

While the density and spatial organization of the network differ between bone types, other properties of the network, such as edge length distribution, degree distribution, and the tree-like topology are surprisingly universal despite the different growth modalities. This shows that the same cellular mechanisms is at work during the outgrowth of cell processes and the formation and pruning of cell–cell contacts, independent of bone type and species. It also suggests that these mechanisms operate in a narrow parameter range that results in an optimized topology with regards to node density and edge length. For example, a certain fraction of high-degree nodes, or closely spaced 3-nodes, might be ideal for sensing damage or integrating signals from remote branches of the network.

One fundamental challenge in biology is to understand how individual cells come together and form higher-order structures that perform complex functions. Multicellular networks such as the OLCN are a common example of such emergent self-organization. Complex network theory offers a well-established framework to robustly quantify and compare such structures using measures from statistical physics and graph theory. Our results provide insight into which principles might be at work when individual bone cells self-organize into a large, interconnected network during bone formation and mineralization. In future work, the definition of small-worldness S via generic random networks (Humphries and Gurney 2008) could be extended to spatially explicit 3D models of network formation to learn more about how different growth mechanisms produce different spatial network topologies (Barthelemy 2011). The workflow that was used in this proof-of-concept study on a limited number of samples (volume imaging followed by segmentation and conversion into a graph) is broadly applicable to many biological network-like structures at different length scales. One remarkable outcome of our analysis is that the osteocyte network shares many characteristics with other multicellular networks, such as neuronal networks, and with higher-order network structures such as the vascular system or the airways (Eguiluz et al 2005, Bullmore and Sporns 2009, Blinder et al 2013, Herriges and Morrisey 2014). The analogy between the OLCN and neuronal networks has previously inspired speculations about the potential signal processing power of the whole skeleton (Turner et al 2002). The approach presented here would allow for a more in-depth quantitative comparative analysis between different networks, potentially revealing unexpected similarities and universal ordering principles.

We developed a strategy for the quantification of the architecture of the osteocyte network using different bone types from mouse and sheep for the proof of principle. Based on these data, we reported the first detailed quantitative characterization of its topological and statistical properties on the cellular level. We defined a number of robust, quantitative measures that are derived from the theory of complex networks, and used these measures to gain insights into how efficient the network is organized with regard to intercellular transport and communication. Our analysis shows that, on the canalicular level, highly organized fibrolamellar bone from sheep, which grows slower but in a more controlled way on top of a template surface, is less interconnected, but more efficiently organized than woven mouse bone, which grows faster but in the absence of an organizing surface. Despite pronounced differences at the tissue level, the topological architecture of the osteocyte canalicular network at the subcellular level is identical across sample types, suggesting a universal growth mechanism during network formation. Our method could be useful for comparing quantitatively the quality of the network of canaliculi in bone of different tissue and individuum age, loading condition and disease state. The results presented here provide insight in which principles might be at work when individual bone cells self-organize into a large, interconnected network during bone formation and mineralization.

We thank H Schell from the Julius Wolff Institute of Charité, Berlin for providing ovine bone samples collected in the SFB 760, and the Institute of Biology of the University of Potsdam for providing the murine bone samples.