Reduced synchroneity of intra-islet Ca2+ oscillations in vivo in Robo-deficient β cells

The spatial architecture of the islets of Langerhans is hypothesized to facilitate synchronized insulin secretion among β cells, yet testing this in vivo in the intact pancreas is challenging. Robo βKO mice, in which the genes Robo1 and Robo2 are deleted selectively in β cells, provide a unique model of altered islet spatial architecture without loss of β cell differentiation or islet damage from diabetes. Combining Robo βKO mice with intravital microscopy, we show here that Robo βKO islets have reduced synchronized intra-islet Ca2+ oscillations among β cells in vivo. We provide evidence that this loss is not due to a β cell-intrinsic function of Robo, mis-expression or mis-localization of Cx36 gap junctions, or changes in islet vascularization or innervation, suggesting that the islet architecture itself is required for synchronized Ca2+ oscillations. These results have implications for understanding structure-function relationships in the islets during progression to diabetes as well as engineering islets from stem cells.


Introduction
The islets of Langerhans, which comprise the endocrine pancreas, are highly organized microorgans responsible for maintaining blood glucose homeostasis. Islets are composed of five endocrine cell types (a, b, d, PP, and e), which, in adult rodents, are arranged such that the b cells cluster in the core of the islet, while other non-b endocrine cells populate the islet periphery (Orci and Unger, 1975;Pfeifer et al., 2015;Erlandsen et al., 1976). This configuration heavily prioritizes homotypic over heterotypic contacts between endocrine cells. Mature human islet architecture is more complex, and, though its exact organization pattern is still debated, it appears to follow a non-random distribution of the different endocrine cell types (Bonner-Weir et al., 2015;Bosco et al., 2010;Brissova et al., 2005;Cabrera et al., 2006;Dybala and Hara, 2019) and prioritize homotypic over heterotypic interactions between endocrine cell types (Hoang et al., 2014;Kilimnik et al., 2012). Though islet architecture is thought to be fairly static in healthy adult humans and rodents, it is quite dynamic during development, pregnancy, and diabetes (Dolensˇek et al., 2015;Miller et al., 2009;Nair and Hebrok, 2015;Sharon et al., 2019;Sznurkowska et al., 2020;Bocian-Sobkowska et al., 1999;Kharouta et al., 2009;Rieck and Kaestner, 2010). In diabetes, rodent and human islet architecture is disrupted, resulting in endocrine cell intermixing and a reduction in the ratios of homotypic versus heterotypic cell contacts (Brereton et al., 2015;Kilimnik et al., 2011;Nir et al., 2007;Striegel et al., 2015;Xiao et al., 2018). How islet architecture and endocrine cell-type sorting within the islet affect islet function remains poorly understood (Steiner et al., 2010).
Homotypic interactions between b cells are important for synchronous glucose-stimulated insulin secretion (GSIS). In GSIS, glucose from the blood enters b cells through glucose transporters and is metabolized, causing an increase in the ratio of intracellular ATP/ADP. This increase in ATP/ADP triggers closure of ATP-sensitive K + channels, resulting in membrane depolarization and the opening of voltage-gated Ca 2+ channels, triggering an influx of Ca 2+ into the cell, which in turn promotes exocytosis of insulin granules (Ammälä et al., 1993;Ashcroft and Rorsman, 2013;Bertram et al., 2010). This process is cyclical and thus oscillatory. In dispersed islets, which have no cell-cell contacts, b cells display heterogenous uncoordinated oscillations in membrane depolarization and Ca 2+ influx at both basal and elevated glucose, which results in high basal insulin secretion and low uncoordinated GSIS (Halban et al., 1982). In intact islets, where homotypic interactions between b cells are present, the oscillations in membrane potential, intracellular Ca 2+ ([Ca 2+ ] i ), and insulin secretion that underlie GSIS are respectively synchronous at elevated glucose and relatively silent at low glucose. This synchronous or 'pulsatile' pattern of insulin secretion from islets is thought to underlie pulsatile insulin levels in circulating blood, a quality important for keeping the peripheral tissues insulin sensitive, for robust liver response to insulin signaling, and for allowing time for replenishment of the readily releasable pool of insulin granules in b cells. Thus, perturbations in pulsatility are thought to contribute to diabetes pathology (Lang et al., 1981;Matveyenko et al., 2012;Satin et al., 2015;Pedersen and Sherman, 2009). Indeed, pulsatile insulin levels in circulating blood of mouse models of diabetes and patients with diabetes, pre-diabetes, and non-diabetic family members of diabetic patients are often disrupted (Satin et al., 2015;O'Rahilly et al., 1988).
One mechanism through which homotypic b cell-b cell contacts can synchronize GSIS is through electrical coupling. Within an islet, b cells are electrically coupled to their homotypic neighbors via varying levels of Connexin36 (Cx36) gap junctions that allow for exchange of cations between neighboring b cells (Benninger et al., 2008;Farnsworth et al., 2014). b cells coupled by these gap junctions display varying levels of excitability and metabolic rates, aspects which themselves are thought to display a non-random spatial organization (Hraha et al., 2014). This coupling of spatially organized heterogenous b cells populations creates an electrical syncytium that responds homogenously to glucose such that, at low glucose levels, insulin secretion is inhibited through repression of [Ca 2+ ] i oscillations across the islet, but at high glucose levels, pulsatile insulin secretion occurs through synchronous [Ca 2+ ] i oscillations that spread in fast waves across the islet from distinct initiation sites (Benninger et al., 2008;Farnsworth et al., 2014;Hraha et al., 2014;Skelin Klemen et al., 2017;Westacott et al., 2017a;Kravets, 2020;Ravier et al., 2005;Speier et al., 2007;Head et al., 2012). Indeed, in silico modeling experiments showed that decreasing the ratio of homotypic b cell-b cell nearest neighbors is predicted to result in perturbation to synchronous [Ca 2+ ] i oscillations (Hoang et al., 2014;Head et al., 2012;Nittala et al., 2007;Hoang et al., 2016). A highly functionally connected subpopulation of b cells, termed 'hub' or 'leader' cells, may also direct synchronous [Ca 2+ ] i , oscillations by harboring a disproportionate amount of functional connections to other b cells, a property allowed for by a high frequency of homotypic interactions between b cells (Johnston et al., 2016;Salem et al., 2019).
It is thus hypothesized that disrupting proper endocrine cell type sorting in the islet in a way that distorts the relative amount of homotypic b cell-b cell contacts, even without affecting any other property of the cell, would be sufficient to disrupt synchronized oscillatory behavior among b cells. However, direct empirical evidence supporting this hypothesis is lacking. This is because of the fact that while most genetic mouse models that show abnormally disorganized islet architecture also display defects in glucose homeostasis (Brereton et al., 2015), the disrupted islet architecture is usually linked to either developmental defects in b cell differentiation or maturation (Hang et al., 2014;Yamagata et al., 2002;Gannon et al., 2000;Gu et al., 2010;Ahlgren et al., 1998;Borden et al., 2013;Doyle and Sussel, 2007;Sinagoga et al., 2017;Huang et al., 2018;Jimenez-Caliani et al., 2017;Bastidas-Ponce et al., 2017;Crawford et al., 2009) or to pathologies related to b cell damage in diabetes (Kilimnik et al., 2011;Nir et al., 2007;Szabat et al., 2016;Brissova et al., 2014;Baetens et al., 1978;Starich et al., 1991). This introduces a strong confounding factor for studying the role of islet architecture in islet function. Therefore, it is difficult to disentangle the relative effect of b cell-intrinsic defects and whole-islet architectural defects, such as reduced ratio of homotypic b cell-b cell contacts, on perturbation of synchronous oscillations that underlie normal GSIS.
Recently, we have described a mouse model in which the cell-surface receptors Robo1 and Robo2 are deleted specifically in b cells (Robo bKO), resulting in disruption of canonical endocrine cell-type sorting within the islets (Adams et al., 2018). Unlike other models of disrupted islet architecture and endocrine cell-type sorting in the islet, the b cells in the islets of Robo bKO express normal levels of markers for b cell differentiation, functional maturity, and regulation of GSIS, and show normal b cell-intrinsic response to glucose. We reasoned that this model would allow us to test the role of islet architecture and endocrine cell-type sorting in regulating synchronous [Ca 2+ ] i oscillations in response to glucose among b cells in a fully differentiated, non-diabetic islet setting.

Results
Robo bKO islets have a decreased ratio of homotypic b cell-b cell contacts but their b cells are functionally mature In silico simulations where the degree of b cell-b cell coupling is changed through a decrease in homotypic nearest neighbors predict that disruption of islet architecture will disrupt synchronous intra-islet [Ca 2+ ] i oscillations and pulsatile hormone secretion (Hoang et al., 2014;Nittala et al., 2007;Hoang et al., 2016). To test whether b cells in Robo bKO islets have a decreased frequency of homotypic b cell-b cell contacts compared to control islets, we performed nearestneighbor analysis on islets from pancreatic sections of Robo bKO and control mice ( Figure 1A). We found that Robo bKO islets possess a significantly decreased frequency of homotypic b cellb cell contacts and a significantly increased frequency of heterotypic b cell contacts with either a or d cells ( Figure 1B, C). a and d cells trended towards a higher ratio of homotypic contacts in Robo bKO compared to controls though this did not reach statistical significance (p=0.06) and showed no difference in the level of heterotypic contacts with each other compared to controls ( Figure 1D, E). Thus, Robo bKO islets have a lower frequency of homotypic b cell-b cell contacts and increased frequency of b cell heterotypic contacts compared to control islets while non-b endocrine cells are less affected.
We have previously shown that genetic deletion of Robo1 and Robo2 selectively in b cells using either Ins1 Cre ;Robo1 D/D Robo2 flx/flx or Tg(Ucn3-Cre);Robo1 D/D Robo2 flx/flx mice (Robo bKO) results in disrupted islet architecture and endocrine cell-type sorting without affecting b cell death or the expression of the b cell maturation markers MafA and Ucn3 (Adams et al., 2018). To verify that b cells in Robo bKO islet are more broadly mature, we expanded the analysis to look at transcript levels of 15 additional b cell maturity markers. RNA sequencing and differential gene expression analysis on FACS-purified b cells from both Robo bKO and control islets revealed no change in transcript levels of any hallmark b cell maturity or differentiation genes (Figure 1-figure supplement 1, Supplementary file 1). Thus, unlike other mouse models with disrupted islet architecture, b cells in Robo bKO islets appear to maintain maturity and differentiation despite loss of normal islet architecture.
Though Robo bKO b cells show normal maturity and differentiation as determined by RNA sequencing and immunostaining for maturity and differentiation factors, experiments have shown that the Slit-Robo signaling pathway is involved in the stimulus secretion cascade linking glucose to insulin secretion in in vitro cultured b cells (Yang et al., 2013). Confoundingly, our RNAseq analysis showed no change in the expression of the major genes responsible for regulation of Ca 2+ dynamics or GSIS between control and Robo bKO b cells (Figure 1-figure supplement  1). Thus to test whether Robo bKO b cells are able to undergo normal [Ca 2+ ] i oscillations in Figure 1. Robo bKO islets have a decreased ratio of homotypic nearest neighbors than controls. (A) Immunofluorescence images (left and middle panels) and cell connectivity maps generated by nearest-neighbor analysis (right panels) of control and Robo bKO islets. b cells (red), a cells (green), and d cells (blue) are denoted by nodes on the connectivity maps. A line the same color as both nodes it connects denotes a homotypic interaction of that corresponding cell type. A white line connecting two nodes denotes a heterotypic interaction between cell types. Scale bars are 100 mm. (B) Figure 1 continued on next page response to stimuli, we performed in vitro [Ca 2+ ] i imaging on single b cells from dissociated Robo bKO and control islets exposed to glucose followed by KCl ( Figure 2). We found no difference in the proportion of b cells that undergo [Ca 2+ ] i oscillations in response to 10 mM glucose between control and Robo bKO b cells (Figure 2A-C). We also saw no significant difference in area under the curve (AUC) of [Ca 2+ ] i in response to 10 mM, in peak height corresponding to the first phase of insulin secretion, or in AUC [Ca 2+ ] i in response to KCl in control and Robo bKO b cells ( Figure 2D-H). Together, this demonstrates that Robo bKO b cells are not defective in their ability to undergo [Ca 2+ ] i oscillations in response to glucose, suggesting that b cell-intrinsic factors mediated by Robo deletions do not have a strong impact on the [Ca 2+ ] i dynamics that underlie GSIS in individual b cells.
We reasoned that the altered degree of homotypic b cell-b cell interactions in Robo bKO islets, together with the retained b cell maturity and intact intrinsic b cell [Ca 2+ ] i dynamic functionality, provides a unique model by which to empirically test the hypothesis that endocrine cell-type organization within the islet affects synchronous oscillatory behavior in intra-islet b cells.

Robo bKO islets display unsynchronized Ca 2+ oscillations in vivo
We thus set out to investigate how the reduced homotypic b cell-b cell interactions in Robo bKO islets affect synchronous b cell behavior by measuring dynamic [Ca 2+ ] i oscillations in response to glucose. Robo bKO islets are fragile in isolation and culture (Adams et al., 2018), making them unsuitable for in vitro analyses of whole-islet [Ca 2+ ] i oscillations. To overcome this limitation, we adopted an intravital [Ca 2+ ] i imaging method that enables imaging of islet [Ca 2+ ] i dynamics in situ within the intact pancreas (Reissaus et al., 2019). In brief, this method employs an intravital microscopy (IVM) platform and adeno-associated viral (AAV) delivery of insulin promoter-driven GCaMP6s, a fluorescent Ca 2+ biosensor, to quantitate b cell [Ca 2+ ] i dynamics in vivo in both Robo bKO and control islets. This method also allows for retention of the islet's in vivo microenvironment, blood flow, and innervation, thus providing a more realistic condition than in vitro approaches to study islet function allow for.
We verified that synchronous [Ca 2+ ] i oscillations are maintained in vivo in islets by measuring GCaMP6s intensity in b cells within AAV8-RIP-GCaMP6-infected islets of control (Robo WT) mice ( Figure 3). As expected, control mice displayed whole-islet synchronous [Ca 2+ ] i oscillations when imaged at anywhere from~0.03 to~1 Hz for at least 10 min after glucose elevation (Figure 3 and Videos 1 and 2). We quantified the level of synchroneity in intra-islet oscillations by analyzing the correlation between GCaMP6s active areas within individual islets ( Figure 3C, D). In brief, this analysis measures the proportion of GCaMP6s active area within an islet where the normalized GCaMP6s intensity over time has a >0.7 Pearson's correlation coefficient. This is demonstrated by synchronous region maps, which show highly correlated areas of the islet in the same color, revealing that almost all of the active [Ca 2+ ] i area is synchronized within control islets ( Figure 3D). While oscillations vary in frequency between islets, the area of highly correlated [Ca 2+ ] i oscillations between b cells within any one islet is very high, confirming that control islets possess highly synchronous intra-islet [Ca 2+ ] i oscillation in response to glucose in vivo ( Figure 3 and Figure 4E).
Conversely, we found that Robo bKO islets on average display less synchronous intra-islet [Ca 2+ ] i oscillations in vivo when imaged at speeds from~0.03 Hz to~1 Hz (Figure 4, Figure 4-figure supplement 1, and Videos 3, 4, and 5). Quantification of this asynchronous behavior through Figure 1 continued Probability of b cell homotypic contact in Robo bKO islets vs. controls. n = 4 mice; control 0.74 ± 0.02 SEM, Robo bKO 0.53 ± 0.04 SEM, p<0.005 t-test. (C) Probability of b cell heterotypic contacts in Robo bKO islets vs. controls, n = 4 mice, control 0.13 ± 0.02 SEM, Robo bKO 0.28 ± 0.02 SEM, p<0.05 MW. (D) Probability of aa or dd homotypic contacts in Robo bKO islets vs. controls, n = 4 mice; control 0.07 ± 0.01 S.EM, Robo bKO 0.13 ± 0.02 SEM, p=0.06 MW. (E) Probability of a-d heterotypic contact in Robo bKO islets vs. controls, n = 4 mice; control 0.03% ± 0.01 SEM, Robo bKO 0.06% ± 0.02 SEM, p=0.22 t-test. (B-E) Similar shaded points in graphs indicate islets from the same mouse, while mean and error bars represent statistics performed on average values from each mouse. Error bars show SEM. 9-11 islets from an individual mouse were measured as technical replicates, and the average values per mouse were used as biological replicates. MW: Mann-Whitney; SEM: standard error of the mean. The online version of this article includes the following source data and figure supplement(s) for figure 1:   , which is similar to levels in control islets. This synchronous population showed no difference in phase lag or wave speed when compared to controls, indicating that they share similar wave dynamics in addition to synchroneity with control islets (Figure 4-figure supplement 2B, C). The existence of this highly synchronous population of Robo bKO islets further suggests that the ability of individual b cells to oscillate [Ca 2+ ] i in response to stimuli is unaffected by deletion of Robo. This supports the idea that Robo bKO b cells do not have intrinsic defects in [Ca 2+ ] i dynamics in response to glucose, and that instead these oscillation defects are due to b cell-extrinsic factors within the islets of Robo bKO mice.

Robo bKO islets have altered functional network properties
To gain a more in-depth understanding of how Ca 2+ dynamics are altered in Robo bKO islets in vivo, we performed networks analysis based on the methods previously described by Stožer et al., 2013. We thus analyzed how oscillations in b cells within an islet correlate with each other and created a network map of single planes from [Ca 2+ ] i recordings that had high enough resolution to identify individual cells. In these network maps, the correlation threshold for functional connections was set to R = 0.95, thus any b cells that had a correlation coefficient R ! 0.95 were considered functionally connected ( Figure 5A). In line with results from the correlated area analysis ( Figure 4E), the average functional connectivity between all b cell pairs within an islet was significantly lower in Robo bKO islets compared to controls, with Robo bKO islets that showed high areas of correlation also showing network connectivity levels similar to controls ( Figure 5A, B). This decrease in functional connectivity was also demonstrated by the difference in the probability distributions of average percent links in control and Robo bKO islets ( Figure 5C). Interestingly, we also found that the standard deviation of average islet connectivity was significantly lower in Robo bKO islets compared to controls, indicating that individual b cells within Robo bKO islets showed a larger variation in connectedness than b cells in control islets ( Figure 5D). We also measured the number of highly connected b cells within the islets, or 'hub' cells, by measuring the percent of b cells that were functionally connected to !25% of the islet. We observed a trend towards a lower percentage of hub cells per islet in Robo bKO compared to controls, though this did not reach the threshold for statistical significance likely because of an outlier in the Robo bKO group (Mann-Whitney [MW] p=0.07). When we performed outlier analysis on these data, it revealed a significant decrease in the percent of hub cells in Robo bKO vs. control islets (MW p<0.05; Figure 5E).
Further, we measured the average clustering coefficient and global efficiency of these networks in order to characterize their 'small world' properties which have been shown to characterize wildtype islet [Ca 2+ ] i oscillations networks (Stožer et al., 2013;Watts and Strogatz, 1998). Small world networks are characterized by a high level of local clustering that is reflected by a high average clustering coefficient (C avg ), but also a high level of integration that is reflected by a high global efficiency (E glob ). In our analysis, we saw a trend toward decrease in both C avg and E glob in Robo bKO islets compared to controls though this did not quite reach statistical significance (p=0.05 for both MW; Figure 4-figure supplement 2F, G), indicating that Robo bKO islets act less like small world networks than control islets.
Altogether, network analysis suggests that Robo bKO islets have higher variability in functional connections per b cell within an islet, are less functionally connected on average, have a lower frequency of hub cells, and are less small world-like when compared to controls. This supports a scenario in which having less homotypic b cell-b cell interactions reduces the level of coupling within an islet.

Robo bKO islets do not show defects in innervation and vascularization
One possible explanation for the observed disruption of intra-islet [Ca 2+ ] i oscillations in Robo bKO islets is through disruptions in innervation or vasculature within the islet as these factors have been implicated in controlling the ability of b cells to synchronize and the level of glucose exposure across the islet (Eberhard and Lammert, 2009). To assess whether the amount of innervation within the islets of Robo bKO differs from that of control islets, we stained pancreatic sections for the pan-neuronal marker Tubb3 and quantified the area of nerves normalized to islet area, and saw no difference in innervation between control and Robo bKO islets (Figure 6figure supplement 1). To determine if there were changes in the amount of intra-islet vasculature in Robo bKO mice, we quantified the amount of matrix components secreted by vessels as a surrogate for vasculature (laminin and collagen IV) in Robo bKO and control islets. We observed no significant difference in area of vessel matrix components between Robo bKO and control islets for either of these proteins ( Figure 6A-D). Despite the similar amount of vasculature in Robo bKO islets, it is possible that the pattern of vascularization may differ in such a way that delivery of glucose to certain regions of the islet is perturbed or delayed due to insufficient vascular coverage throughout the islet. To test whether glucose perfusion from the blood vessels was similar across the islet in Robo bKO compared to controls, we performed intravital imaging of islets during intravenous injection of the fluorescent glucose analog 2-NBDG ( Figure 6E-K and Videos 8 and 9). We then assessed the timing of glucose arrival to areas <10 mm and >10 mm from the closest vessel to assess whether glucose can reach all portions of Robo bKO islets within a similar time frame as controls ( Figure 6F-H and Videos 8 and 9). No   Figure 4 continued on next page significant differences in timing of glucose arrival between randomly sampled regions of the islet adjacent to blood vessels and regions >10 mm from the closest vessel in control and Robo bKO were seen ( Figure 6F, G). We also saw no change in the AUC of average normalized 2-NBDG intensity in regions close to vessels within the same islet in controls compared to Robo bKO and in regions far from vessels between controls and Robo bKO ( Figure 6H, I). Further, to measure the variability of 2-NBDG perfusion across individual islets, we measured the standard deviation of areas <10 mm from vessels within the same islet and areas >10 mm within the same islets in control and Robo bKO and saw no difference ( Figure 6J, K). Altogether, this indicates that all areas of the islet are exposed to glucose at essentially the same time in both control and Robo bKO islets.

Robo bKO islets retain normal expression levels and localization of Cx36 gap junctions
Another possible explanation for the loss of synchronized whole-islet [Ca 2+ ] i oscillations in Robo bKO is that the gap junction protein Cx36 is mislocalized or mis-expressed when Robo is deleted, and thus cell coupling is inhibited. Indeed, the phenotype described above for Robo bKO islets is reminiscent of that observed in mice heterozygous for a Gjd2-null allele (Benninger et al., 2008;Ravier et al., 2005). To test whether Cx36 is mis-localized in Robo bKO b cells, we stained for Cx36 along with F-actin to visualize the cell borders in b cell lineage-traced Robo bKO and control tissue sections ( Figure 7A) and saw that in Robo bKO islets Cx36 still localized normally to the b cell borders. Further, to test whether Robo bKO mice downregulate Cx36 in islets, we measured the area of Cx36 protein immunofluorescence normalized to islet area in Robo bKO and control islets ( Figure 7B, C). We found no difference in Cx36 area between Robo bKO islets and controls ( Figure 7C). We also verified the specificity of our Cx36 antibody by staining Cx36 -/tissue   sections, which showed no staining within the islet ( Figure 7B). This is consistent with our RNAseq analysis that showed no change in Cx36 expression in Robo bKO compared to controls (Figure 1figure supplement 1). Overall, this suggests that loss of synchronous intra-islet [Ca 2+ ] i oscillations is not due to decreased expression or mis-localization of Cx36 in Robo bKO b cells. Instead, this is consistent with a scenario in which coupling within Robo bKO islets is decreased due to a decrease in the ratio of b-b cell contacts rather than a loss of the gap junction machinery itself.

Discussion
In this study, we provide evidence for the importance of islet architecture for proper islet function in vivo. We show that in Robo bKO islets, which possess disorganized islets architecture that results in a lower frequency of homotypic b cell-b cell interaction, synchronized [Ca 2+ ] i oscillations in the islet are disrupted. This disruption of synchronized [Ca 2+ ] i oscillations is not due to loss of functional b cell maturation and identity, altered Robo-mediated b cell-intrinsic defects in glucose-stimulated Ca 2+ dynamics, loss of Cx36 expression or localization to the b cell membrane, loss of innervation, or change in amount or pattern of islet vascularization, suggesting that endocrine cell-type sorting within the islet is, by itself, important for islet function, at least in the aspect of coordinated [Ca 2+ ] i oscillations among b cells.
Robo, and its ligand Slit, have been previously shown to affect [Ca 2+ ] i oscillations in b cells in vitro (Yang et al., 2013). However, Robo-mediated b cell-intrinsic effects are likely not the cause of asynchronous in vivo [Ca 2+ ] i oscillations in Robo bKO islets. This is supported by two independent observations: (1)   lineage-tracing reporter (Adams et al., 2018), which uses the same Cre deriver that is used to delete Robo in those b cells. Thus, high expression of H2B-mCherry suggests efficient recombination of the Robo floxed allele. All this together support the conclusion that loss of synchronized [Ca 2+ ] i oscillations is caused by a b cell-extrinsic mechanism, such as the reduction in frequency of homotypic interactions between b cells observed in Robo bKO islets.
Robo receptors have known roles in angiogenesis and axon guidance, and thus could affect precisely how the islet is innervated and vascularized (Blockus and Chédotal, 2016). However, we found that the amount of innervation and vascularization between Robo bKO islets and controls is similar, suggesting that abnormal vascularization or innervation are likely not the cause of disruption in [Ca 2+ ] i oscillations found in Robo bKO islets. Additionally, our intravital imaging experiments show that the fluorescent glucose analog 2-NDBG delivered intravenously reaches all portions of Robo bKO islets within the same time frame as control islets. Thus, incomplete or delayed glucose perfusion across the islet due to improper islet vascularization is likely not the cause of asynchronous [Ca 2+ ] i oscillations in Robo bKO islets.
One possible cause for asynchronous [Ca 2+ ] i oscillations in Robo bKO islets is through disrupted coupling between b cells. Interestingly, the observed [Ca 2+ ] i oscillation phenotype in Robo bKO islets is reminiscent of the phenotype seen in heterozygous Cx36 mutants where the level of synchroneity in whole-islet [Ca 2+ ] i is highly variable from islet to islet. This suggested that Cx36-mediated gap junction coupling may be reduced in Robo bKO islets, possibly through reduction or mis-localization of Cx36. However, we found that the area of Cx36 was not different between Robo bKO and controls, and that Cx36 gap junctions still localized to the correct intracellular domain at b cell borders in Robo bKO islets. Thus, we hypothesized that coupling may be reduced instead by a reduction in the number of homotypic b cell-b cell interactions across the islet, which could reduce the overall functional connectivity. In line with this hypothesis, we found a significant reduction in b cell homotypic interactions in Robo bKO islets compared to controls, and a significant reduction in islet    Figure 6 continued on next page connectivity as determined by network analysis. We also observed a significant increase in the standard deviation of average islet connectivity in our network analysis, which indicates that cells within Robo bKO islets show more variability in their level of functional connections than controls. We further hypothesized that if homotypic b cell coupling is disrupted across the islet by a decrease in frequency of homotypic b cell interactions, it could possibly affect the number of hub cells within Robo bKO islets. Since hub cells are thought to direct synchronous [Ca 2+ ] i oscillations among b cells within an islet, it is conceivable that disruption of this subpopulation would affect synchronous GSIS. Indeed, our network analysis revealed a trend towards a lower percent of hub cells found in Robo bKO islets compared to controls that reached statistical significance after outlier analysis. Further, we found that the average clustering coefficient (C avg ) and global efficiency (E glob ) of the networks in Robo bKO islets trended towards being lower than controls (p=0.05 for both), which is indicative of a reduction in small world properties of the network. Taken together, these results suggest that reduction in frequency of homotypic interactions reduces overall coupling and increases the variability of coupling between b cells in the islet, thus disrupting synchronous [Ca 2+ ] i oscillations, possibly in part due to a reduction in hub cells and loss of small world network properties. However, a direct comparison between amount of homotypic interactions and the degree of gap junction coupling or synchroneity in [Ca 2+ ] i oscillations within the same islet remains to be tested as this task is beyond the reach of the intravital imaging system used in our experiments. With the advent of the acute pancreatic slice method, this should be feasible as live measurements like [Ca 2+ ] i imaging and fluorescence recovery after photobleaching (FRAP) analysis to test gap junction coupling, and post hoc staining of fixed slices for architectural analysis of the same slice are possible (Marciniak et al., 2014). An important future experiment would thus be to test the level of gap junction coupling and synchroneity of [Ca 2+ ] i oscillations, along with the frequency of homotypic interactions within the same islet to directly Figure 6 continued mice, control 0.17 ± 0.02 SEM, Robo bKO 0.17 ± 0.02 SEM, p=0.99 MW. >12 islets from an individual mouse were measured as technical replicates, and the average values per mouse were used as biological replicates. (D) Quantification of area of laminin staining normalized to islet area showing no difference in amounts of basement membrane marking blood vessels in Robo bKO compared to control islets. Control n = 8 mice, Robo bKO n = 8 mice, control 0.31 ± 0.05 SEM, Robo bKO 0.38 ± 0.05 SEM, p=0.39 t-test. >12 islets from an individual mouse were measured as technical replicates, and the average values per mouse were used as biological replicates. (E) 2-NBDG perfused through the islet vasculature in a control and Robo bKO islet. (F) Average normalized 2-NBDG intensity (F/F max ) over time in areas < and >10 mm from nearest vessel in control islets. (G) Average normalized 2-NBDG intensity (F/F max ) over time in areas < and >10 mm from nearest vessel in Robo bKO islets. (H) Area under the curve (AUC) of 2-NBDG intensity over time in areas <10 mm from vessels in control and Robo bKO islets. Control 6.9 ± 0.26 SEM, Robo bKO 6.9 ± 0.14 SEM, p>0.90 t-test. (I) AUC of 2-NBDG intensity over time in areas >10 mm from vessels in control and Robo bKO islets. Control 6.8 ± 0.24 SEM, Robo bKO 6.8 ± 0.17 SEM, p>0. https://elifesciences.org/articles/61308#video8 correlate frequency of b cell homotypic contacts to the level of coupling and synchronous oscillatory behavior using this method.
Another possible factor that contributes to the aberrant [Ca 2+ ] i oscillations in Robo bKO islets is disruption of paracrine and autocrine signaling between endocrine cells, caused by the increased heterotypic and decreased homotypic b cell interactions in these islets. The spatial islet architecture dictates the types and amounts of paracrine and autocrine interactions among the different cell types, and these interactions are thought to regulate synchronous GSIS (Benninger and Hodson, 2018). For example, a cells and d cells, which secrete glucagon and somatostatin respectively, were shown to affect cAMP oscillations in b cells, which are important for synchronizing b cell Ca 2+ oscillations during the second phase of insulin secretion (Hodson et al., 2014;Tian et al., 2011;Grapengiesser et al., 2003). b cells may also be able to electrically couple to neighboring d cells, and depolarization of b cells may thus trigger secretion of somatostatin from connected d cells, which in turn could affect islet Ca 2+ and electrical dynamics (Briant et al., 2018). Additionally, ATP, which is co-secreted with insulin from b cells, has an autocrine effect such that it binds purinergic receptors on b cells to help synchronize electrical and Ca 2+ dynamics (Hellman et al., 2004;Tudurí et al., 2008). Moreover, receptor-ligand interactions occurring between neighboring b cell have also been shown to affect GSIS independent of electrical coupling. This is demonstrated by the fact that dissociated wildtype b cells have worse GSIS than intact islets lacking Cx36 gap junctions (Benninger et al., 2011). EphA-ephrinA binding between b cells has been identified as a mechanism that can govern this coupling-independent enhancement of GSIS (Konstantinova et al., 2007). These interactions rely on b cell-b cell contact and thus could be disrupted if the ratio of homotypic interactions in the islet were reduced; however, this particular effect is likely downstream of [Ca 2+ ] i oscillations and thus less likely to contribute to the oscillation phenotype observed in Robo bKO islets (Benninger et al., 2011). Ultimately, the spatial arrangement of endocrine cell types within the islet dictates the extent to which these paracrine, autocrine, and cell surface receptorligand interactions occur within the microenvironment, and thus changing this architecture would likely affect emergent electrical and Ca 2+ dynamics to some degree. Thus, while we propose that defects in synchronous Ca 2+ oscillation in response to glucose in Robo bKO islets are largely due to inefficient b cell coupling as a result of decreased homotypic b cell interactions in these islets (Farnsworth et al., 2014), it is also possible that changes in paracrine and autocrine signaling, as well as cell surface receptor-ligand interactions, contribute to this phenomenon.
Several patterns of [Ca +2 ] i oscillations have been described in islets in vitro, in both isolated islets and acute pancreas tissue slices: slow oscillations with a period between~3 and 5 min, fast oscillations with periods <1 min, and mixed oscillations consisting of fast oscillations superimposed on the plateau fraction of slow oscillations (Ravier et al., 2002). The Dual Oscillator Model proposes that these slow Ca 2+ oscillations are controlled by intrinsic glycolytic oscillations and set the period for pulsatile insulin secretion, whereas the fast oscillations result from feedback of ions on their respective channels and control the amount of insulin secreted by dictating the plateau fraction of the slow oscillations. Thus it is of interest if these same oscillation patterns observed in vitro in wildtype islets also occur in vivo. We observed fast, slow, and mixed oscillations in vivo in control islets, with the most common pattern being mixed and the least common being fast alone, demonstrating that these patterns are not an artifact of in vitro Ca 2+ imaging but occur in vivo in response to IP glucose  (Supplementary file 2). Slow oscillations in control islets with >80% correlated area had periods ranging from~1 min to as long as~10.3 min. For fast oscillations, periods ranged from~4 s to at least~20 s and possibly as slow as~31 s though this high end of the range may be inflated due to possible undersampling from islets imaged at~0.1 Hz. However, we chose to include islets imaged at~0.1 Hz in this table because the slow oscillation speed is accurately measured and this imaging speed still confirms the presence of fast oscillations and thus provides information about oscillation pattern (i.e., mixed vs. slow). Interestingly, islets coming from the same mouse tended to have more similar oscillation periods, a property which has been observed in vitro as well (Nunemaker et al., 2005). Unfortunately, due to the restraints imposed by our current intravital imaging protocol, we cannot know the exact duration of and/or concentration of glucose exposure within the islet during imaging.
Though the unsynchronized phenotype in Robo bKO islets is clear from the videos we captured at approximately 1 Hz, 0.2 Hz, 0.1 Hz, and 0.03 Hz, we note that other recent experiments that analyzed network dynamics in islets were performed at an imaging speed of 1 Hz or faster. That said, we believe that it is important to capture [Ca 2+ ] i dynamics in 3D in order to ensure we are not biased by examining certain 'synchronized regions' that would not be representative of the islet dynamics as a whole. We would further argue that the slow oscillations we observe are on a~1-10 min time scale, and therefore sampling these oscillations even at a frequency of 1 frame per 10s of seconds is sufficient for capturing synchronicity. Additionally, when we imaged at 1 Hz, we observed the same pattern of oscillations in Robo bKO islets and controls as was observed in slower imaging speeds (Figure 4-figure supplement 2). Furthermore, our analyses of network dynamics in simulated islets showed that functional connections and phase lag can still be accurately measured with time resolution as low as 0.1 Hz (the slowest we simulated), with only a modest loss in sensitivity (data not shown). Thus, though our imaging speeds may cause underestimation of loss in coordinated area or functional connections, it is still sensitive enough to observe the robust difference in network dynamics between Robo bKO islets and controls at these imaging speeds.
We also note that the mosaic GCaMP6s expression, which results from infection of AAV8-RIP-GCaMP6s, is variable between islets and does not allow monitoring of Ca 2+ dynamics in all b cells across the islet. However, the distribution of GCaMP6s in these islets is similar if not broader than what is achievable with incubation of commonly used Ca 2+ dyes that have proven successful in measuring Ca 2+ dynamics in the islet. Specifically, the AAV8 we used penetrates into the core of the islets instead of just the superficial layers as is often seen with Ca 2+ dye incubation, allowing for equal if not greater representation of the islet in imaging. Thus, though GCaMP6s expression is mosaic and variable between islets, the expression is comparable to other common methods that have been useful in assessing Ca 2+ dynamics. However, future experiments would benefit from using the genetically encoded GCaMP6s floxed allele combined with a b cell-specific Cre driver to increase whole-islet resolution of Ca 2+ wave dynamics (Madisen et al., 2015).
Taken together, our data thus supports the hypothesis that disrupting proper endocrine celltype sorting in the islet in a way that distorts the relative amount of homotypic b cell-b cell contacts, without affecting b cell-intrinsic [Ca 2+ ] i dynamics, is sufficient to disrupt synchronized [Ca 2+ ] i oscillations among b cells. These results lend evidence to the idea that correct spatial islet architecture is important for islet function, and may have an important impact on understanding islet dysfunction in diabetes, and on approaches to generate functionally better islets from stem cells in vitro. bKO islets and controls showing no significant difference, n = 4 mice for each genotype, control 0.46 ± 0.06 SEM, Robo bKO 0.69 ± 0.11 SEM, p = 0.11 t-test. Similar shaded points in graphs indicate islets from the same mouse, while mean and error bars represent statistics performed on average values from each mouse. Error bars shown are SEM. 10-14 islets were measured per mouse as technical replicates, and the average for each mouse was considered a biological replicate. SEM: standard error of the mean.  (Branchfield et al., 2016), Ins1-Cre (Thorens et al., 2015), Urocortin3-Cre (van der Meulen et al., 2017), and Rosa26-Lox-Stop-Lox-H2BmCherry (Blum et al., 2014) mice were previously described.

Key resources table
All mouse strains were maintained on a mixed genetic background. Control colony mates in all analyses were Robo +/+ with either Ins1-Cre or Ucn3-Cre.

RNA sequencing
RNA was isolated from FACS sorted lineage-traced b cells (Adams et al., 2018) from control and Robo bKO mice using phenol chloroform extraction (TRIzol). DNA libraries were generated using Takara's SMART-Seq v4 Low Input RNA Kit for Sequencing (Takara, Mountain View, CA) for cDNA synthesis and the Illumina NexteraXT DNA Library Preparation (Illumina, San Diego, CA) kit for cDNA dual indexing. Full-length cDNA fragments were generated from 1 to 10 ng total RNA by SMART (Switching Mechanism at 5 0 End of RNA Template) technology. cDNA fragments were fragmented and dual indexed in a single step using the Nextera kit's simultaneous transposon and tagmentation step. Quality and quantity of completed libraries were assessed using Agilent DNA series chip assay (Agilent Technologies, Santa Clara, CA) and Invitrogen Qubit ds DNA HS Kit (Invitrogen, Carlsbad, CA), respectively. Each library was standardized to 2 nM. Cluster generation was performed on Illumina cBot, with libraries multiplexed for 1 Â 100 bp sequencing using TruSeq 100 bp SBS kit (v4) on an Illumina HiSeq2500. Images were analyzed using standard Illumina Pipeline, version 1.8.2.

Intravital imaging for dynamic [Ca 2+ ] i
Mouse pancreata were exposed in anesthetized mice by making a small incision on the right side of the mouse and externalizing the tip of the pancreas. A glass dish was placed over the exposed pancreas and the mouse was placed on microscope stage over an inverted objective with isoflurane anesthesia for the remainder of imaging. Islets were identified on the surface of the pancreas by detecting Histone H2B-mCherry fluorescent nuclei labeled by a b cell-specific lineage-tracing reporter (Adams et al., 2018). Once islets were identified, mice were given injections of 1 g/kg body weight glucose (30% in saline) intraperitoneally. Blood glucose levels were monitored through tail vein bleeds. Once the blood glucose reached at least 200 mg/dL, GCaMP6s activity was identified using the microscope eye piece. Time courses of GCaMP6s were imaged using a custom-built inverted multiphoton microscope (Bruker Fluorescence Microscopy, Middleton, WI), as described previously (Lugo-Cintró n et al., 2020). The system consisted of a titanium:sapphire laser (Spectra-Physics, Insight DS-Dual), an inverted microscope (Nikon, Eclipse Ti, Melville, New York, NY), and a Nikon Apo 20x objective. For recordings, z-stacks were set to 1, 3, 8, or 12 slices each 8 mm apart and images were captured at 1 Hz, 0.2 Hz, 0.1 Hz, or 0.03 Hz, respectively, over at least 10 min at a resolution of 512 Â 512 pixels. GCamP6s and mCherry signal were excited with a laser tuned to 890 or 910 nm. After time courses were recorded, high-resolution image z-stacks were taken with 60 zplanes taken 1 mm apart or 8 z-planes taken 8 mm apart at 1024 Â 1024 pixel resolution. For some images, 70 kDA rhodamine-dextran was injected retroorbitally to mark the vasculature of the islets in vivo.

Gap junction, vasculature, and innervation quantification
Cx36 levels were quantified from images of islets co-stained with rabbit anti-Cx36 (Invitrogen) and guinea pig anti-insulin (Dako) antibody. Vasculature levels were quantified from images co-stained with rat anti-laminin b1 (Invitrogen) or rabbit anti-col IV and guinea pig anti-insulin (Dako). Innervation was quantified from images of tissue sections co-stained with rabbit anti-Tubb3 (BioLegend) and guinea pig anti-insulin (Dako). For vasculature and Cx36, eight Z-planes were taken 1 mm apart, and for innervation, five z-planes were taken 1 mm apart, on a Leica SP8 Scanning Confocal microscope using a 40Â oil immersion objective (Cx36) or 20Â objective (vasculature and innervation). Filtering and thresholding were applied to both channels for each islet, and the area of each staining was measured using Fiji's analyze particles functions. The area of gap junctions, blood vessels, or innervation (marked by their respective antibody) was divided by the area of DAPI (for Cx36) or insulin (vasculature and innervation) for each islet.

Intravital imaging and analysis of 2-NBDG perfusion
Mouse pancreata were exposed in anesthetized mice by making a small incision on the right side of the mouse, and externalizing the tip of the pancreas, and a tail vein catheter was inserted. A glass dish was placed over the exposed pancreas and the mouse was placed on a microscope stage over an inverted objective with isoflurane anesthesia for the remainder of imaging. Intravital imaging was performed on a custom-built multi-photon microscope as described above. Islets were identified on the surface of the pancreas by detecting Histone H2B-mCherry fluorescent nuclei labeled by b cellspecific lineage-tracing reporter (Adams et al., 2018). Recordings were taken of a single plane within the islet at 1 Hz imaging speed and 512 Â 512 pixel resolution using a 20Â objective and the laser tuned to 820 nm. During recording, 50 mL of 5 mg/mL 2-NBDG diluted in PBS (Invitrogen N13195) were injected through tail vein catheter, with recordings continuing for 5 min after injection. Perfusion of 2-NBDG was measured in the pancreas using the ImageJ Time Series Analyzer Plugin with the following workflow. A 10 Â 10 mm nondestructive grid was overlayed on each movie and the time point with the most inter-vascular 2-NBDG intensity was used to determine blood vessel location in the islet. 5 mm diameter circular regions of interest (ROIs) were added randomly across the islet in areas within 10 mm of vessels and areas >10 mm away from the closest vessel as determined by overlayed grid. Average intensities of ROIs were then measured over the time span 5 s before and 10 s after the frame with the first >5 unit increase of intensity in the 2-NBDG channel averaged over the entire image, marking the entry of 2-NBDG into the islet. The average intensity for each ROI was then normalized to basal-peak intensity (F/F max ) by subtracting the minimum intensity over the À5 to 10 s time frame from each time point and then multiplying those values by 1/ (max intensity À min intensity). The normalized average intensities for all ROIs within the same islet within 10 mm and >10 mm from vessel ROIs were then averaged respectively to give a value for average 2-NBDG perfusion close to and far away from vessels for each islet at each time point. The AUC was measured using Prism (GraphPad) for the normalized average 2-NBDG intensity over À5 to 10 s for regions within 10 mm and >10 mm for each islet.
Nearest-neighbor analysis b cells were identified using the lineage tracer Rosa26-Lox-Stop-Lox-H2BmCherry crossed to Ucn3-Cre, and tissue sections were stained with antibodies against glucagon and somatostatin to identify a and d cells, respectively. The 3D Tissue Spatial Analysis Toolbox for Fiji (Tran Thi Nhu et al., 2017) was used to identify specific cell types using the above markers and calculate the number of celltype-specific nearest neighbors from all identified endocrine cells.
Correlated area, wave speed, and phase lag analysis All images were analyzed using previously published methods (Westacott et al., 2017b) with custom MATLAB (MathWorks) scripts. Islets were defined as all mCherry+ cells that were direct neighbors with multiple other mCherry+ cells within a cluster. Islets were initially selected as showing activity through blinded selection by two separate researchers who are knowledgeable about islet [Ca 2+ ] i dynamics and who did not perform the imaging. If either researcher determined any islet to be inactive, the islet was not used for analysis. From this assessment, 16 out of 22 control islets and 18 out of 23 Robo bKO islets were determined to be active. Subsequent analysis was also performed blinded to genotype. For activity analysis, images were smoothed using a 5 Â 5 pixel averaging filter. Areas without significant fluorescence were removed. Saturated areas were also removed by limiting the area to intensity below the maximum value. Photobleaching was low, and as such could be approximated as a linear decline. Thus any linear trend was removed to correct for minimal photobleaching. Any islets with significant motion artifacts were removed or time courses were shortened to the time over which no significant movement occurred (displacement of <0.5 cell width). For the time course of each pixel in the image with significant fluorescence, a peak detection algorithm was used to determine if the areas had peak amplitudes significantly above background (Westacott et al., 2017b). A region was considered 'active' if the corresponding time course for each pixel had a peak amplitude >1.5Â background. The fraction of active area was calculated as the number of pixels detected as 'active' across all z-planes, normalized to the total number of pixels that showed significant fluorescence across all z-planes that were not saturated. Coordination was determined based on coincident timing of identified peaks, where areas were segmented by identified peaks occurring at similar time points (Westacott et al., 2017b). The cross-correlation of the time courses for two 5 Â 5 pixel subregion was taken. If the correlation coefficient was >0.7, then the two subregions were considered highly coordinated and merged into a larger region. The coordinated area was calculated as the number of pixels in the largest area of coordination across all z-planes normalized to the total number of pixels of the islet that were determined to be 'active' for all planes. This analysis is based on previous analysis (Westacott et al., 2017b), but adjusted for three-dimensional data.
Phase lag and wave propagation speed was determined, as in Benninger et al., 2008. For every 5 Â 5 region, the phase lag was calculated from a Fourier transform of each time course. Only regions with a correlation coefficient of >0.7 when compared to the islet average were used. First, the peak frequency was identified from the power spectrum, as generated from a Fourier transform of the islet average time course. The phase lag was only calculated for the peak frequency. The phase lag was calculated from the difference in phase between the maximum phase region and the minimum phase region and converted into a time lag according to (dt = (1/f) * tan^(À1)(f1 À f2)), where j1Àj2 is the phase difference, and f is the sampling frequency of the time courses analyzed. Speed was calculated by dividing the distance between the maximum phase lag and minimum phase lag regions by the time lag. For phase and speed analysis, only islets with >80% coordinated area was used. If the islets did not show a clear phase transition across the islet, they were excluded from the analysis.
All statistical analysis was performed in Prism (GraphPad) or MATLAB. If the data were determined to follow a normal distribution, an F-test was used to determine if variances were equal, then a Student's t-test or Welch t-test (for unequal variance) were utilized. If a normal distribution could not be verified for the data, a Mann-Whitney test was utilized for determining significance. p<0.05 was considered significant.

Network and hub analysis
The network analysis is based upon that previously applied to study islet Ca 2+ dynamics (Johnston et al., 2016;Stozˇer et al., 2013). First, a single plane in the z-stack was chosen that showed the greatest coverage of GCaMP6s signal. Visible cells were manually identified and the Ca 2+ time courses for each cell were derived. To remove bias, an equal number of islets from each temporal resolution were used for both the control and knockout. There was no difference in results between islets imaged at a speed of 0.1, 0.2, and 1 Hz. Adjacency matrices for an islet were created using the corr() function in MATLAB to compare the entire time course between each cell and calculate the Pearson's correlation coefficient between cells in a pairwise manner. Diagonals were set to 0 as cells are not considered synchronized with themselves. The control group was used to set the R_th threshold (as set forth by Stozˇer et al., 2013). All control time courses showed a power law with R_th = 0.98. However, the threshold for analysis was set to 0.95 because most Robo bKO islets did not show any connectivity with R_th at 0.98, which thus prevented comparisons. This threshold was statistically significant for p<0.001 for all islets.
To identify b cell hubs, a probability distribution was created by counting the number of cells that had X number of links, for all values of X. Cells with zero links were not plotted in the histogram but were considered in the subsequent analysis. Hubs were identified based on the percent of cells that are synchronized (linked) being more than 25% of the islet. This allows for size invariance and consistency for comparison between control and KO groups. In contrast in the methods set forth in Johnston et al., 2016, hubs are considered cells with more than 60% of the islet's connections. The number of hubs was then normalized by islet size. This is dependent on the individual islet probability distribution function, and inter-islet comparison is contingent on the fact that the islets follow the same network distribution. This is not necessarily the case when comparing control and Robo bKO islets.
The clustering coefficient was calculated as the number of connections for cell(i) divided by the number of possible connections. Because no constraints were placed on possible connections, this was the total number of cells in the islet -1. The average clustering coefficient presents the mean of C for all cells.
Global efficiency was calculated as follows. First, shortest path length was calculated with MAT-LAB function graphallshortestpaths(). This function uses the Johnston et al. algorithm (Johnston et al., 2016) to find the shortest path between every pair of cell in the islet. For example, the path length between cell(i) and cell(j) is 1 if they are directly synchronized, or 2 if cell(i) is not synchronized with cell(j), but each is synchronized with cell(k). Cells without any path connecting them were not considered. Finally, the characteristic path length (L) was calculated by summing up all nonzero path lengths normalized to total possible connections (size * (size-1)). The global efficiency is related to the inverse of global path length (Latora and Marchiori, 2001): In vitro single-cell Ca 2+ imaging Islets were isolated according to standard protocol from adult Robo bKO and control mice. For islet dispersion, 12 mm round No. 1.5 coverslips contained in a 24-well plate were pre-coated overnight with 50 mL 1:15000 PEI (Sigma P3143) overnight. Groups of 50-100 mouse islets were dispersed into single cells in 3 mL Accutase (Thermo Fisher A1110501) at 37˚C for 10 min. During the incubation, PEI was replaced with 100 mL Geltrex (Thermo Fisher A1413302) and centrifuged at 500 g for 5 min at 4˚C, followed by removal of excess Geltrex. The cells were washed once with islet culture medium (RPMI1640 supplemented with 10% FBS [v/v], 100 units/mL penicillin, and 100 mg/mL streptomycin [Invitrogen]) and resuspended in 1 mL medium before plating 500 mL per coverslip. The plate was centrifuged for 5 min at 500 g and cultured overnight before imaging. For measurements of cytosolic Ca 2+ , dispersed islet cells were pre-incubated in 5 mM Fura2-AM (Thermo Fisher F1201) in islet media containing 11.1 mM glucose for 45 min at 37˚C, followed by 15 min incubation in islet media containing 2.7 mM glucose. Coverslips were transferred to a RC-48LP imaging chamber (Warner Instruments) mounted on a Nikon Ti-Eclipse inverted microscope equipped with a 20Â/0.75N.A. SuperFluor objective and PerfectFocus (Nikon Instruments). The chamber was perfused with a standard external solution containing 135 mM NaCl, 4.8 mM KCl, 2.5 mM CaCl 2 , 1.2 mM MgCl 2 , 20 mM HEPES, and glucose as indicated (pH 7.35). The flow rate was set to 0.4 mL/min (Fluigent MCFS-EZ) and temperature was maintained at 33˚C using solution and chamber heaters (Warner Instruments). Excitation was provided by a SOLA SE II 365 (Lumencor) set to 10% output and an inline neutral density filter (Nikon ND4). Fluorescence emission was collected with a Hamamatsu ORCA-Flash4.0 V2 Digital CMOS camera at 0.1 Hz. Excitation (x) and emission (m) filters were used in combination with a ET FURA2/GFP C164605 dichroic (Chroma): Fura2, ET365/20x, ET535/30m; mCherry ET572/35x, and ET632/60m. b cells were identified by the expression of mCherry. Baseline-normalized [Ca 2+ ] i was quantified using Nikon Elements and GraphPad Prism software.

Statistical analysis
All analyses were performed in Prism (GraphPad) unless otherwise stated. For data with normal distribution (determined visually and by Shapiro-Wilk test) and equal variance (determined by F-test), a t-test was used. For normally distributed data with unequal variance between groups, a t-test with Welch's correction was used. For non-normally distributed data, a Mann-Whitney test was used to determine significance. p>0.05 was considered significant.