Open-source analysis and visualization of segmented vasculature datasets with VesselVio

Summary Vascular networks are fundamental components of biological systems. Quantitative analysis and observation of the features of these networks can improve our understanding of their roles in health and disease. Recent advancements in imaging technologies have enabled the generation of large-scale vasculature datasets, but barriers to analyzing these datasets remain. Modern analysis options are mainly limited to paid applications or open-source terminal-based software that requires programming knowledge with high learning curves. Here, we describe VesselVio, an open-source application developed to analyze and visualize pre-binarized vasculature datasets and pre-constructed vascular graphs. Vasculature datasets and graphs can be loaded with annotations and processed with custom parameters. Here, the program is tested on ground-truth datasets and is compared with current pipelines. The utility of VesselVio is demonstrated by the analysis of multiple formats of 2D and 3D datasets acquired with several imaging modalities, including annotated mouse whole-brain vasculature volumes.


INTRODUCTION
The acquisition of high-resolution and large-scale 3D vasculature datasets has been facilitated by recent and continual developments of powerful imaging technologies, including light-sheet fluorescence microscopy (LSM) (Di Giovanna et al., 2018;Kirst et al., 2020;Todorov et al., 2020) and micro-computed tomography (mCT) (Quintana et al., 2019;Schaad et al., 2017). Simultaneously, constant improvements in computational power and the widespread availability of powerful programming languages provide great means to process and extract detailed features from the resulting datasets. Use of micron-scaled imaging technologies and powerful computational tools to characterize the intricate details of microvascular networks will improve understanding of microvascular structure, function, and remodeling in health and disease.
Despite ongoing advancements in imaging technologies and open-source image segmentation software for 3D vasculature datasets (Kirst et al., 2020;Todorov et al., 2020;Haft-Javaherian et al., 2019), the publication of open-source and stand-alone analysis software applications that do not require prior programming experience is scarce. Existing standalone analysis tools that are accessible and freely available are unable to extract features from 3D datasets (Zudaire et al., 2011;Niemisto et al., 2005). Other current open-access, terminal-based analysis packages have similar 2D limitations (Mazzaferri et al., 2018;Montoya-Zegarra et al., 2019;Rust et al., 2020). Open-source analysis tools and software packages capable of analyzing 3D datasets are limited. Several of these tools depend on uncorrected centerline analyses, provide limited feature output, and markedly over-label branch points (Arganda-Carreras et al., 2010;Todorov et al., 2020). Other modern analysis packages that extract more accurate (Kirst et al., 2020;Tetteh et al., 2020) and detailed features (Hahn et al., 2019;Chapman et al., 2015) from 3D vascular networks require considerable programming skills or heavy interaction with terminals, potentially leading to unwelcome usage barriers and steep learning curves for researchers. Many publications also make use of private code or proprietary software for feature extractions (Jafarnejad et al., 2019;Epah et al., 2018;Kelch et al., 2015). In these MOTIVATION Modern imaging techniques and open-source image segmentation packages have led to the generation of large-scale and high-resolution vasculature datasets. However, open-source analysis packages for these datasets require programming knowledge, leading to an access barrier for some researchers. To overcome this barrier, we present VesselVio, a standalone and open-source application developed for the analysis and visualization of vasculature datasets.
instances, limited code availability and high software prices can hinder widespread analysis access. Thus, there is an apparent need for an easily accessible tool for 3D vasculature-dataset analysis.
Here, we present VesselVio, an open-source application for the analysis and visualization of vasculature datasets and vascular graphs. The back-end pipelines leverage custom feature-extraction techniques, high-level Python libraries, justin-time compilation, and parallel processing for rapid, detailed feature extraction and visualization of vasculature datasets. To make these analysis pipelines easily accessible, we also developed a downloadable (or, alternatively, single-line executable) front-end application. To test the performance and utility of VesselVio, we analyzed ground-truth synthetic vascular datasets, annotated mouse whole-brain datasets imaged using LSM, mouse inferior colliculus segmentations imaged using mCT, 2D retinography images, and pre-constructed mouse whole-brain vascular graphs.

RESULTS
Graphical user interface enables dataset loading, analysis, and visualization We sought to build an open-source application that allows users to extract and visualize numerous quantitative features from vascular networks (Figures S1A-S1C). VesselVio was developed for vascular datasets that have already been binarized, making it compatible with datasets of any imaging origin ( Figure 1A). The application is also compatible with pre-constructed graphs from other programs. Further, binarized volumes can be loaded alongside annotations, such as the p56 mouse brain atlas from the Allen Brain Institute, to analyze region-specific features ( Figure 1C). Centerline extraction and graph construction enable detailed network feature extraction The features of a vascular network can be extracted by identifying centerlines of the network and creating undirected graphs from the centerline points (Selle et al., 2002;Czech et al., 2011;Lee et al., 1994). Local neighborhood connectivity of the skeleton points is identified and used to create an undirected graph, G = (V, E). Graphs are constructed using the Python package igraph, which was selected because of its efficiency with graph loading and shortest-path finding, as well as its native compatibility with Windows and MacOS. Following initial graph construction, the datasets undergo a series of correction processes that eliminate spurious branchpoints, smooth centerlines, and remove isolated/endpoint segments at user-defined lengths ( Figure 1B). These initial processing stages enable downstream quantifications of network and segment features ( Figure 1C).

Spurious branchpoint clique filtering improves branchpoint quantification
After graph construction, network endpoints and branchpoints are identified by examining the degree of connectivity of the centerline vertices ( Figure 2A), a commonly used approach (Arganda-Carreras et al., 2010;Hahn et al., 2019;Todorov et al., 2020). Although centerline neighbor 26-connectivity can be used to accurately detect endpoints, this approach spuriously over-labels branch points, leading to artificially inflated counts ( Figure 2B). Some programs allow for interactive user input to correct mislabeled segments, inaccurate centerlines, or spurious branchpoints (e.g., Imaris [Quintana et al., 2019]). However, manually correcting branchpoints in 3D datasets can become tedious and time intensive, particularly when datasets are many gigabytes large. Presented here is an automated approach to correct these spuriously labeled branchpoints. Spuriously labeled branchpoints form small, all-to-all connected subgraphs in the constructed graphs, otherwise known as cliques. Previous techniques used to eliminate spurious branchpoint cliques rely on parent vessel orientations (Palá gyi et al., 2006) or generate candidate weights based on 26-neighborhood connectivity (Chen et al., 2009). VesselVio implements a two-pass filter that eliminates spurious branchpoint cliques through weighting based on the candidate vertex radius and the radii of neighboring vertices. Weighting candidates based on radius in addition to connectivity also mimics parent vessel hierarchy and improves segment radii calculations (Figure 2A).
To test our approach, a set of synthetically generated vasculature datasets with ground-truth branchpoint labels were acquired and analyzed (Tetteh et al., 2020), and results were compared with previous programs (Arganda-Carreras et al., 2010;Todorov et al., 2020). VesselVio clique filtering resulted in a mean 97.4% accuracy of branchpoint labeling (n = 20). Further, the hierarchical branchpoint filtering of VesselVio results in a mean 2.6% error that outperforms the mean 122.3% error of previous techniques ( Figure 2B).
Modified Euclidean distance calculations for segment radius estimation Several techniques exist for identifying vessel radii. One technique involves recording the largest maximally inscribed spheres that can rest within mesh vessel centerline points (Antiga et al., 2003;Antiga and Steinman, 2004), but this technique often depends on the creation of directed graphs (i.e., manually directed vessel hierarchy) and thus was not suitable for an automated pipeline. A similar method identifies the Euclidian distance between a vessel centerline and the center of the nearest non-vessel neighbor (Kirst et al., 2020;Mouches and Forkert, 2019). However, one pitfall to this approach is that vessels with near-resolution or at-resolution radii are incorrectly sized when their closest non-vessel neighbor is connected by a voxel face or a pixel edge ( Figure 2C), leading to oversized single-voxel/-pixel vessels ( Figure 2D). This issue is not as apparent for non-vessel neighbors connected to voxel edges/corners or pixel corners (Figures S1D-S1F). As such, a simple half-unit correction is implemented along specific orientations to preserve small-segment radii measurements for nonvessel neighbors connected to voxel faces and pixel edges ( Figures 2C-2F).
VesselVio provides detailed feature reporting Following the identification of branchpoints, endpoints, and centerline radii, individual vessel segments and their features can be isolated from vascular networks. Segment features extracted in the analysis pipeline include average length, radius, tortuosity, surface area, and volume. Next, this information is used to identify network features, including network volume and skeleton length, vessel segment partitioning, segment counts, and averages of segment features. Graph file exports from VesselVio enable users to independently identify relevant network characterization metrics, such as cohesion, network diameter, or clustering coefficient measurements (Hahn et al., 2019). Additionally, the centerline coordinate information extracted from the network are stored in the graph under a coordinates attribute, which can be accessed for custom geometric analyses. Lastly, by leveraging the same back-end analysis pipeline with an additional mesh-construction stage, VesselVio can render vascular datasets for visualization and inspection of results ( Figure 1C). 3D dataset compatibility and utility testing Analysis of annotated mouse cerebrovasculature imaged with LSM To demonstrate compatibility with 3D volumes, multiple datasets were analyzed. First, publicly available BALB/c mouse whole-brain vasculature datasets were acquired; these datasets were imaged with LSM (Todorov et al., 2020). To examine inter-regional network characteristics, 483 structures representative of the whole brain were selected using the Annotation Processing page ( Figure S1C). The brains were subsequently loaded and analyzed (n = 3; Figures Figure 3D). Further, an average of 3.5 million branchpoints were identified across the brain (Table S1A). Network surface area and length were observed to be largely Analysis of male and female mouse inferior colliculus segmentations imaged with mCT Next, cerebrovascular resin casts of female and male CFW mice (n = 5) were created, as described previously (Quintana et al., 2019). Following cast creation, imaging, and volume reconstruction, the inferior colliculi (IC) were segmented and analyzed (Video S1; Video S2). Network feature analysis revealed two primary sex differences in these datasets: female mice had greater average segment length and segment surface area compared with males ( Figures  . Furthermore, the runtime scaling of the pipeline was tested; processing speed was observed to scale linearly with voxel counts ( Figure S2H).

Compatibility testing with additional 3D datasets
Several other 3D vasculature datasets were also examined. First, an isotropic rodent lymph node microvascular network imaged with synchrotron mCT was visualized to reveal complex capillary networks (Video S3) (Jafarnejad et al., 2019). Next, the arterial and venous cerebrovasculature of an individual human brain from a public dataset was composited and visualized (Video S4) (Bernier et al., 2018). Finally, an anisotropic extraction of cerebrovasculature from C57 mice was examined (Video S5) (Di Giovanna et al., 2018). Together, the analysis and visualization of rodent cerebrovasculature, rodent lymph node, and human cerebrovasculature datasets demonstrate VesselVio's compatibility with 3D datasets generated with various imaging techniques.

Analysis of human 2D retinography datasets highlights disease-specific vascular phenotypes
To demonstrate the utility of VesselVio with 2D datasets, retinographs sourced from the HRF Image Database were analyzed and compared (Budai et al., 2013). This database contains high-resolution images of healthy control patients, patients with diabetic retinopathy, and patients with glaucomatous eyes. Analyses revealed several topological and segment-based differences among the groups. Relative to healthy controls, diabetic patients present with reductions in vessel-percentage-area fraction and surface area ( Figures 4A and 4C), reduced average segment radius ( Figure 4K), and increased vessel tortuosity (Figures 4L and 4O). Next, relative to healthy controls, glaucomatous retinas had reduced vascular surface area ( Figures 4A and C), increased branchpoint counts ( Figure 4D), increased segment counts ( Figure 4F), and reduced average segment length (Figure 4J). Finally, the mean radius of segments was reduced in both groups compared with healthy controls, as well as the distribution of segments across small-and medium-radius bins ( Figures 4K and 4M).

DISCUSSION
We developed an open-source application, VesselVio, for the analysis and visualization of segmented 2D and 3D vasculature datasets. Focus was centered on constructing an analysis pipeline capable of producing quantitative characterizations of numerous whole-network and small-scale segment features. Ground-truth result comparisons and three sets of analyses were conducted to demonstrate the accuracy of VesselVio for vasculature analysis.
The first two sets of analyses examined whole-brain vascular characteristics of BALB/c mice and sex-differences in the cerebrovasculature of female and male CFW mice. There are many sex differences in the cardiovascular system, including vascular tone (Boese et al., 2017), microcirculation (Huxley and Kemp, 2018), and blood-brain barrier function (Robison et al., 2019;Humphries et al., 2017). In the latter analysis, we sought to examine how the cerebrovascular network in a specific nucleus, the IC, may differ by sex. IC analyses revealed differences in the average length and surface area of segments, particularly large segments. The resulting whole-brain analyses provided continued insight into inter-regional differences in vascular density (Kirst et al., 2020;Ji et al., 2021). These results serve to demonstrate the 3D network analysis utility of VesselVio.
The third set of analyses examined retinal vascular networks of healthy controls, patients with diabetic retinopathy, and patients with glaucoma. Alterations in retinal microvascular networks are associated with vision loss (Carmeliet, 2005), and structural changes associated with excessive or reduced angiogenesis can also serve as indicators for underlying disease states (Lee et al., 2015;Selvam et al., 2018). Our analyses recapitulated previously reported reduced vessel area and diameter in glaucomatous eyes (Chan et al., 2017), as well as increased segment tortuosity in diabetic retinopathy eyes (Sasongko et al., 2011). Contrasting results in diabetic retinopathy eyes were noted in previous studies that reported increased vessel diameters (Bek, 2017;Klein et al., 2012), whereas we and others observed decreased diameters (Adhi et al., 2013;Fryczkowski et al., 1988). However, these differences may be due to the vessel type (the high-resolution fundus [HRF] vessels are not separated by arterial/venous hierarchy) or the vascular region analyzed. Although VesselVio is not intended to be used in any clinical diagnostic contexts, this analysis demonstrates the ability of this application to identify and characterize vascular network alterations in pathological conditions.
Future studies seeking to pinpoint detailed microvasculature alterations or broadly characterize regional vascular network features can employ the analytical and visualization components of VesselVio to obtain and enhance understanding of their results. In sum, VesselVio bridges the gap between modern vasculature imaging and analysis techniques by providing an optimized, open-source analysis pipeline and user-friendly application free for use.

Limitations of the study
In common with most analytical programs, the functionality of VesselVio is limited by the quality or resolution of the images that are loaded into the program. For example, if anisotropic datasets are loaded without pre-smoothing filters, then skeletonization of these datasets can produce erroneous segments. To adjust for these errors, the option to prune small, connected end segments is included (Video S6), but this pruning is unbiased and could hinder some types of analyses, such as angiogenesis assessments.
Rather than analyzing whole-brain datasets at once without subregion selection (e.g., brainstem analysis), in this pipeline it is optimal to segment the brain into smaller regions prior to analysis (e.g., interbrain, midbrain, and hindbrain analysis). This is because region-based analysis occurs after region segmentation. A limitation of this approach is that region annotation prior to skeletonization can lead to boundary effects, where vessels that cross annotation boundaries become disconnected during graph construction. The most notable effect of this disconnection would be altered endpoint counts. This approach was selected because of the offered speed optimizations and drastic memory usage reductions, as specific research questions often analyze subregional vasculature rather than whole-brain vasculature. In contrast, because of these boundary effects, it may be argued that VesselVio is most suitable for the analysis of smallto medium-sized datasets where subregion analysis is not necessary. Future approaches may seek to find a balance between this pipeline and other headless open-source pipelines, such as TubeMap (Kirst et al., 2020), that avoid boundary effects by requiring greater memory availability.
The presence of spurious branchpoints is unavoidable in the process of the implemented skeletonization algorithm and 26connectivity graph construction. When holes or imaging artifacts are present in datasets, the skeletonization algorithm can also produce erroneous centerlines. VesselVio takes an automated approach to correct these spurious branchpoints and erroneous centerlines. A limitation to the automated branchpoint-filtering algorithms may mean that in some instances, the automated correction of spurious branchpoints or erroneous centerlines may be less optimal than manual corrections. However, given that there can be hundreds to hundreds of thousands of spurious branchpoints in skeletonized datasets, manual corrections are not feasible. Future contributions may seek to implement userguided graph construction that can match the efficiency of our automated approach.
Because VesselVio requires that binarized datasets be loaded for analysis, our radius calculations are based on the binary volume and its centerline. Though our modified Euclidean distance transform (mEDT) radius calculations enable more detailed radius assessments than the traditional EDT calculations, this radius calculation method is likely not as effective as another recent approach that uses vessel-filled datasets, which are created with corrosion casting or fluorescent gel perfusion techniques. This approach adjusts uncentered centerlines and calculates the radius of centerline voxels based on an iterative assessment of vessel fluorescence intensity and angle (Ji et al., 2021). As such, future contributions to this application may seek to support radius calculations for vessel-filled datasets that are based on this described approach (Ji et al., 2021).
Finally, because VesselVio is a graphical user interface (GUI)based application primarily intended for standalone downloads, bundling the application for distribution acts as a limitation to speed. Bundling leads to slower performance in comparison to the non-bundled terminal-based version of our application or other headless analysis pipelines. Although functionality is equivalent between the bundled application and the application built from terminals at runtime, terminal builds offer more efficiency on startup and analysis. Because of this, single-line executable usage of the app is detailed in our documentation.

STAR+METHODS
Detailed methods are provided in the online version of this paper and include the following:

EXPERIMENTAL MODEL AND SUBJECT DETAILS Animals
All experiments were approved by the West Virginia University Institutional Animal Care and Use Committee, and animals were maintained in accordance with NIH Animal Welfare guidelines. Adult female and male CFW mice (7-to 8-weeks of age; strain #024; Charles River Laboratories, USA) were obtained and maintained under 14:10 light-dark cycles (150 ± 25:0 lux light:dark; lights on from 0500-1900 h). Following arrival, animals were given 1 week to acclimate to vivarium conditions before tissue collection. Food (2018 Teklad; Envigo, USA) and reverse osmosis water were provided ad libitum throughout the entire duration of the experiment.

METHOD DETAILS
Vascular corrosion casting and mCT imaging Vascular casts of mouse brains were created using a resin corrosion casting method described previously (Quintana et al., 2019). Perfusions occurred in the light phase between 1200-1600 h. Prior to perfusion, mice were injected i.p. with 25 U of heparin (63739-931-28; McKesson Corporation, USA) in 250 mL of saline and then deeply anesthetized with an isoflurane (4% induction, 1.5% maintenance) and 0.4 L/min oxygen flow mixture. Following confirmation of complete anesthetization with a pedal withdrawal reflex test, mice were transcardially perfused at a flow rate of 4 mL/min first with 15 mL of 25 U/mL heparin in saline, followed by 15 mL of 4% paraformaldehyde (#416780030; Acros Organics, Belgium) in saline with a pH of 7.4, followed lastly with PU4ii resin (VasQTec, Switzerland) formulated exactly as directed by the manufacturer with methyl ethyl ketone dilution (M209, Fisher Chemical, USA). Five days after perfusion, the craniums were decalcified with a 12-h wash of 5% formic acid (BDH4554; VWR International, USA), the brains were dissected, and remaining tissue was removed from the casts with two 12-h washes of 7.5% KOH (BDH7622; VWR International, USA) at 50 C. Casts were then rinsed with three 1-hour Milli-Q water washes, and the cleaned casts were osmicated in a 1% solution of osmium tetroxide (#75632; Sigma Aldrich, USA) for 12-hours to allow for optimal x-ray diffraction during mCT scans. Casts were imaged on a SkyScan 1272 (Bruker, USA) at 50 kV/ 200 mA with 360 rotations in step sizes of 0.17 , no filter, 900 ms frame exposures, and 4 frame averages/step to produce an isotropic voxel resolution of 2.7 mm 3 . Scan parameters were determined based on the manufacturer's guidance to achieve optimal x-ray transmission through the sample. Scans were then reconstructed using NRecon (Bruker) with beam hardening corrections at 15%, ring artefacts reduction at 3, smoothing at 0, custom alignment compensations set for each sample, and 0.02-0.40 dynamic image ranges. Following reconstruction, the volumes were resliced coronally for inferior colliculi segmentations. Using interpolated polygonal tracing with CTAn (Bruker), inferior colliculi were manually segmented from bregma À4.9 to À5.4 using the lobule 2 of the cerebellar vermis as a landmark for bregma À4.9 (Video S1) (Paxinos and Franklin, 2019). Scans, reconstructions, reslicing, and inferior colliculi segmentations were all conducted using the Bruker SkyScan analysis software suite. The resulting segmentatinos were analyzed at a 2.7 mm 3 isotropic resolution with a 10 mm isolated segment filter length and a 5mm endpoint segment prune length (Table S1B).

Dataset details
BalbC mouse brain vasculature datasets The BalbC whole-brain vasculature datasets were sourced from the VesSAP repository (Todorov et al., 2020). To analyze the whole brain datasets, 483 regions were selected using the Annotation Processing page. The downloaded brains were loaded directly into the program without alteration and were analyzed at a 3 mm 3 isotropic resolution with a 10 mm isolated segment filter length and a 5mm endpoint segment prune length (Table S1A).

HRF dataset
The HRF image dataset was downloaded directly from the HRF Image Database (Budai et al., 2013). The images were loaded directly into the program and analyzed with a 10 px isolated segment filter length and a 5 px endpoint segment prune length (Table  S1C). Synthetic dataset Twenty synthetic vasculature datasets with corresponding branch point label datasets were downloaded and holes in the vasculature were manually filled (Tetteh et al., 2020). The datasets were then loaded into the program and analyzed with a 10 voxel isolated segment filter length and a 5.5 voxel endpoint segment prune length. Branch point labels from the labeled datasets were quantified by skeletonizing the spheres around the branch points, labeling the components, and counting the components (Table S1E).

VesselVio pipeline
Dataset input and processing preparation Images are loaded into VesselVio using the Simple-ITK image reader or nibabel, depending on the filetype (Lowekamp et al., 2013;Brett et al., 2020). It is important to note that VesselVio is only compatible with vascular datasets that have been pre-segmented. To ensure that datasets loaded into the program are prepared appropriately for subsequent analysis, all loaded volumes are binarized and prepared as n-dimensional contiguous arrays of 0-value background and 1-value foreground elements. All array processing in VesselVio is conducted using NumPy and Numba (Harris et al., 2020;Lam et al., 2015). Volume skeletonization and centerline extraction VesselVio employs a custom implementation of a widely used medial axis parallel thinning algorithm to locate vessel centerlines from 2D and 3D datasets (Lee et al., 1994), with speed optimizations and parallel processing enabled by Numba. This algorithm was selected because it produces comparatively few erroneous branch point extensions, particularly when used with high-resolution datasets. This algorithm is also capable of thinning 2D and 3D datasets, making it optimal for this pipeline. Following skeletonization, (n, 3) or (n, 2) arrays are created containing Cartesian coordinate information of the location of all vessel centerlines.

Radii calculations
Before the creation of the undirected graph, vessel centerline radii measurements are conducted. Previous publications focused on voxel/pixel vasculature analysis have utilized Euclidean distance transforms (EDT) on segmented vasculature images to find centerline radii (Tsai et al., 2009;Todorov et al., 2020). These blanketed EDT methods find the Euclidean distance (ED) from the centerline point to the center of the nearest non-vessel neighbor point. However, an apparent unconsidered flaw in this traditional method is that finding the ED to the nearest non-vessel neighbor can overestimate the radius of the vessel if its nearest non-vessel neighbors are 6-connected in 3D space (N 6 (v)) or 4-connected in 2D space (N 4 (p)) ( Table S1D). For example, given an isotropic voxel resolution of 1 mm 3 , a straight vessel along the X-axis that is 1 voxel thick should have a putative radius of 0.5 mm ( Figure S1D). However, the traditional EDT method will record this vessel as having a radius of 1 mm 3 , doubling the apparent diameter. This is because ED measurements are calculated between the coordinates of the center of the centerline voxel and the center of its nearest non-vessel neighbor. Radii inflations in traditional EDT measurements are less pronounced as vessels become larger, but small-diameter vessels can be binned incorrectly using this blanketed technique.
To avoid this overestimation of vessel size, distance calculations can be made to the face of N 6 (v) or the edge of N 4 (p) non-vessel neighbors. To achieve this, a modified EDT (mEDT) calculation was implemented with 0.5-unit corrections for all N 6 (v) or N 4 (p) nonvessel neighbors (Figures S1D-S1F). When the nearest non-vessel neighbors are located beyond the immediate N 26 (v) or N 8 (p) space, the mEDT is applied to points that are connected to the centerline point only via face-to-face steps in 3D space or edgeto-edge steps in 2D space (i.e., not to diagonally connected points). In other words, given the phenotypical tubular shape of vessels, the mEDT approach is not applied to non-vessel neighbors of the sets N 26 (v) -N 6 (v) and N 8 (p) -N 4 (p) or non-vessel points beyond N 26 (v) or N 8 (p) neighborhoods that are situated diagonally to the centerline point ( Figures S1D-S1F).
In previous publications, EDTs are applied to the entire image, and the centerline values are then extracted. However, this process is computationally expensive and leads to difficulty with making the intended 1D corrections. As such, the mEDT method is only applied to centerline points, rather than the entire image. Points in the previously constructed centerline coordinate arrays are used as seeds for the placement of expanding boxes that search centerline point neighborhoods for non-vessel neighbors. This box expands around a centerline point until at least four non-vessel neighbors in the binarized volume are identified, rather than a single neighbor as with previous implementations. Then, rather than directly computing the ED (square root of sum of squared deltas) between the centerline point and neighbor (1), the absolute deltas of the centerline and neighbor coordinates are computed and loaded into a lookup table with precomputed mEDT values (2). This lookup table is constructed based on the input resolution to provide the half unit corrections along the described orientations (e.g., an XYZ delta array of [0,0,5] returns a 4.5 voxel distance).
After the distance for each non-vessel neighbor is calculated, the four lowest values are averaged, and this average is defined as the radius for the centerline point. Four values are averaged to improve the accuracy of radii measurements and account for potential small vessel surface divots, bubbles, or imaging artifacts common in various vascular preparation techniques, such as corrosion casting. This process is parallelly repeated for all centerline points. Graph construction Following radii measurements, undirected graphs are constructed to represent the vasculature skeleton. Graph creation and processing in VesselVio uses the Python igraph package (Csardi and Nepusz, 2006). To construct the vascular graphs, the number of centerline points is identified, and an equal number of isolated vertices are added to the graph. At this point, 2D images are padded as 3D arrays. Centerline points are assigned an index based on their order of appearance in the point coordinate array. Each vertex is given a coordinate and radius attribute based on its corresponding point values. Next, the [x max , y max , z max ] oriented 13-connectivity neighborhood of all centerline points is scanned to identify N 26 (v) neighbors. Upon neighbor identification, an edge is created between the vertices. Unidirectional scanning in the edge identification process prevents parallel edges from being created (Kirst et al., 2020

Mesh visualization and application interface
To construct meshes for visualization, we leverage the high-level Python package PyVista (Sullivan and Kaszynski, 2019) that wraps The Visualization Toolkit. We create individual polydata datasets from our segment splines, apply tube filters to create centerline and scaled network meshes, and assign each tube a radius, length, tortuosity, surface area, and volume scalar for visualization. All scalars and scaled segment sizes are based on the mean of segment features. These segments are combined into an undirected grid for surface extraction and subsequent rendering. Additional branch point and endpoint meshes are created. Then, original (voxel/pixel based) and smoothed surface meshes (marching cubes based-based) are created from the filtered input volume so researchers can visually validate the output features from the program by comparing simple/scaled networks to the original/smoothed meshes. Lastly, construction of the front-end application for VesselVio was accomplished with PyQt5 under GNU GPLv3 licensing (Figures S1A-S1C).

Computational resources
VesselVio was developed with Python 3.8 and was tested locally on a 2019 16 00 MacBook Pro (2.6 GHz 6-Core Intel Core i7, Intel UHD Graphics 630 16GB) and a Windows 10 computer with dual 3.39 GHz Intel Xeon Gold 6128 processors, 128 GB of Sk Hyinx Ddr4 Sdram Memory Module RAM, and dual 8GB NVIDIA Quadro P4000 graphics cards.

QUANTIFICATION AND STATISTICAL ANALYSIS
All statistical analyses were conducted using Prism 9 (GraphPad; USA). Synthetic vasculature results were analyzed using one-way ANOVA tests and repeated-measures two-way ANOVA tests. HRF datasets were analyzed using one-way ANOVA tests. Mice IC data were analyzed using two-tailed student's t tests. Distributions of segment counts per radii bin, average segment length per radii bin, and segment tortuosity per radii bin were all analyzed using two-way ANOVA tests. Following main effects observations in one-or two-way ANOVA tests, multiple comparisons were made using Sidak's test. P-values below 0.05 were considered statistically significant.