A method for three-dimensional quantification of vascular smooth muscle orientation: application in viable murine carotid arteries

When studying in vivo arterial mechanical behaviour using constitutive models, smooth muscle cells (SMCs) should be considered, while they play an important role in regulating arterial vessel tone. Current constitutive models assume a strictly circumferential SMC orientation, without any dispersion. We hypothesised that SMC orientation would show considerable dispersion in three dimensions and that helical dispersion would be greater than transversal dispersion. To test these hypotheses, we developed a method to quantify the 3D orientation of arterial SMCs. Fluorescently labelled SMC nuclei of left and right carotid arteries of ten mice were imaged using two-photon laser scanning microscopy. Arteries were imaged at a range of luminal pressures. 3D image processing was used to identify individual nuclei and their orientations. SMCs showed to be arranged in two distinct layers. Orientations were quantified by fitting a Bingham distribution to the observed orientations. As hypothesised, orientation dispersion was much larger helically than transversally. With increasing luminal pressure, transversal dispersion decreased significantly, whereas helical dispersion remained unaltered. Additionally, SMC orientations showed a statistically significant (p<0.05\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$p < 0.05$$\end{document}) mean right-handed helix angle in both left and right arteries and in both layers, which is a relevant finding from a developmental biology perspective. In conclusion, vascular SMC orientation (1) can be quantified in 3D; (2) shows considerable dispersion, predominantly in the helical direction; and (3) has a distinct right-handed helical component in both left and right carotid arteries. The obtained quantitative distribution data are instrumental for constitutive modelling of the artery wall and illustrate the merit of our method.


Introduction
Smooth muscle cells (SMCs) play a crucial role in regulating arterial vessel tone. When an SMC contracts, it exerts a force along its long axis. Therefore, the orientation of SMCs within the artery wall is mechanically of importance. Several mechanical models of the artery wall include a smooth muscle component, e.g. Masson et al. (2011), Spronck et al. (2015), Zulliger et al. (2004). To our knowledge, all currently available constitutive models assume SMC orientation to be strictly circumferential and ignore any dispersion in SMC orientation. In order to develop constitutive models that describe SMC orientation more realistically, knowledge on SMC orientation and its dispersion is essential.
The assumption in current constitutive models of strictly circumferentially oriented SMCs may have arisen from experimental studies (Peters et al. 1983;Walmsley 1983), stating that SMC orientation is circumferential. Other studies have shown a more disperse SMC orientation with two main orientations (Holzapfel et al. 2002) and have shown that SMCs are enveloped by collagen bundles, leading to a parallel orientation of collagen and SMCs in the media of rat aortas (O'Connell et al. 2008). In the aforementioned studies on SMC orientation, arteries were first fixed and subsequently histologically sectioned, potentially causing artefacts. To overcome this problem, a method using intact arteries is required.
Previously, we developed a two-dimensional (2D) method to quantify SMC orientation in vessels using two-photon laser scanning microscopy (TPLSM) (Spronck et al. 2014b). This method, which was used in rat ureters, had the advantage over previous methods that it did not require sectioning. However, it only yielded 2D orientations, neglecting out-of-plane orientation. Clearly, using such a method, it is impossible to quantify transverse angles. In addition, due to the curvature of the wall, only a very limited part of the wall could be analysed, leading to the inclusion of only a few SMCs per acquisition. To parameterise a constitutive model, ideally, three-dimensional (3D) orientation and dispersion information are available.
We hypothesised that SMC orientation would show considerable dispersion, based on recent published work assessing porcine aortic SMC orientation in 2D (Tonar et al. 2015). As we expected SMCs to be mainly oriented in the plane of the arterial wall (Peters et al. 1983), we further hypothesised SMC orientation dispersion in the helical direction to be greater than in the transversal direction. To test these hypotheses, we developed a method for the quantification of 3D SMC orientation, applicable to viable, intact murine carotid arteries. This 3D method does allow for the quantification of transverse angles and is suited for the analysis of a much larger part of the acquired vessel wall than the previous 2D method (Spronck et al. 2014a, b). In the present study, we used this 3D method to assess orientation and dispersion differences between left and right carotid arteries (potentially relevant from a developmental biology point of view). Additionally, we investigated the effect of increasing luminal pressure on 3D SMC orientation.

Animals and staining
Mice were killed with an overdose of isoflurane (Forane, Baxter, Deerfield, IL, USA). Left and right common carotid arteries of ten male C57BL/6JRj mice were excised, labelled at the proximal end, and carefully mounted between micropipettes (Fig. 1a) to maintain viability for a prolonged period of time (>2 h) as described previously (Megens et al. 2007). Carotid artery sections of 5 mm were mounted, yield- b Example x-y image slices acquired at increasing imaging depth, i.e. at decreasing z-coordinate (1: z = 93 µm, 2: z = 81 µm). Scale bar 200µm. c Examples of reconstructed x-z image slices. Figure scale equal to scale in b ing a length of 3 mm of artery usable for imaging (i.e. without imaging one of the pipettes). Calcium-free HEPES buffer was used as submersion medium. SMC nuclei were stained using 2.7 µM SYTO 13 (Life Technologies , Gaithersburg, MD, USA) for 30 min. Arteries were maximally vasodilated by the addition of nitroprusside (Sigma-Aldrich, St. Louis, MO, USA) to a concentration of 10 µM and were imaged at room temperature. Animals were checked for situs inversus totalis (SIT) by the inspection of the location of the visceral organs; none of the animals showed SIT. All animal experiments were approved by the local authority.

Image acquisition
Viable mounted arteries were imaged using an upright TPLSM system (Leica TCS SP5II MP, Leica Microsystems , Mannheim, Germany) equipped with a Leica HCX APO L 20×/1.00 W water dipping objective (Fig. 1a). Two-photon excitation was achieved using a pulsed Ti-Sapphire laser (Spectra Physics Mai Tai DeepSee, Newport, USA) tuned at 825 nm. Signal detection was performed using a descanned hybrid diode detector (GaAsP) between wavelengths of 510 and 555 nm. Image slices of 1480 × 1480 pixels were acquired with a pixel size of 0.5×0.5 µm 2 , resulting in a field of view of 0.74×0.74 mm 2 . Pixel dwell time was 1.1 µs, and acquisition was performed at 12-bit precision. Acquisition of one slice took 2.6 s. Slice spacing was 0.5 µm, resulting in an effective voxel size of 0.5 × 0.5 × 0.5 µm 3 . All subsequent image processing (see following section) was performed at this resolution; no downsampling was performed. Image slices formed stacks of 334 ± 30 slices (mean ± SD). Left carotid arteries were imaged at luminal pressures of 40, 80, and 100 mmHg, respectively. Right carotid arteries were imaged at 40 mmHg. At each luminal pressure, one image stack was acquired for each animal. Representative image slices are shown in Fig.  1b,c. All subsequent image processing was performed at full resolution, i.e. no downsampling was performed. All arteries were imaged directly after excision. A maximum of 2 h was required to complete the imaging of the left and right arteries of one mouse. Previously, we have shown that in this time span, mounted arteries remain viable (Megens et al. 2007).

Deconvolution and vesselness filtering
Raw image stacks were deconvolved using 10-pass adaptive blind 3D deconvolution (AutoQuant X2, MediaCybernetics , Bethesda, MD, USA). Deconvolution of one stack took approximately 35 min on a modern computer [Intel i7-2670QM quad-core central processing unit . Further data processing and analysis were performed using MATLAB R2014b (MathWorks , Natick, MA). After deconvolution, the deconvolved data stack was 3D vesselness filtered (Frangi et al. 1998) to enhance elongated structures. Note that in contrast to the approach described in our previous paper (Spronck et al. 2014b), filtering is performed in 3D on the image stack as a whole. A 3 × 3 Hessian matrix was calculated for each voxel using a Gaussian width of σ x = σ y = σ z = 1 µm. Hessian components were calculated by convolving the image stack with discretised Gaussian derivative kernels. The eigenvalues of the Hessian matrix (λ 1 ≤ λ 2 ≤ λ 3 ) were used for performing the vesselness filtering. Vesselness parameters were α = β = 0.5 (Frangi et al. 1998). c was determined for each image stack as half the value of the maximum Hessian norm, as suggested by Frangi et al. (1998). Vesselness filtering of one stack took approximately 65 min on a modern computer (see above). All further image processing steps were completed in <1 min per stack.

Clustering and cell orientation calculation
To separate cell nuclei from background, the vesselness stack was thresholded at a cut-off value of 0.01. Voxels were then clustered based on their 3D 6-connected neighbourhood (Gonzalez and Woods 2008). Clusters were included based on their volume using upper (320.1 µm 3 ) and lower (38.5 µm 3 ) thresholds calculated as mean ± 2SD from aortic SMC nucleus sizes (O'Connell et al. 2008). An example of a stack of clustered nuclei is given in Fig. 2. For each cluster, an inertia matrix was calculated (Jähne 1993;Vader et al. 2009). The eigenvector corresponding to the largest eigenvalue of this matrix represents the principal SMC orientation. Coordinate and orientation definitions. a Definition of cylindrical coordinate system using radial (ρ), angular (θ), and axial (z c ) coordinates. b, c Orientations of SMCs can be visualised as points on a unit hemisphere. b The orientation of such a hemisphere with respect to the vessel. c The two spherical coordinates that are used to quantify SMC orientation. θ h and θ t are an SMC's helix and transverse angles, respectively. d Exact definitions of θ h and θ t . θ h is the angle that the orientation's projection on the axial-circumferential plane makes with the circumferential direction. θ t is the angle that the orientation makes with the axial-circumferential plane. SMC smooth muscle cell

Coordinate transformation
To assess SMC orientations with respect to the vessel wall, a cylinder was least-squares fitted through the detected centroids as follows: first, a rough fit was obtained by fitting a cylinder through all centroids. Using the axis of this cylinder, the inner layer of centroids was detected by selecting the centroids closest to the axis in a moving 3D search window of 50 µm in axial and 10 • in angular directions. Second, a cylinder was fitted through the detected inner layer of centroids. The axis of this cylinder was used as a reference to convert all centroid locations to cylindrical coordinates (ρ, θ, z c , Fig. 3a). Subscript c in z c is used to distinguish the axial cylindrical coordinate (z c , Fig. 3a) from the vertical Cartesian coordinate (z, Figs. 1, 2). The orientations of the SMCs were expressed locally with respect to the vessel wall (Fig. 3b) and are represented as a point on a unit hemisphere (Fig. 3c). Orientations were quantified by a helix angle (−90 • ≤ θ h < 90 • ) and a transverse angle Fig. 3c, d). Visual inspection of the detected cell nuclei (as exemplified in Fig. 2) shows that, from a certain imaging depth, the number of detected nuclei decreased. In order to assure that we imaged a region of sufficient image quality, centroids were selected from an angular region of interest of −22.5 • ≤ θ ≤ 22.5 • . In order to prevent border effects, an axial window slightly smaller than the z crange was manually selected (704 ± 23 µm, mean ± SD). In six stacks, using this large axial window led to SMC layers getting merged during SMC layer detection (see following section) due to inhomogeneities in vessel diameter along the axial coordinate. In these cases, a smaller window of 379 ± 74 µm (mean ± SD) was used.

Smooth muscle layer separation
Cell density as a function of the radial coordinate (ρ) was calculated using kernel smoothing density estimation (KDE) (Silverman 1986). KDE allows for estimation and visualisation of the probability density from a set of data points. In our case, the KDE shows the density of SMCs at a certain radial coordinate (ρ). A Gaussian kernel with σ = 1 µm was used on the list of radial coordinates (ρ) of all included nuclei of one stack, yielding the ini- Fig. 4 Representative example of SMC density as a function of the radial coordinate (ρ, Fig. 3a). Inner and outer layers of SMCs are clearly visible. Layers were detected by finding the two largest density maxima. The distance between those maxima was defined as 2 . SMCs are assumed to be part of the inner/outer layer if they are within a distance of from their respective maximum. SMC smooth muscle cell tial kernel density estimate f (ρ) [µm −1 ]. This estimate was converted to an absolute density d(ρ) in cells per µm 3 via with N the total number of detected nuclei, r θ the angular region of interest in radians, and r z c the axial region of interest in µm. An example density plot is shown in Fig. 4a. SMCs showed to be concentrated in two distinct layers, as visible in the example. We analysed SMC orientation separately for these (inner and outer) layers.
To separate the layers, we detected the two density maxima, the distance between which we defined as 2 . Each layer is taken to be centred around its maximum, with boundaries of ± on each side (Fig. 4a). The artery diameter (as reported in "Results" section) is taken to be the middle between the inner and outer layer density maxima.

Orientations of SMC nuclei in detected layers
After smooth muscle layer separation, SMC orientations can be depicted on a projected hemisphere, as defined in Fig. 3c. A representative example of such a visualisation is given in Fig. 5a-c. The need for using a projection is elaborated on in "Appendix 1".

Orientation statistics
As our method yields orientations in 3D, appropriate, spherical statistics are required to accurately quantify our results (Fisher et al. 1993). In particular, the obtained orientations are undirected (i.e. axial data), requiring a 3D probability distribution that is antipodally symmetric (Mardia and Jupp 2000). One of the simplest distributions of this kind would be a (rotationally symmetric) Watson distribution (Mardia and Jupp 2000). However, this distribution requires Orientations of detected nuclei, visualised using Lambert equal-area projection a before layer separation, b for the inner layer, and c for the outer layer. The Lambert equal-area projection is elaborated on in "Appendix 1". d-f Corresponding fitted Bingham distributions to orientations in a-c. The obtained Bingham distributions provide a quantification of the distribution of the individual orientations. The distribution parameters are explained in Fig. 6. SMC smooth muscle cell data to be rotationally symmetric around the main direction. In the vessel wall, the largest spread in SMC direction will be in the axial-circumferential plane, whereas the spread in the radial-circumferential plane will be much smaller. The probability distribution used should be able to accommodate this non-rotationally symmetric geometry. The Bingham distribution fulfills this requirement, while still being antipodally symmetric (Bingham 1964(Bingham , 1974. The Bingham density function has the form (Bingham 1964): : where κ 1 and κ 2 are two concentration parameters, μ 1 :: and μ 2 :: are two orthogonal unit vectors indicating the two spread directions of the distribution (corresponding to κ 1 and κ 2 ), and x : the current local (Cartesian) coordinate on the unit sphere. The normalisation function d(κ 1 , κ 2 ) was evaluated by numerically solving the integral (Bingham 1964;Onstott 1980) d(κ 1 , κ 2 ) = 2π 0 π 0 exp κ 1 cos 2 φ + κ 2 sin 2 φ sin 2 θ sin θ dθ dφ. (3) Fitting the Bingham distribution to our data was performed in two steps by using the moment method (Bingham 1964;Onstott 1980;Tanaka 1999). First, the local centroid orientations were converted to Cartesian coordinates and the following matrix was computed: with N the number of detected centroids. μ 1 :: and μ 2 :: can now be obtained as the eigenvectors of X X T corresponding to the smallest and middle eigenvalues, respectively. μ 3 :: , corresponding to the largest eigenvalue, represents the distribution's centre. It is converted back to spherical coordinates (Θ h , Θ t , Fig. 6a) using Parameters describing the distribution's dispersion. κ 1 and κ 2 describe the distribution's concentration along two orthogonal directions. The more negative a κ is, the narrower the distribution is in the respective direction (Onstott 1980). Θ i is the angle that the principal spread direction (corresponding to κ 2 ) makes with the axial-circumferential plane. In this example, where μ 3,x , μ 3,y , and μ 3,z represent the x, y, and z components of μ 3 :: . Note the use of capital thetas, representing the orientation's centre. Θ i is defined as the angle that μ 2 :: makes with the x y plane (Fig. 6b). The second step in fitting the Bingham distribution is finding κ 1 and κ 2 . This is performed through maximisation of the log-likelihood function (Bingham 1964;Onstott 1980;Tanaka 1999) where τ 1 and τ 2 are the smallest and middle eigenvalues of X X T . Maximisation was performed by the minimisation of −F using the MATLAB Optimization Toolbox fmincon function and the interior-point algorithm, while constraining κ 1 and κ 2 to the interval (−∞, 0). Bingham distributions of the representative examples in Fig. 5a-c are shown in Fig. 5d-f.

Qualitative results
Orientation plots for all arteries studied show a cloud of points that is approximately centred around the circumferential direction, as also visible in the representative example . Θ h , Θ t , Θ i , κ 1 , and κ 2 are defined in Fig. 6. SD standard deviation, SMC smooth muscle cell ( Fig. 5a-c). Although these plots provide a good qualitative overview of the orientations of the SMC nuclei, it is hard to derive quantitative conclusions from them. Therefore, Bingham distributions were fitted to the acquired orientations, as exemplified in Fig. 5d-f, providing quantitative orientation data. These quantitative (Bingham parameter) results are discussed below. Figure 7 shows Bingham parameters and diameters for left (n = 10) and right (n = 10) carotid arteries at a luminal pressure of 40 mmHg. Left and right artery diameters are nearly equal. As Θ i is relatively close to zero, κ 1 represents dispersion in transversal direction and κ 2 represents dispersion in helical direction. The fact that κ 1 is more negative than κ 2 indicates that dispersion in helical direction is larger than in transversal direction. This holds for left as well as right arteries. Θ h is positively different from zero in the right as well as the left arteries, in both the inner and outer layers (all p < 0.05), which is indicative of a right-handed helical pattern. No statistically significant differences were found between left and right arteries, except that only in the right arteries, Θ h was significantly larger in the outer layer than in the inner layer. Θ t shows the same pattern, being significantly different from zero in all cases. Please refer to "Discussion" section for an interpretation of the latter finding.

Dependence of smooth muscle orientation on luminal pressure
As expected, artery diameter increased significantly with increasing luminal pressure, illustrating vessel deformation at macro-level (Fig. 8). Figure 8 further shows Bingham parameters for the left carotid arteries (n = 10) at three luminal pressures (40, 80, and 100 mmHg). Θ h and Θ t are both again significantly different from zero in all cases. κ 1 shows a very consistent pattern in the inner layer, becoming more negative with increasing pressure, thereby indicating a decrease in SMC orientation dispersion. As mentioned above, because Θ i is relatively close to zero, this decrease is approximately along the transversal direction. κ 2 increases with increasing luminal pressure from 40 to 100 mmHg in the inner layer, reflecting an increase in dispersion in helical direction. κ 1 and κ 2 show the same pattern with increasing luminal pressure in the outer layer, although it is not statistically significant here.

Reproducibility
Reproducibility of our findings was assessed and is described in "Appendix 2".

Discussion
In this study, we quantified SMC orientation and dispersion in viable murine carotid arteries in full 3D. We found that SMCs are arranged in two layers, and we quantified SMC orientation and dispersion for each layer separately. Orientation dispersion in the helical direction was larger than in the transversal direction. We found a statistically significant right-handed helix in both left and right arteries, in both layers. With an increase in luminal pressure, the dispersion in transversal orientation was found to decrease. The quantification method we developed in this study has three key advantages over previous methods: (1) arteries are mounted in their in vivo geometry and are kept viable (Megens et al. 2007); (2) by using TPLSM, 3D microscopic image stacks at large imaging depths (as compared to e.g. confocal laser scanning microscopy) can be acquired (Oheim et al. 2001); and (3) the present processing method analyses , and 100 mmHg (circular mean ± circular SD for Θ h , Θ t , and Θ i ; arithmetic mean ± SD for κ 1 , κ 2 , and diameter). * p < 0.05 versus angle zero (Rayleigh test). # p < 0.05 inner versus outer layer (Rayleigh or paired Student's t test). + p < 0.05 for parameter difference between pressures (Rayleigh or paired Student's t test). Θ h , Θ t , Θ i , κ 1 , and κ 2 are defined in Fig. 6. SD standard deviation, SMC smooth muscle cell the full 3D image stack as a whole, instead of analysing data on a slice-by-slice basis as most other methods do [including our previous method (Spronck et al. 2014b)]. As arteries are mounted in their in vivo geometry (point 1), artefacts that could potentially arise due to fixation (e.g. shrinking) and/or sectioning are absent. It has been shown that the present vessel mounting and TPLSM imaging set-up (points 1,2) can be used to visualise both structural and functional properties of intact viable murine arteries. Its physiological relevance is demonstrated in various studies. Examples include the study of the endothelial glycocalyx (Reitsma et al. 2011;Slaaf and Zandvoort 2011), endothelial junctional adhesion molecule A (JAM-A) and its effect on monocyte or platelet recruitment (Karshovska et al. 2014;Schmitt et al. 2014), or the presence of smooth muscle progenitor cells (Subramanian et al. 2010) and proliferating endothelial cells ) in healthy and atherosclerosis-prone arteries. Consequently, we believe that the results obtained in our study are relevant and well translatable to an in vivo situation, albeit that we studied mouse and not human arteries. In a sliceby-slice analysis, flat image slices are acquired of a curved object, causing cells at various depths of the wall to end up in one image slice (Spronck et al. 2014b). This limitation is overcome in our present approach by analysing the acquired image stacks as a whole and using a curved (cylindrical) coordinate system (point 3), omitting the need for a small region of interest (Spronck et al. 2014b). Hence, our 3D method allows a larger number of nuclei to be analysed and thus provides more robust results. Furthermore, transverse angles and dispersion therein can be analysed, whereas in a slice-by-slice analysis, out-of-plane information is lost.
Our results display statistically significant mean transverse and helix angles (Θ t and Θ h ,Figs. 7,8). SMC orientation shows considerable dispersion in both helical and transversal directions. To critically evaluate our findings, we performed the following two additional experiments for five left arteries at a luminal pressure of 40 mmHg: (1) we digitally rotated our image stack (after deconvolution) by 180 • around the z-axis.
(2) We physically rotated the mounting chamber by 180 • around the z-axis (i.e. around the objective axis). As compared to the non-rotated situation, both experiments should lead to a negation of Θ t and Θ i . Experiment 1 indeed led to these negations (data not shown) and did not affect Θ h , κ 1 , and κ 2 . Experiment 2, on the contrary, still yielded positive Θ t values (Fig. 9). The observed significant average transverse angle (Θ t ) therefore seems to be a methodological artefact and not a true anatomical feature. We believe that this artefact results from inhomogeneities in image field brightness, slightly shifting the fitted cylinder axis and causing an artificially significant Θ t value. The mean helix angle (Θ h ) observed in these two experiments remained positive, as would be theoretically expected. For the physically rotated acquisitions, we used the same imaging region as used in these arteries for the non-rotated acquisitions. We achieved this by measuring the distance between the nonrotated imaging region and the proximal pipette (measuring the movement of the calibrated microscope table). After the rotation of the mounting chamber, imaging was performed at the same distance from the proximal pipette.
We further verified our helix angle findings by a third experiment, digitally flipping image stacks (after deconvolution) along the x-axis. Results were as expected: a negation of Θ h and Θ t , and Θ i and both κ parameters remaining unaltered. Therefore, it is unlikely that our statistically significant average helix angle findings are artefactual.
Historically, SMC orientation has been studied in a variety of different species and arteries. Some authors described a spiralling orientation of medial smooth muscle (Benning- Fig. 9 SMC orientations and arterial diameters in normal and rotated left carotid arteries at a luminal pressure of 40 mmHg (circular mean ± circular SD for Θ h , Θ t , and Θ i ; arithmetic mean ± SD for κ 1 , κ 2 , and diameter). * p < 0.05 versus angle zero (Rayleigh test). Θ h , Θ t , Θ i , κ 1 , and κ 2 are defined in Fig. 6. SD standard deviation, SMC smooth muscle cell hoff 1927; Pflieger and Goerttler 1970;Schultze-Jena 1939;Strong 1938), whereas others described smooth muscle orientation as circumferential (Pichler et al. 1953;Ushiwata and Ushiki 1990). Benninghoff described already in 1927, a spiralling smooth muscle orientation with a helix angle increase towards the periphery (Benninghoff 1927). In 1939, Schultze-Jena assessed the muscle spiral in a variety of arteries and found opposing-handed helices in left and right arm arteries (Schultze-Jena 1939). A subsequent work by Pichler et al. (1953) assessed the potential mechanical consequences of a helical orientation. Pichler concludes, based on measurements in helically cut strips of carotid artery, cut at helix angles of +30 • and −30 • that a single-helical structure is unlikely. Although these historical studies provide excellent qualitative anatomical insight, only limited quantitative conclusions can be drawn.
The results we obtained, indicating a significant righthanded helix, contradict the conclusion of most quantitative studies that SMC orientation is strictly circumferential. In human cerebral arteries, Walmsley (1983) found a statistically non-significant mean helix angle of −0.2 ± 2.4 • , in contrast to our angles of 2.4 ± 2.2 • (inner layer) and 5.0 ± 5.0 • (outer layer; mean ± SD; left carotid arteries at 40 mmHg). Peters et al. reported a "truly circumferential" smooth muscle orientation and an "inconsequential" left-or right-handedness of the helix (Peters et al. 1983) in human cerebral arteries, whereas we found a consistent, significant helix, however, in viable carotid arteries of mice. Some studies take a slightly different approach and consider the arterial (aortic) smooth muscle orientation to be a mixture of two or more distributions. Holzapfel et al., in their modelling paper in 2002, analysed SMC orientation in human aortas (Holzapfel et al. 2002). They implemented a method similar to our previous 2D method (Spronck et al. 2014b), but used histologically fixed instead of fresh tissue. Interestingly, Holzapfel et al. assumed the SMC orientation distribution to be centred around the circumferential direction, with two equally weighted helical parts. Using this assumption, they argued for two helices with angles of ±8.4 • . Tonar et al. (2015) assessed potential segmental differences in SMC orientation in porcine aortas, using a mixture of one to five Von Mises distributions. Their paper concludes that "The orientation of vascular smooth muscle was successfully fitted using two von Mises distributions in most of the samples…" and that "Only a minor fraction of samples required a tertiary von Mises component to describe the orientation…". We carefully studied Supplemental Material 1 of their paper. When assessing the N c = 2 columns (i.e. the columns that correspond to two Von Mises distributions), one can count that in 55 out of 82 samples, one of the two Von Mises components was weak and/or redundant. Tonar et al. (2015) state this in section 4.2 as "From a mathematical modelling point of view, this indicates that a structural model with a single vascular SMC system would be sufficient for most of our data". Another finding of Tonar et al. is a segmental difference in SMC orientation dispersion (not in main orientation) along the aortic tree. This may also be the case along the carotid artery, potentially explaining the spread in κ 1 and κ 2 parameters we measured (whiskers in Figs. 7 and 8), since we may not have imaged the exact same part of the common carotid artery in all mice.
In our study, we analysed both left and right arteries and found a right-handed helix at both sites, violating leftright symmetry. Although very speculative, this may be explained by the intrinsic asymmetry of molecules, especially the proteins actin and myosin (Delhaas et al. 2008). As both molecules are arranged in a right-handed helix, when stretched, they will exhibit torsional motion. If a larger structure is built from such smaller helical components, an overall helicity may result. A more detailed discussion of Fig. 10 SMC orientations and arterial diameters in left carotid arteries at a luminal pressure of 40 mmHg in young (age 8 weeks) and old (age 23 weeks) mice (circular mean ± circular SD for Θ h , Θ t , and Θ i ; arithmetic mean ± SD for κ 1 , κ 2 , and diameter). * p < 0.05 versus angle zero (Rayleigh test). Θ h , Θ t , Θ i , κ 1 , and κ 2 are defined in Fig. 6. SD standard deviation, SMC smooth muscle cell this hypothesis can be found in Delhaas et al. (2008). In our case, when this hypothesis would hold, one would expect left and right arteries to show the same helical handedness, as we observed. It may also be the case that this molecular asymmetry effect acts through the extracellular matrix. As smooth muscle cells have been shown to align with the extracellular matrix collagen (O'Connell et al. 2008), asymmetry of the extracellular matrix will presumably translate to asymmetry in SMC orientations as well.
Our study population of ten male C57BL/6JRj mice consisted of five young (age 8 weeks) and five older (age 23 weeks) individuals. We assessed whether age had an influence on orientations by comparing the results for the two age groups for left arteries at 40 mmHg (Fig. 10) and found no statistically significant orientation parameter difference. We therefore decided to pool data from both age groups, yielding one group of n = 10.
Our study involved the use of lower and upper volume thresholds to exclude clusters that were too large or too small. The effect of this thresholding was analysed in a subset (n = 5) of samples at 40 mmHg (Fig. 11). Thresholding excluded 26.8±3.9 % (mean ± SD) of the clusters. The locations of these clusters are shown in Fig. 11. The combination of vesselness filtering and subsequent cluster volume thresholding was aimed at excluding non-SMCs. It should be noted that the number of nuclei excluded by vesselness filtering is not visible in Fig. 11. Calculating this number would require running the image analysis without vesselness filtering. However, vesselness filtering additionally functions as a means to convert the images from 12-bit greyscale ([0, 4095]) values to continuous ([0, 1]) values. Running the analysis without vesselness filtering, therefore, would require the definition of a new image threshold in the 12-bit domain, making comparison of the vesselness-filtered and vesselness-unfiltered results problematic.

Fig. 11
Effect of cluster volume thresholding. Curves show SMC density as a function of the relative radial coordinate (ρ − r ), where ρ is the radial coordinate as defined in Fig. 3a and r is half of the detected vessel diameter ( 1 2 d). Pooled data of five samples at a luminal pressure of 40 mmHg are shown [those (non-rotated) samples that were also used for assessing reproducibility (Table 1) and the effects of rotation ( Fig.  9)]. The lower volume threshold (V threshold,L , blue, dashed line) removed clusters mainly at the inner side of the wall, whereas the upper volume threshold (V threshold,U , red, dotted line) removed clusters mainly at the outer side of the wall. Percentages pertain to the number of clusters in each stack that were excluded/included, respectively (mean ± SD). SD standard deviation, SMC smooth muscle cell Our finding of an arterial helix could have implications for the constitutive modelling of the artery wall. SMCs oriented in a slight net helix will create a slight torsion in the artery wall, which changes with SMC contraction. Such effects can and should be investigated using a constitutive modelling approach. Additionally, the fact that SMC orientation shows dispersion may have mechanical implications. An excellent example of the mechanical consequence of the inclusion of structural dispersion in constitutive models is given by the collagen modelling performed by Gasser et al. (2006). They found that, due to the dispersion of collagen, the recruitment of collagen fibres occurs at much lower stretches. In addition, when they excluded dispersion, collagen fibres showed larger than physiological rotations before bearing any load. Clearly, such effects also play a role in the modelling of SMC contraction in the artery wall.
Our study was limited in the sense that we assumed the orientation of the SMC nucleus to be representative of the orientation of the entire SMC. This assumption was also made in Holzapfel et al. (2002), Peters et al. (1983), Walmsley (1983). Staining nuclei instead of entire cells has the advantage that it leads to more robust clustering results. Canham et al. (1982) studied the coalignment of the muscle cell and nucleus and found the upper limit on misalignment in three dimensions to be 2.4 • . We used an aspecific nuclear dye, staining nuclei of all cells in the artery wall. Therefore, non-SMCs (e.g. fibroblasts) may have been included in our analyses. However, by applying vesselness filtering, the narrow layer separation, and nucleus volume thresholds, we assume to have filtered out the majority of non-SMCs.

Conclusions
We conclude that (1) vascular SMC orientation can be quantified in 3D; (2) SMC orientation shows considerable dispersion, predominantly in the helical direction, which decreases transversally with increasing luminal pressure; and (3) 3D quantification of SMC orientation reveals a distinct right-handed helical component in both left and right carotid arteries. These quantitative distribution data are essential to improve constitutive modelling of the artery wall.
1. Consider the surface element on the unit hemisphere spanning from θ h to θ h + dθ h and θ t to θ t + dθ t . The area of this element is not dθ h dθ t , but dθ h dθ t cos θ t . Consequently, visualising orientations by using a rectangular θ h versus θ t plot over-emphasises orientations with large θ t . 2. A related error occurs when trying to represent an orientation of θ t = ±90 • . At what θ h should this orientation be drawn? Clearly, for θ t = ±90 • , θ h is undefined.
These problems are illustrated in Fig. 12a. A solution to this problem could be to display the hemisphere in 3D, as performed by Alastrué et al. (2010). Although correct, such 3D figures are, to our opinion, hard to interpret from a (2D) document.
The problem of displaying spherical information on a flat surface is not new and occurs, for example, when creating a 2D map of the 3D world. The transformation from latitudes and longitudes of locations on the surface of a sphere into locations of a plane is called a map projection (Snyder and Voxland 1989). A large number of map projections exist, all suited for a specific goal. Importantly, there is not one general "best" projection (Snyder 1987). For our application, we require the projection to be equal area, tackling the aforementioned point 1 (Snyder 1987) and illustrated in Fig. 12b. Many projections possess this property, but the one that is most often used to plot directional data is the Lambert azimuthal equal-area projection (Borradaile 2003). This is the projection used in this paper. Its geometrical construction is explained in Fig. 12c. Mathematically, it is obtained via the following equations: x = 2 1 + cos θ t cos θ h cos θ t sin θ h , and y = 2 1 + cos θ t cos θ h sin θ t . Fig. 12 Considerations for using map projections. a Typical problem when representing a three-dimensional orientation without using a projection. Consider four points A, B, C, and D with indicated spherical coordinates. Points A and B clearly represent a different orientation that is 90 • apart. Points C and D, however, appear to also represent a different orientation, whereas in reality, they are equal (they are both at the "top" of the hemisphere). b The problem in a can be solved using a map projection. By representing the orientations using a projection, points C and D coincide. The projection used is the Lambert azimuthal equalarea projection. c A Lambert azimuthal equal-area projection is created by using a plane tangential to the centre point on the unit hemisphere (point B in this case). All points will be projected on this plane. From all points on the hemisphere, a circular arc is drawn to the projected plane that (1) ends perpendicular to the projected plane and (2) has point B as its centre μ, overall mean; σ μ , between-animal standard deviation; σ r , reproducibility standard deviation. See "Appendix 2" for details on computing these values. Θ h , κ 1 , and κ 2 are defined in Fig. 6. d, diameter; I and O, inner and outer layer, respectively

Appendix 2: Reproducibility
To evaluate our findings, in five arteries, we acquired a second image stack after physically rotated the mounting chamber by 180 • around the z-axis (see "Discussion" section and Fig. 9). As these rotated stacks were acquired at the same location as the corresponding non-rotated stacks, these data can be used to analyse the repeatability of our method. In order to assess repeatability, we calculated three summary values of parameters P i j ∈ {Θ h,i j , κ 1,i j , κ 2,i j , d i j }, with d i j representing diameter. Subscript i ∈ {N, R} denotes non-rotated (N)/rotated (R) stacks; subscript j ∈ {1, 2, 3, 4, 5} denotes animal number.
The following summary values were computed: 1. The overall mean (μ): 3. The reproducibility standard deviation (σ r ): In these equations, M and S denote the circular mean and circular standard deviation operators for P i j = Θ h,ij ; and the arithmetic mean and standard deviation operators for P i j ∈ {κ 1,i j , κ 2,i j , d i j }, respectively. m j is used as a shorthand notation of M(P N j , P R j ).
Results are shown in Table 1. In all cases, the reproducibility standard deviations were smaller than the between-animal standard deviations. Reproducibility of the helix angle (σ r ) was good at 0.76 • (inner layer) and 1.37 • (outer layer), respectively.