## Abstract

In this work, we outline an automated method for the extraction and quantification of material parameters characterizing collagen fibre orientations from two-dimensional images. Morphological collagen data among different length scales were obtained by combining the established methods of Fourier power spectrum analysis, wedge filtering and progressive regions of interest splitting. Our proposed method yields data from which we can determine parameters for computational modelling of soft biological tissues using fibre-reinforced constitutive models and gauge the length scales most appropriate for obtaining a physically meaningful measure of fibre orientations, which is representative of the true tissue morphology of the two-dimensional image. Specifically, we focus on three parameters quantifying different aspects of the collagen morphology: first, using maximum-likelihood estimation, we extract location parameters that accurately determine the angle of the principal directions of the fibre reinforcement (i.e. the preferred fibre directions); second, using a dispersion model, we obtain dispersion parameters quantifying the collagen fibre dispersion about these principal directions; third, we calculate the weighted error entropy as a measure of changes in the entire fibre distributions at different length scales, as opposed to their average behaviour. With fully automated imaging techniques (such as multiphoton microscopy) becoming increasingly popular (which often yield large numbers of images to analyse), our method provides an ideal tool for quickly extracting mechanically relevant tissue parameters which have implications for computational modelling (e.g. on the mesh density) and can also be used for the inhomogeneous modelling of tissues.

## 1. Introduction

The mechanical behaviour of the arterial wall is mainly governed by the organization and composition of the three major microstructural components: collagen, elastin and smooth muscle cells [1–3]. The influence of these components on the cardiovascular function in health and disease has been the subject of extensive research [4–10]. While elastin is load-bearing at low and (to a smaller extent) high strains, it is collagen that endows the arterial wall with strength and load resistance, thus making it the most relevant mechanical tissue constituent [11–14]. Research indicates that changes in the mechanical properties of the healthy arterial wall play a role in arterial disease and degeneration [7]. For example, increased stiffening of the vessel wall with age may (among other factors) be related to an increased ratio of collagen to elastin as well as increased collagen cross-linking [4,15]. Higher wall stiffness may also be a contributing factor for atherosclerosis [6]. Furthermore, changes in the mechanical environment can lead to growth and rearrangement of collagen, which can cause enlargement of intracranial aneurysms, with an associated elevated rupture risk and mortality rate of 35–50% [8]. Collagen remodelling is also believed to play a role in the healing of dissecting aortic aneurysms [16]. Therefore, evaluating and monitoring morphological data on collagen within the arterial wall is essential to facilitate a better understanding of the underlying mechanical principles governing the behaviour of the vessel wall. Additionally, such data can be used to improve modelling of the cardiovascular system and increase our understanding of disease progression.

Many problems related to the mechanical function of arteries can be studied in the framework of finite element (FE) analysis. FE-based constitutive models for arterial tissues are available, even in commercial codes such as Abaqus (SIMULIA, Providence, RI, USA), and numerical modelling is a well-accepted means by which to gain insights into the functional relationships between structural and mechanical properties within arterial tissues, as well as studying tissue- and organ-level deformations or stresses [17]. Collagen fibres generally display wavy patterns in the unstressed arterial wall. When strain is applied, it leads to a progressive recruitment of collagen fibres, which align themselves in preferred (principal) directions, causing the characteristic nonlinear mechanical response of arterial tissues. Higher strains effect the orientations of collagen fibres [18], resulting in an improved fibre orientation coherence and smaller dispersions of the entire fibre families [1,19]. Anisotropic, fibre-reinforced material laws have been developed for modelling such tissues [2,3], and some of these models also include a measure of the fibre dispersion [20–23].

To visualize collagen, one can make use of either stained histological sections (common stains: picrosirius red (PSR), haematoxylin and eosin, Masson's trichrome, Movat's pentachrome), or different microscopy techniques, for example polarized microscopy [19,24]; electron microscopy (Xia & Elder [25] and references therein); fluorescence microscopy [26]; multiphoton microscopy (MPM) [27–30], featuring enhanced penetration depth in soft biological tissues, good optical sectioning and good resolution. Both, fluorescence microscopy and MPM use collagen as a source of second harmonic generation (SHG) [31–33] and autofluorescence [34–36], which allows for direct observation without staining. For a brief overview of structural quantification of collagen fibres in arterial walls, see Holzapfel [37].

To assess morphological collagen data from two-dimensional images in a simple, fast and automated fashion, we use a Fourier-based image analysis approach in combination with wedge filtering and progressive regions of interest (ROI) splitting. Fourier power spectrum analysis and wedge filtering are both well-established techniques for characterizing collagen organization, based on the autocorrelation of intensity fluctuations within a given image [25,38–43]. ROI splitting is a tool for extracting morphological data from sub-images, allowing us to assess changes of such data among varying length scales.

Our approach differs from aforementioned works, in that we focus on extracting relevant biomechanical material parameters that can be used directly in numerical modelling. We aim to provide an automated method to span the work of an experimentalist focused on obtaining high-quality images from soft tissues, and a computationalist focused on modelling of such tissues using material parameters extracted from two-dimensional images. Specifically, we demonstrate the following: first, the accurate determination of the principal directions of fibre reinforcement (i.e. the preferred fibre directions); second, the quantification of the fibre dispersion about this principal direction; third, the gauging of appropriate length scales for characterizing important features of the tissue morphology (which has implications on the FE mesh density of the geometrical model); finally, a comparison of the entire fibre distributions (instead of averaged measures) among different length scales using higher-order statistics.

## 2. Methods

### 2.1. Images

To validate and demonstrate our methods in a controlled setting, we generated two binary test images, as shown in figure 1. Both test images are made of 1511×1511 pixels and contain 100 straight white lines, randomly distributed so that no two lines touch each other. The lengths of the lines vary according to a Gaussian distribution with a standard deviation of 35 per cent of the mean length. In figure 1*a*, all lines are aligned parallel at 25° (0° and 90° correspond to the horizontal and vertical axes of the image, respectively). In figure 1*b*, 30 per cent of the lines have been randomly rotated to yield a Gaussian distribution with a mean of 45° and a standard deviation of 10°. All random numbers were generated using Matlab (MathWorks Inc., MA, USA).

Figure 2 gives an overview of the three microscopy images that we chose for our studies. Figure 2*a* shows a greyscale version of a representative micrograph of a PSR-stained, in-plane tissue section of the human intima (image size: 1936 × 1936 pixels). We obtained the image via polarized microscopy through the course of a different research study; details regarding the tissue preparation and histological procedures can be found in Schriefl *et al.* [44]. The image features two almost symmetrically arranged collagen fibre families with respect to the main cylinder axes, indicated by the white dashed lines in figure 2*a*). The circumferential and axial directions of the vessel correspond to 0° and 90°, respectively. Figure 2*b* shows an MPM image of one collagen fibre family taken from an in-plane section of the human media (image size: 1237×1237 pixels). The SHG signal of collagen was detected using a commercially available coherent anti-Stokes Raman scattering and SHG microscope system based on a Leica SP5 (Leica Microsystems, Inc.). The device is equipped with a picosecond solid-state laser-based light source (picoEmerald; APE, Germany; HighQ Laser, Austria) with an integrated optical parametric oscillator (OPO). The SHG signal of collagen was generated by tuning the OPO to 830 nm. Backscattered SHG signal was collected in epi-mode using a non-descanned detector and a BP 465/170 filter. Imaging was performed using a Leica HCX PL APO CS 40.0×0.25 oil objective. Figure 2*c* shows a SHG image of collagen fibres from a transversal section of chicken cartilage (image size: 1443 × 1443 pixels). For details regarding tissue preparation and imaging, see Lilledahl *et al.* [43].

### 2.2. Fourier transformation and wedge filtering

To obtain material parameters for applications in computational modelling, our first task is to extract quantitative data regarding collagen fibre orientations from two-dimensional images. Towards this end, we represent the original greyscale image by a distribution function *f*(*x,y*), where (*x,y*) defines a point in the real two-dimensional image. Using a two-dimensional fast Fourier transformation (FFT), we obtain where (*u,v*) is a point in the Fourier space. To avoid frequency-domain effects originating from periodic discontinuities at the boundaries of the original image, we first apply a window using a raised cosine function that reduces the greyscale values to 0 at the image periphery (figure 3*a*). After FFT, we perform a coordinate shift, transforming the lowest spatial frequency to the origin, from where it increases as we move towards the image edges. The FT is then multiplied with its complex conjugate *F**, yielding the power spectrum *P*, namely
2.1

Figure 3*b* displays the Fourier power spectrum shown with a logarithmic intensity scale of the windowed image of figure 3*a*. Note the 90° shift of the axes owing to the Fourier transformation. Frequency components of the collagen fibres along different orientations are represented as changes in amplitude, *I*(*Φ*), along a specific angle *Φ*. The two fibre families are visible as two white, elongated clouds oriented between the major axes.

The collagen fibres are now discriminated by spatial frequency and orientation. We use wedge-shaped orientation filters to extract the fibre angles and their corresponding amplitudes from the power spectrum [38]. The amplitudes *I*(*Φ*) are obtained by summation of every *P*(*u,v*) within individual wedges (1° wedge width), yielding a discrete distribution of relative amplitudes as a function of the angle *Φ*. For these angles to correlate with the real image, we shift them back by 90°, as shown in figure 3*c*. For subsequent fitting, we smooth the relative amplitudes using a moving average filter with a span range of 5°. The effect of data smoothing is illustrated in figure 3*d*.

### 2.3. Distribution fitting

To describe the angular distribution of fibres, say , we use a *π*-periodic von Mises distribution similar to [3,21], given by
2.2where *Φ* denotes the angle, the concentration parameter determining the shape of the distribution, and the location parameter describing the mean (or principal) fibre orientation (figure 4*a*). Both parameters (*b* and *μ*) are determined during the fitting process. *I*_{0}(*b*) in (2.2) is the modified Bessel's function of the first kind of order zero and the function of order *n* is defined by
2.3

Note that the von Mises distribution is an angular distribution and a close approximation to the wrapped normal distribution (which is the normal distribution wrapped around the circle) [45]. The von Mises distribution is normalized according to Holzapfel & Ogden [3], such that 2.4which yields the density function 2.5

Note that the location parameter *μ* becomes zero in the model because the *x*_{1}-axis of the coordinate system and the mean fibre direction are made to coincide (figure 4*b*). From equation (2.5), we can compute the dispersion parameter
2.6with . For an isotropic fibre distribution (*b* = 0), in the case of complete fibre alignment, i.e. no dispersion, , as illustrated in figure 4*a*.

We obtained an analytical solution for *κ* using the trigonometric identity , yielding the result shown in (2.6). Because is *π*-periodic, (2.3) also holds for the integration limits . Note that *κ* is half of the circular variance [46]. The dispersion measure *κ* can be used to construct a symmetric structure tensor **H** that describes the fibre distribution in a continuum mechanical framework. Note that this approach assumes a planar, symmetric fibre distribution and is appropriate only for modelling thin lamellar structures or a subset of three-dimensional problems under the assumption of an in-plane arrangement of the fibres (see §4 for more details).

To fit the distribution from figure 3*d*, we use maximum-likelihood estimation (MLE). The fundamental properties of the FT (see §4) allow us to generate the fibre angle dataset (required for MLE) from the amplitude distribution. For this purpose, we treat the distribution of relative amplitudes as a histogram, where the number of fibre angles equates to the value of the corresponding amplitude. For instance, from a relative amplitude value of 20 per cent at 0° we generate 20 angles with 0°. This approach ensures that the created angular dataset can reproduce the original amplitude distribution and provides us with a sufficient number of angles for the MLE.

As can be seen in figure 3*d*, the distribution is composed of two fibre families, with a summation of the amplitudes at the overlapping region. Therefore, we use a mixture of two von Mises distributions, given by
2.7

Thus, four parameters are fitted, namely the concentration parameters *b*_{1} and *b*_{2} and the location parameters *μ*_{1} and *μ*_{2} of the two distinct distributions. Note that the mixture of the two von Mises distributions in (2.7) does not need to be normalized to *π* by a constant factor, because it is only used in the fitting process for determining the four fitting parameters [47], and is not used further in the model.

### 2.4. Varying length-scale analysis

To determine meaningful measures of fibre orientations at different length scales, we divide the original image continuously into ever smaller sub-images, denoted as ROIs. With each subsequent dividing step *n*, the number of ROIs grows exponentially by *n*^{2}, as illustrated for *n* = 1 − 4 in figure 5*a*. Every ROI is windowed (see earlier text) and the Fourier power spectrum is calculated according to (2.1). We then fit a line through the centre of the power spectrum in a least-squares sense [25], from which the overall fibre angle *φ* is determined for each ROI (figure 5*b*). Note that *φ* could also be determined by means of wedge filtering, but because we are interested only in one parameter for each ROI (instead of the entire angular distribution of amplitudes for each ROI), it is computationally (much) more efficient to fit a line through the power spectrum. This yields *n*^{2} angles *φ* within the original image for each dividing step *n*, illustrated in figure 5*c* as a histogram for *n* = 30 (0° and 90° correspond to the horizontal and vertical direction of the image, respectively). From this angular distribution, we determine two measures describing the average behaviour of the calculated fibre orientations: the mean angle and the median angle . We provide both measures in the illustrative example (*n* = 30) in figure 5*c*. Normalizing the histogram to yields the computed probability mass function that we use in the subsequent analysis. Using both and , we aim to determine an appropriate range for *n*, where the resulting distribution of *n*^{2} fibre angles reflects the inhomogeneous morphology of the collagen fibres in the image. We achieved this goal by determining a range where and display (relatively) stable behaviour among subsequent partitions (*n*'s).

To test whether the calculated distributions from dividing images into ROIs reflect the inhomogeneity of the morphology, we plot statistical measures across different subdivisions (i.e. different *n*'s) of images. In this context, the weighted error entropy is the most relevant statistical measure [48,49]. Error entropy has been used in the context of stochastic learning as a reliable metric, and has the ability to capture arbitrary statistics [50]. The error entropy *E* is defined as
2.8where
2.9

While *p*_{d} quantifies the difference between the computed (*p*_{c}) and actual (*p*_{a}) PMF, the error entropy *E*(*n*) quantifies the uncertainty in the difference between the two distributions.

Hence, when the two distributions are identical, the error entropy is zero. However, this assumes that both distributions *p*_{c} and *p*_{a} are centred around the same value (a difference distribution that is Dirac delta at a different mean value has less error entropy than a difference distribution, which is centred around the correct mean). To prevent this occurrence, we weight the error entropy with difference in mean values. The weighted error entropy (*E*_{w}) is given by
2.10where *k* is a weighing constant, and and are the mean values of the given distributions. In order to determine the *E*_{w}, the actual PMF (*p*_{a} in equation (2.9)) must be computed if it is not known. For the first two test images (figure 1*a,b*), *p*_{a} is the input PMF (the normalized angular distribution of the lines). Because *p*_{a} is usually unknown for experimental images (figure 2*b,c*), we choose a specific PMF to be *p*_{a}, where the gradient of absolute entropy of the PMF stabilizes and reaches a plateau. To ensure that the chosen *p*_{a} is robust, we perform a perturbation analysis. If *n* yields the ROI size when the gradient of absolute entropy stabilizes, we also choose theoretical distributions corresponding to (*n* − 3) and (*n* + 3). A robust choice would imply that our results (the plots of *E*_{w}) do not change for such small perturbations in the chosen theoretical PDF *p*_{a}.

## 3. Results

We have developed automated methods to extract collagen fibre orientations and associated dispersions at different length scales from two-dimensional images that can then be used in computational modelling. To this end, we focused on three key parameters: first, the dispersion parameter, *κ*, which is a well established and used measure of anisotropy [3,21,23]; second, the fibre angles and , as a measure for the average orientation of a given fibre distribution; third, the *E*_{w}, which analyses and compares entire fibre distributions rather than the averaged fibre orientations, which is, therefore, useful when average measures are not sufficient. Furthermore, the behaviour of the average fibre angles and *E*_{w} at different length scales allows us to determine an appropriate range of ROI sizes, which will yield fibre distributions representing the tissue morphology from the image, which has implications for computational modelling (e.g. on the mesh density). Our method is based on transforming an image to Fourier space where the power spectrum is calculated, from which the earlier-mentioned parameters can be extracted (figures 3 and 5).

### 3.1. Distribution fitting

To obtain the dispersion parameter *κ*, we fit a given angular distribution using MLE. As an illustrative example, we choose the (more challenging) case of an image featuring two (rather than one) collagen fibre families (figure 2*a*). The result of fitting a mixture of two von Mises distributions (equation (2.7)) is shown in figure 6. The fitting parameters are *b*_{1} = 2.503, *μ*_{1} = −39.6° and *b*_{2} = 2.149, *μ*_{2} = 39.4° for the two distributions, respectively (*b* denotes the concentration parameter and *μ* the location parameter). Using (2.6) and the shape parameters *b*_{1} and *b*_{2}, the dispersion parameters and were calculated. To quantify the goodness of the fit, we use the Pearson's *χ*^{2} test (yielding a *p*-value) and the coefficient of determination *R*^{2}. For the fit in figure 6, we obtained *p* = 0.8535 and *R*^{2} = 0.935.

### 3.2. Varying length-scale analysis

Figure 7*a* shows the behaviour of the average fibre angles and with increasing (decreasing ROI size) for the first test image from figure 1*a*. Initially, both measures yield the correct angle of 25° in accordance with the actual input PMF , which is a Dirac delta distribution located at 25° (figure 1*a*). The mean angle starts to deviate from the correct angle around *n* = 13 (left arrow in figure 7*a*) marking the upper cut-off point for , whereas the median angle shows a more stable behaviour yielding the appropriate angle until the upper cut-off point for around *n* = 36 (right arrow).

For the *E*_{w} in figure 7*b*, *p*_{a} was chosen to be a Dirac delta distribution located at 25°, corresponding to the actual PMF of the first test image. Note that if, instead of the Dirac delta distribution, we choose *p*_{a}, based on the results of the perturbation analysis (not plotted), where the gradient of absolute entropy stabilizes (rather than the known Dirac delta distribution), it results in the theoretical distribution, verifying our scheme for computing the unknown actual PMFs for images of soft tissues. The *E*_{w} is quite stable until around *n* = 13. While the reliable estimate for the median holds until *n* = 36 (upper cut-off point for in figure 7*a*), the entire PMF shows only a relatively stable behaviour for a third of that range (until *n* = 13), similar to .

The results in figure 8*a* for and of the second test image, which includes some dispersion (figure 1*b*), show a very different behaviour with respect to the first test image. Initially, the mean angle yields increased values around 30° (owing to the fact that 30% of the fibre angles are distributed with a mean of 45°) and starts to drift off at *n* = 17. Overall, the values for fluctuate with increasing *n*, and a stable behaviour is never observed. The median angle on the other hand fluctuates only in the beginning of the analysis (*n* = 1–10), and stabilizes at a lower cut-off point of about *n* = 11 (see arrow); is then stable up to *n* = 24 () before it starts to drift away. The *E*_{w} in figure 8*b* indicates that after an initial increase, the PMFs do not change much until *n* = 14, after which shows a steeper increase.

Figure 9 displays the results of the average fibre angles in (*a*), and the *E*_{w} in (*b*) of a SHG image featuring collagen fibres in an in-plane section of the human media (shown in figure 2*b*). Both and have a lower cut-off point at *n* = 3. After *n* = 19, drifts away sharply, whereas displays a (relatively) stable range until the upper cut-off point at approximately *n* = 40. The weighted error entropy curve features two regimes: (i) a very stable region until *n* = 18, and (ii) a monotonic, almost linear increase henceforth. This is similar to the drift in the curve for , whereas the drift is less steep for .

The results from the analysis of the SHG image featuring collagen fibres of the chicken cartilage (shown in figure 2*c*) are displayed in figure 10. Both, and in panel (*a*) have a lower cut-off point at *n* = 4, followed by a stable domain until the upper cut-off point at *n* = 25. In the stable domain, yields values approximately 6° higher than . The *E*_{w} in figure 10*b* has a trend similar to that in figure 9*b*, with a stable regime followed by an almost linearly increasing curve. At the beginning of figure 10*b*, up to around *n* = 6, *E*_{w} reveals changes in the PMF owing to aforementioned artefacts. This is not the case for *E*_{w} in figure 9*b*, because the analysed image (shown in figure 2*b*) does not contain regions with very different fibre orientations (compare figure 2*b* with 2*c*).

## 4. Discussion

In this study, we outline an automated Fourier-based imaging analysis method to extract and quantify material parameters characterizing collagen fibre orientations from two-dimensional images. Our method combines Fourier power spectrum analysis and wedge filtering (two well-established image analysis tools for morphological data extraction; [25,38,40]). Additionally, we show how ROIs (sub-images) can be used to obtain useful orientation information among different length scales.

We aim to extract data from which we can (i) determine material parameters for computational modelling of soft biological tissues using fibre-reinforced constitutive models and (ii) gauge which length scales are most appropriate to capture tissue morphology, which can also have implications for computational modelling (e.g. mesh density). Specifically, we focus on three key parameters: first, we compute the dispersion parameter *κ*, a measure of anisotropy [3,21,23]. Mainly, *κ* is based on a three-dimensional (three-dimensional), rotationally symmetric von Mises distribution. Because we deal with two-dimensional images, we make use of the two-dimensional equivalent as discussed in Holzapfel & Ogden [3]. Second, we introduce two measures for the average orientation of a given distribution , and show that is the more appropriate measure to use if the collagen fibres are dispersed. In the presence of fibre families, the average orientation corresponds to the preferred (or principal) fibre direction of one fibre family [44], and hence can be used as a parameter for continuum mechanics-based constitutive models (see earlier studies [2,51–54] and references therein). Third, we calculate the *E*_{w} as a measure of changes in the entire fibre distributions at different length scales, as opposed to their average behaviour.

### 4.1. Analysis validation

In our approach, three fundamental properties of the FT are pertinent, ensuring that the distribution of relative amplitudes in the Fourier space is a representation of the collagen fibre distribution. They are (i) the rotation of a collagen fibre results in an equal rotation in the Fourier space (*the rotational property*), (ii) the Fourier transform of a region containing several fibres is equal to the sum of the individual Fourier transforms of the same fibres (*the addition theorem*), and (iii) a spatial shift of a fibre does not affect the amplitude of its Fourier transform (*the shift theorem*). Hence, Fourier space provides a proper representation of the fibre orientations of a two-dimensional image.

Generally, the dimensions of any original image used for our analysis should be much larger than the structures of interest shown in the image (in our case, collagen fibres) to ensure that they can still be contained within ROIs at increasingly smaller length scales. To validate the functionality of our method, we created two test images shown in figure 1*a,b*. The purpose of the first test image (figure 1*a*) is to verify the capability of the method to extract correct angular data among various length scales. The distribution of all lines in the first test image corresponds to a Dirac delta function at 25°. Both average measures , yield the correct angles of 25° in the range of *n* = 1–12 (figure 7*a*). As expected, is more stable at increasingly smaller ROI sizes (increasing *n*) than .

To investigate which of the two average measures is more appropriate if the fibres are dispersed, we modified the first test image to include some dispersion, resulting in the second test image shown in figure 1*b*. From the results in figure 8*a*, it becomes evident that the mean average measure does not stabilize. It is, therefore, not suited for identifying a range of ROI sizes that would yield angular distributions that are representative for the orientations of the lines (simulated fibres). On the other hand, the behaviour of the median average measure shows three interesting characteristics: (i) it yields a lower cut-off point, which represents a maximum ROI size above which the average angles of the entire distributions fluctuate highly as *n* changes; (ii) it yields an upper cut-off point which marks a minimum ROI size below which the angular values start to drift away; and (iii) it identifies a stable range of ROI sizes between both cut-off points. Such stable ranges can be observed only if the underlying distributions (from which is determined) reflect the orientations of the features in the image. If, for example, a ROI size is too large, details within individual ROIs are smeared out and information is lost. On the other hand, if a ROI size is too small, the information within individual ROIs is compromised by, for example, image artefacts, resolution limits or edge effects.

### 4.2. Distribution fitting

We use MLE to fit a mixture of two von Mises distributions to the angular distribution obtained from the image in figure 2*a*, featuring two collagen fibre families in the human intima of the thoracic aorta [44]. The results of this fitting are shown in figure 6, and as a measure of the goodness-of-fit we determined the *p* and *R*^{2} values. Because the Pearson's *χ*^{2}-test depends on the size of the dataset (i.e. the number of generated angles), we also calculate the coefficient of determination *R*^{2}. Because *R*^{2} is computed by comparing the original distribution with the estimated PDF, it does not depend on the number of angles and is, therefore, an additional measure for the goodness-of-fit independent of *p*. We also emphasize that *R*^{2} is not optimized during the fitting procedure.

We would like to note that although it is quite common to fit a curve to a histogram using the least-squares method, there are some drawbacks to this approach that should be kept in mind. For one, least-squares fitting depends on the bin size of the histogram. Another potential pitfall arises from the *normality assumption*, which states that the errors are normally distributed with mean zero. Because bin counts in a histogram are non-negative this assumption does not hold. Also, the *constant variance assumption* and the *independent-errors assumption* are not justified when fitting distributions [55].

### 4.3. Fibre angles and dispersion

In the literature, the importance of collagen fibre orientations on the mechanical behaviour of arterial walls has been well established [3,37,44,56]. To account for the effect of strain on changes in the fibre orientations and the angular distributions, to approximate the *in vivo* strain state of arteries and to ensure more straightened fibres necessary for angular measurements, the investigated samples were usually either pressurized or pre-stretched biaxially beforehand [8,19,57,58]. Our approach allows us for a fast and automated determination of the angular distribution from two-dimensional images. By fitting a given distribution (in our example a distribution of two fibre families, see figures 2*a* and 3) with a mixture of two von Mises distributions, we can compute the principal fibre orientations , and the corresponding dispersion parameters , , one for each fibre family (figure 6).

Note that this approach assumes a planar, symmetric fibre dispersion for each fibre family individually. Such two-dimensional data can be directly applied to mechanical modelling of thin lamellar tissue structures on the basis of membrane or thin shell theory, for example, cerebral arteries, using appropriate constitutive models (see earlier studies [3,21] for further discussion). However, many arteries behave as thick-walled cylindrical tubes (e.g. the human aorta), and hence three-dimensional modelling approaches are required. In this case, for a subset of fully three-dimensional problems, one can make use of the (mainly) in-plane arrangement of the collagen fibres in arterial walls. For example, in the case of the human descending aorta and common iliac arteries, the collagen fibres lie in the plane of the tissue (angular deviation ) [44]. The assumption of in-plane arranged collagen fibres is not necessarily appropriate for other soft tissues. For example, in the articular cartilage, the fibres display very different orientations depending on the specific zone. In the tangential zone, they lie within the transversal plane, a variety of orientations are observed for the middle zone, while in the deep zone the fibre arrangement is perpendicular to the transversal plane [59]. To ensure that our analysis yields correct material parameters, most of the collagen fibres should lie within the sectioning (imaging) plane in order to avoid projections to this plane, which would lead to deviations in length and (possibly) angles from the true morphological state. Therefore, in the case of an in-plane arrangement of the collagen fibres, the two-dimensional collagen fibre data obtained through our method may also be used for three-dimensional modelling (using the two-dimensional *κ* parameter according to equation (2.6)) under the assumption of strictly in-plane organized collagen fibres (neglecting the out-of-plane deviations is a simplification of the problem with regards to computational modelling, motivated by experimental observations [44]).

### 4.4. Varying length-scale analysis

We developed an analysis method that continuously splits an image into increasingly smaller ROIs to obtain a physically meaningful measure of fibre orientations across length scales. In most applications, mean and median measures suffice when we are interested in the average behaviour of fibre families. Using two test images, we have shown that the median angle is the more appropriate measure to use for this purpose.

To demonstrate the applicability and usefulness of our method on microscopy images of soft tissues containing collagen fibres, we ran our analysis on SHG images of the human media and the chicken cartilage (figure 2*b* and 2*c*, respectively). Results for both images, as shown in figures 9*a* and 10*a*, reveal a stable range of ROI sizes, located between the lower and upper cut-off point, where shows only small fluctuations (). The advantage of the length-scale analysis can be seen in figure 10*a*, where a ROI size in the range of leads to highly fluctuating values for , followed by a wide range (*n* = 5−25) of stable behaviour before slowly drifts away. Using our method, a finite element analyst can, for example, determine an appropriate mesh density (based on a ROI size where shows stable behaviour) that also captures the relevant tissue structure. For example, the cartilage tissue in figure 2*c* features very different fibre orientations at the perimeter compared with the centre of the tissue. If a ROI size is too large, this specific tissue inhomogeneity will be lost. Only if the ROI size is small enough, can such spatial changes in orientation be resolved. Analysis of the cartilage tissue shown in figure 10*a* tells us that a mesh density corresponding to a ROI size between is appropriate to capture important structural features of the image.

However, if average measures are not sufficient (e.g. quantification of the risk of a pathological condition such as an abdominal aortic aneurysm [60]), and to depict how higher-order statistical measures vary across different values for *n*, we compare weighted error entropies across different values *n*. There are different measures that can be used to compute the distance between two distributions [61–63]. For example, the Kullback–Leibler divergence [64] is a measure of distance between two distributions, but it is not symmetric and is only defined for two distributions whose range of non-zero values is the same (which is not the case in our examples). Mutual information, as another example, quantifies how much the knowledge of one variable reveals about another variable. However, correlation information between the two variables is needed, and such information is not available in our case. The *E*_{w} is an appropriate measure for our purpose. It compares entire fibre distributions (rather than the averaged behaviour) and reveals how much of the tail information is lost when averaged behaviour is considered. For the microscopy images in our study, the *E*_{w} shows that the entire fibre distributions start to change with increasing *n* before a change in can be observed (cf. (*a*) and (*b*) in both figure 9 and figure 10). Therefore, if one is interested in the behaviour of the entire fibre distributions, then a study of *E*_{w} helps to identify ranges of *n* where calculated fibre distributions are similar, and at that length scale, they start to differ.

The usefulness of both parameters and *E*_{w} is illustrated in the results for the second test image (figure 1*b*) in figure 8. For example, if one is interested in determining an appropriate ROI size capturing local changes in tissue morphology, then a corresponding ROI size ranging from *n* = 11 to *n* = 20 is suitable, because in this range shows a (relatively) stable behaviour. On the other hand, if the behaviour of the entire PMF is of interest, then the result for *E*_{w} tells us that we are losing information starting around *n* = 14, which is insight that cannot be gained from average measures.

Another valuable aspect of our methods pertains inhomogeneous computational modelling of tissues. For this purpose, the location-specific structural parameters (fibre angle *κ*) obtained for every ROI can be incorporated element-wise into an FE mesh. Such patient-specific data allow for a more precise modelling of the cardiovascular system and could help improve our understanding of the interaction between collagen fibre morphology and arterial wall mechanics.

### 4.5. Limitations

The method we present extracts morphological data from two-dimensional images and is, therefore, not able to capture true three-dimensional fibre orientations and dispersions. As it is the case with other methods that aim to quantify fibre angles, it is a prerequisite for the collagen fibres to be (mostly) straightened (see Schriefl *et al.* [44] and references therein). Such data can be used for modelling within the framework of membrane or thin shell theory, but are limited to a subset of three-dimensional physical problems and soft biological tissues (see §4). To fit the fibre distribution, we use MLE that requires a sufficient number of fibre angles. Our method does not measure fibre angles directly, but extracts them from Fourier power spectra. While this automated approach for extracting orientational information from images is well established, we want to emphasize that the resulting fibre angle distributions do not represent accurate fibre angle counts. The choice for the upper and lower cut-off points to determine a (relatively) stable dividing range (*n* range), with the aim of obtaining physically meaningful fibre orientations which reflect the tissue morphology, must be made manually. Depending on the input image, this choice can either be very evident or in the case of slow drifts with increasing *n*-values less defined.

## Acknowledgments

The use of autopsy material from human subjects for this study was approved by the Ethics Committee of Medical University Graz (no. 21-288 ex 09/10).

This work was partly supported by the European Commission under the seventh Framework Programme in the scope of the project SCATh: Smart Catheterization (grant agreement no. 248782), and by an AHA postdoctoral fellowship grant to S.S., University of California, USA. We also thank Dr Heimo Wolinski, University of Graz, Austria, and Dr Magnus B. Lilledahl, Norwegian University of Science and Technology, Norway, for generously providing us with MPM images (figure 2*b*,*c*).

- Received April 26, 2012.
- Accepted June 12, 2012.

- This journal is © 2012 The Royal Society