Applying limiting entropy to quantify the alignment of collagen fibers by polarized light imaging

1.   Introduction

Collagen is a kind of vital and abundant structural protein in the human body. So far, more than 20 types of collagen have been identified, which are characterized by a triple helical domain ^[1]. However, type I, II, III, and IV collagen consist over 90% of the collagen in the body. In particular, type-I collagen is the most abundant collagen type in many tissues, e.g., skin, cornea, or bone, which consists of nanometer-scale fibrils. Collagen fibers are bundles of fibrils with a diameter of 2µm to 5µm ^[2,3,4]. These fibers are the prominent components of extracellular matrix (ECM), which normally provides the underlying structural framework of biological tissues ^[5]. Alterations of collagen fibers organization parameters, specifically fiber orientation and alignment, play an integral role in many diseases including cancer ^[6,10,11], scaring ^[7], scleroderma ^[9], and glaucoma ^[8]. Thus, quantifying these parameters may be important to study the properties of healthy and diseased tissues and to develop a potential innovative diagnostic tool.

Lichen Sclerosus (LS) is a chronic disorder of skin characterized by hypopigmented fibrotic plaques. It is most commonly found in women, but also found in men and children. Though any skin site may be affected, LS has a preference for the anogenital area ^[12,13]. For example, Vulvar Lichen Sclerosus (VLS) is mainly found in the anogenital area of pre-pubertal and post-menopausal females ^{[14,15,16,17]}. It is among one of the most common referrals for vulvar pruritus and is the most common form of structural change to the vulvar region ^[17,18,19]. According to an educational bulletin by American College of Obstetricians & Gynecologists (ACOG) ^[20], VLS accounts for approximately 10% of patients seen in vulvar clinic visits. Developing lesions lead to vulvar pain and sexual dysfunction, with a huge impact on the patients' quality of life ^[21,22,23]. Moreover, affected women have an 4-5% risk of developing squamous cell carcinoma of the vulva ^[22,24]. However, general practitioners may overlook this skin disorder despite the characteristic clinical and histological appearance, which can lead to considerable delays in diagnosis ^[25,26]. Currently, histopathological examination is the gold standard for differentiation from other inflammatory or tumorous vulvar diseases ^[27]. Although the etiology of VLS is still unknown, the histopathology of VLS is characterized by alterations in the epidermis and dermis as well as inflammatory cell infiltration ^[27]. Among these histopathologic features, homogenization of collagen (hyalinization) has been recognized as the most prominent pathological change ^[27,28,29]. Yet, there were few quantitative studies about this histological diagnostic criterion to date. In most literatures, it was defined as "a zone of thickened, homogenous collagen fibrils" ^[27,28,29]. It is subjective and may lead to misdiagnosis in microscopic examinations of the lesions. However, in a recent study by Hiram et al., they performed scanning electron microscopy of the dermis from a lesion of LS and found that collagen fibers were adherent to each other as well as packed together, and round structures were also observed in the homogenization zone ^[30], as shown in Figure 1. Since collagen fiber alignment is involved with the pathological processes of LS (including VLS), accurate characterization of the structure information of collagen fibers in VLS tissues may be of great interest to understand the pathology and possibly to improve the histological diagnostic criteria.

A variety of imaging methods are available to retrieve structural information of collagen fibers in tissue ^{[33,35,36,38]}. For examples, scanning and transmission electron microscopy has an advantage of visualizing individual collagen fibers at nanometer scale. However, it requires complicated sample preparation procedures, and is expensive ^[37]. Second harmonic generation (SHG) imaging provides label-free detection of collagen fibers with unequalled specificity and contrast but is relatively expensive as it requires multiphoton microscopy ^[34,39,40]. Recently, due to the system simplicity, low cost and effectiveness, polarized imaging methods have been widely used to characterize collagen fiber organization and to detect the pathological changes of biological tissues ^{[31,41,42,43]}, especially those related with fibrous structural changes ^[6,7,43]. It has been found that linear polarization images with single incident polarization can be very different when the incident polarization or sample orientation distribution is changed ^[44,45,46]. Yang et al presented a polarized light spatial frequency domain imaging method to quantify collagen fiber orientation on fresh bovine tendon and porcine aortic valve leaflet tissue (PVL). They found that collagen fibers, when placed between two synchronously rotating polarizers, showed maximum brightness if the polarization direction of polarizers corresponds to preferred local fiber orientation ^[44]. Also, Bancelin et al determined collagen fiber orientation in histological slides of uterine cervical tissue using Mueller polarimetry (MP) and validated the results by SHG. In their work, the orientation of collagen fiber is expected to be collinear with the azimuth of the slow axis of its retardation ^[32]. Therefore, the local fiber orientation can be determined by rotating the two polarizers synchronously and finding the rotation angle in individual pixels in an image stack, at which collagen fiber shows maximum/minimum brightness.

To detect the anomaly of fiber alignment, a measure of disorder for local fiber orientation map is desired. Such measure should satisfy two requirements to have a high potential for clinical application, including (i) invariant to the angular translation and (ii) insensitive to the resolution and accuracy of the system or equipment. Both requirements ensure that the characterization of malignant tissue is generic and the results from labs can be easily put into practice. The traditional methods to quantify disorder involves use of spread, kurtosis, Chi-square test and Dykstra-Parsons coefficient ^[33,47,48]. However, all those indexes are not well designed for assessment of distributions. For instance, spread (standard deviation) ^[33], one of the most commonly adopted statistic indices, is heavily dependent on the relative angles with respect to the reference axis. Illustrated by Figure 2, the change of reference axis can alter the spread, as the spread is not invariant to periodic translation. Nonetheless, the reference axis is linked with the specimen's placement under microscope. The variation of positioning, which is impossible to eliminate, probably contributes to large variation of the quantification result. A more common and advanced way for quantifying disorder is to employ cyclic statistics ^[49]. It connects angle value with complex number and the statistical indexes are evaluated by manipulating complex numbers. Though effective cyclic statistics is, it is sensitive to the accuracy of rotor motors in our system, which inspires proposition of another viable index.

Given the success in communication theory, information entropy, which deals with the probability distribution and is irrelevant to the coding of variable, is supposed to be a solution to the placement problem. The discrete information entropy was proposed by Claude Shannon and is widely adopted in communication theory to quantify the information conveyed in the system ^[50]. There are various interpretations of discrete information entropy in communication field, including a measure of uncertainty, information of a system and even minimal expected number of queries to determine a code ^[51]. Entropy is also an important concept in statistical mechanics, where a system is comprised of enormous number of elements. Since each element occupies a certain microstate, the system is characterized by the probability of each microstate occupied. Therefore, the entropy can be adopted to quantify the level of disorder in such a system ^[52]. In our study, the local fiber orientation map by signal processing over a stack of polarized images, can be deemed as a physical system in statistical mechanics. The orientation of each pixel represents a microstate and the map is thus characterized by the distribution of fiber orientations. Therefore, discrete information entropy should be an effective tool to quantify the disorder of fiber orientations in the map. Some significant works have been done. For example, Ducourthial et al. proposed a revised discrete entropy to measure the orientation distribution of collagen fibers generated by SHG system ^[53]. They normalized the discrete information entropy to reduce the influence of bin size. Chiverton et al. applied entropy to analyze the steel fibers in X-ray Microtomography (XCT) data ^[54]. Their method successfully quantified 3D XCT images in multiscale. Though powerful the abovementioned methods are, the potential in clinical application is still limited. Both methods are based on discrete entropy and relies on discretization of the system. To tackle this problem, an index based on limiting entropy is proposed and described in detail in this manuscript.

In this paper, we established a rotating linear polarization imaging system combined with a novel parameter, homogeneity index (HI), which is based on limiting entropy, to quantify the alignment quality of collagen fibers. Compared with other parameters, i.e., such as cyclic variance, discrete entropy and normalized entropy, HI is more robust to the sample orientation at different placement and systematic sampling increment (angular step size). Since there is a trade-off between total number of bins and accuracy of probability estimation in each bin, being insensitive to number of bins indicate lower requirement for system resolution while detaining high estimation accuracy. Verification was performed on a set of Poly Lactic Acid (PLA) filament phantoms with controlled alignment quality of fibers to validate the ability of our proposed method to measure the level of disorder of the local fiber orientation. Then, we applied this method to assess the collagen fiber alignment in both normal and VLS biopsies for further understanding of the pathology of VLS and possible improvement of the histological diagnostic criteria for VLS. The structure of this paper is as follows: Section 1 is the introduction. Section 2 is a detailed account of theoretical foundation, the experimental system, data acquisition and analysis procedures. The experiment process and results are presented in Section 3, followed by the discussion in Section 4; Finally, Section 5 is the conclusion.

2.   Materials and methods

2.1.   Theory background

This section first presents a novel quantification method, homogeneity index (HI), to characterize images of spatially distributed orientations. Then the proof that HI is robust to variation of the histogram binning size and the precision of measurements (angular step size), is given. Moreover, the traditional cyclic statistics and other entropy-based indexes are introduced briefly and the comparisons between HI and those indices are conducted. Finally, the computational verification is presented.

2.1.1.   Derivation of Homogeneity Index

The discrete entropy was proposed by Claude Shannon and widely used in communication theory to quantify the information conveyed in the system ^[50]. The name originated from the entropy in thermal mechanics and statistical mechanics. Discrete information entropy is independent of the coding of variables and can reflect the complexity and disorder of the investigated system. In this work, it can be written as:

where $P (\theta_i)$ represents the possibility of the entity staying in state $i$, or the portion of occurrences in state $i$. Subscript $N$ is the total number of discrete points. If $P (\theta_i)$ = 1 for certain $i$, meaning that all the fibers of in the image are oriented towards $\theta_i$, such a system can be highly ordered, or has zero information entropy. It can be proven that the system has the highest information entropy when the angles of fibers are evenly distributed (i.e. totally randomized).

Despite the success in communication system, the discrete information entropy is highly dependent on the total number of discrete points ($N$). When $N$ goes to infinity, $H_N$ also goes to infinity, indicating that the discrete entropy cannot be considered as a discrete approximation of a limiting distribution. To address this problem, Jaynes proposed limiting density information entropy ^[55]:

In this formula, $p(\theta)$ is the probability density, and $m(\theta)$ refers to the limiting of sampling distribution. In our application, all the angles are evenly distributed, i.e. $m(\theta)$ = 1/180. As a result,

It can be further proven that the below important relationship holds:

When $N$ is large enough, $O(1/N)$ is negligible.

To show this, let's consider the evenly sampling case and assume that the distribution function (probability density function) $p(\theta)$ is strictly positive and differentiable. Then, for any partition $a = x_0 < x_1 < ... < x_N = b$ of the interval [a, b], the Eq. 2.2 can be written as:

Hence,

Using intermediate value theorem,

And given the fact that

We have,

When $N$ is large, $p(\theta)-p\left(\xi_{i}\right)$ goes to 0. Moreover, because $p(\theta)$ is differentiable, there exists a constant $K$, such that for any $\xi, \eta \in$ [a, b],

Therefore,

where $M$ is a constant.

Then we have

Hence,

The proof is completed. In fact, the condition that the probability density function is differentiable, can be replaced with Lipschitz continuous, which is very natural for physical phenomenons.

According to this formula, $H$ is indeed the limiting density entropy after removing the influence of number of sampling points. Though $H$ eliminates the sampling number effect, it is very hard to calculate $H$ directly from the distribution. Take into account that when $N$ goes to infinity, the right-hand side of Eq. 2.4 will converge to $H$, enabling the approximation of $H$ by $H_N$. Based on this, we propose the following index by taking the right hand side of Eq. 2.4 into 10-based exponential function:

We name this index as homogeneity index. It can be shown that the homogeneity index has following properties:

1) Relatively stable under change of total number of sampling angles (number of bins for sampling angles).

2) Converges as $N$ gets to infinity.

3) Lies on [0, 1] due to concavity of $\log_{10}(x)$.

4) Rises as the alignment quality of fibers decreases. Especially, when the fibers are aligned, HI = 0; when fibers are totally disordered, HI = 1.

2.1.2.   Cyclic Variance and Normalized Entropy

It is popular to adopt cyclic statistics to characterize the distribution of angles. In cyclic statistics, the angles are represented in complex numbers and the related attributes are calculated in the form of complex expressions ^[49]. One advantage of cyclic statistics is that it is inherently invariant to the periodic translation of 360$^{\circ}$. In general, each angle $\theta$ is transformed to $z$ such that,

Then the sample moment is defined by:

The circular mean angle reflecting the location of average angles, can be calculated by:

The circular variance (CVar), in contrast, is defined by:

CVar is ranging from 0 to 1. By definition, CVar is invariant to the periodic translation of 360$^{\circ}$. However, our system requires 180$^{\circ}$ as the period. Therefore, the definition of $z$ is revised as the following form to satisfy the requirement:

CVar computed based on Eq. (2.20) is adopted in this paper. Similar with HI, it can be proved that CVar converges as number of sampling angles increases.

Let

where $\theta_{k} = \frac{180^{\circ}}{n} k$ and

Now we can prove $|M_{n}-M_{0}|$ goes to $0$ when $n$ goes to infinity. In fact,

Therefore, $M_n$ converges to $M_0$. Consequently, CVar($M_n$) converges to CVar($M_0$).

Another measure employed in literature for angular distribution utilizes normalizing discrete entropy to reduce the influence of bin size ^[53]. Namely,

However, this index is only effective when $N$ is small. It fails to distinguish the difference between two angular distributions with large N, since:

Then,

$H$ is independent of $N$, hence a constant for given distribution. Therefore, $S$ is in fact a function of $N$ and will converge to 1 as $N$ rises for all the distributions except for completely oriented or uniform distribution (in these cases, $S$ is always 0 or 1). When number of bins is large, $S$ will always be close to 1.

As mentioned before, being insensitive to the choice of number of bins reduces the requirement of system resolution for high estimation accuracy. This is because for a given system resolution, the total number of pixels in the image is constant. The average counts of occurrences in each bin is inversely proportional to total number of bins. If the index is invariant to the number of bins, one can purposely select large bin number, hence higher average number of occurrences or a better estimation of probability. Figure 3 shows a computation experiment result for different disorder measures, where a normal distribution between 0$^\circ$ and 180$^\circ$ with average 90$^\circ$ and standard deviation 10$^\circ$ is employed. Discrete information entropy grows logarithmically. $S$ increases fast at smaller bin number and will finally reach 1, while HI and CVar are relatively unchanged. This figure depicts the superiority of HI and CVar. Therefore, HI, comparing with CVar will be adopted in our work as a measure of homogeneity.

2.2.   Polarized light imaging system and acquisition of polarized images

In this manuscript, we used a rotating linear polarized light imaging method to measure the distribution of local fiber orientation in both filament phantom slices and unstained human biopsies. The polarized light imaging system was schematically depicted in Figure 4(a). Non-polarized light from a halogen lamp (OSL2, Thorlabs, US) with a broad spectrum and adjustable intensity was first roughly collimated and passed through a filter (FB69010, Thorlabs, US) of 690 nm center wavelength, 10 nm bandpass. Then, a linear polarizer (LPVISE100-A, Thorlabs, US) was used to polarize the light along the polarizer's axis prior to illuminating the sample. A second linear polarizer was introduced behind the sample to measure the change in light polarization upon propagation through the sample. The two co-polarized polarizers are automatically and synchronously rotated over 180$^{\circ}$ polarization range with a step size of 10$^{\circ}$, driven by two stepper motors. At each polarization angle, a monochrome charge-coupled device camera coupled with a lens (BM-141-GE, JAI, Denmark) was used to record the image of sample. Both stepper motors and the camera were controlled by a customized LabVIEW program.

In this study, we defined the direction parallel to the optical table as the initial zero angle. Then, the whole process of acquiring the polarized images of the test sample was as follows.

● Step 1: two polarizers were co-polarized to zero angle so that their transmission axes were aligned to zero angle.

● Step 2: the test sample was placed in the sample holder perpendicular to the optical path.

● Step 3: an image stack of the test sample with 18 frames was captured by the camera as two polarizers rotated synchronously from 0 to 180$^{\circ}$ with an increment of 10$^{\circ}$, presented as Figure 4(a).

2.3.   Quantification of alignment quality

2.3.1.   Local fiber orientation

To better understand the method, note that ignoring other optical effects, if there is a birefringent sample, such as collagen fibers in a thin tissue section, between the two co-polarized polarizers, the transmitted light will be as:

where $I_0$ is the intensity of the input light before the polarizer, $\theta$ is the relative angle between the polarizer and the birefringent axis of the sample, $\Delta n$ is the sample birefringence, d is the sample thickness and $\lambda$ is the light wavelength ^[56,57]. Considering all values are kept constant over different measurement angles, the transmitted light will be as:

with the maximum achieved at $\theta = 0^{\circ}$ or $90^{\circ}$, i.e., $\theta_{polarizer}$ = $\theta_{fiber}$ or $\theta_{fiber} + 90^{\circ}$. In real tissues, although birefringence is the dominant optical effect due to small thickness of test specimens, other effects, e.g. cylinder scattering, are not negligible. Consequently, the intensity–$\theta_{polarizer}$ curve is characterized by one higher peak followed by one lower peak, which corresponds to $\theta_{fiber}$ and $\theta_{fiber} + 90^{\circ}$ respectively. By locating the global maximum of this curve, the fiber orientation angle $\theta_{fiber}$ can be inferred. Thus, the above-mentioned method in section 2.2 can be used to determine fiber orientation angles in thin sample slices. Since the 2-D fiber orientation is determined by locating the global maximum of the intensity–$\theta_{polarizer}$ curve, multiple images are required to have enough angular resolution to find the peak. Though more data is better for curve fitting, 18 is a convenient number and enough for curve fitting. Therefore, we acquired 18 frames for every test sample.

There were three steps to obtain the local fiber orientation results.

● Step 1: three random regions of interest (ROIs) at a size of 50$\times$50 pixels were selected for the image stack of one sample.

● Step 2: a plot of all the intensity values from the 18 frames of the stack was given for each pixel of the randomly selected ROIs, presented as the red points in Figure 4(b).

● Step 3: the local fiber orientation for each pixel was extracted by tracking the location of the maximum intensity of the fitted intensity response curve, presented as Figure 4(b).

To minimize the influence of systematic measurement error on the intensity response and get reliable values of the fiber orientation, fibers with peak-to-valley height less than 2% of the peak value were excluded in the map.

2.3.2.   Calculation of Homogeneity Index (HI)

With the local orientation thus determined for each pixel of the ROI, a local fiber orientation map was obtained. For further quantification, homogeneity index (HI) was used to assess the level of disorder of the local fiber orientation map, which also corresponds to the homogenization degree of fibers. Then HI is calculated by:

where $\theta_i$ represents the fiber orientation and its probability is calculated by counting the number of orientated pixels in the image; The definition of HI index is inspired by the limiting entropy and has the advantage of being invariant to the sample's placement and the polarizers' rotation spacing ^[55].

2.3.3.   Statistic

Data were analyzed using MATLAB software. The classification ability of HI and CVar to quantify the alignment quality of collagen fibers were evaluated using a linear classifier with leave-one-out cross-validation to reduce the probability of a type I error (false positive). The statistical significance of any observed difference was analyzed by Student's t-test. Statistical significance was defined as a P < 0.01.

3.   Experiments and results

3.1.   Verification experiments with PLA filament phantoms

We utilized Poly Lactic Acid (PLA) filament phantoms to test the performance of the polarized light imaging system and validate that the change in polarization state can be used to compute the alignment quality of fibers within a sample. An electrostatic spinning method was applied to produce PLA filament phantoms with designated orientation and polarization characteristics ^[58,59], which simulated collagen and elastic fiber orientation in biological tissue. The alignment quality of the PLA filaments was adjusted by controlling the speed of the roller of the electrostatic spinning device ^{[60,61,62,63]}. As the rotation speed of the roller increased from 500 revolutions per minute (rpm) to 3000 rpm, the filament alignment became ever more ordered ^[57]. Here, we defined the homogenization of fibers as a state of being homogeneous, where fibers distribute randomly with no preference ^[30]. As the quality of alignment increases, the fibers became more aligned with a preference of local orientation. Thus, an increase in order caused a decrease in homogeneity. In our experiments, a set of PLA samples with increasing alignment quality (decreasing homogeneity) of the identical thickness (5.3$\pm$ 0.7 µm) were produced by electrospinning at the same feeding rate, the same collection time (7 min), but with different rotational speeds of 500rpm, 1000rpm, 2000rpm and 3000rpm, respectively ^[57]. And the diameters of the fibers mainly lie in the range of 1-2.2 µm, which is similar to the diameters of collagen fibers. More details about PLA filements in this study can be found in our previous work ^[57].

3.2.   Results of PLA filament phantoms

3.2.1.   System error

Since the measurement bias are introduced by the light sources that emit partially polarized light and by the polarizer that has synchronous rotation error, calibration is a necessity. To do that, we directly acquired images of a "blank" sample (a PLA phantom embed with no fibers obtained at 0rpm) as the two polarizers synchronously rotate ^[57]. The sum intensity of the images was calculated and plotted with the rotation angle of the two polarizers. As expected, the green curve in Figure 5(a) was flat and its intensity value was close to zero over 180$^{\circ}$ range, which indicated that the overall system error is tiny and thus negligible.

3.2.2.   HI of PLA phantoms

The validation of our system combined with HI to successfully quantify the disorder is implemented on PLA phantoms described in section 3.1. Following imaging of the phantoms, we plotted the intensity over 180$^{\circ}$ polarization range of each pixel and constructed the local fiber orientation maps for each ROI. The averaged intensity plots over the three random ROIs for each PLA phantom were plotted in Figure 5(a). The peak value of the blue fitted curve of the PLA phantom obtained at 3000 rpm was much larger than the red fitted curve of the PLA phantom obtained at 500rpm, indicating that the peak value correlates with the alignment quality of fibers, where the bigger peak value means a more organized structure. Figure 5(b) shows the local fiber orientation map within one ROI of the 3000rpm PLA phantom, while Figure 5(d) shows the map of the 500rpm PLA phantom. From these data, we can also extract the corresponding histograms as indicated in Figure 5(c) and Figure 5(e). Then, HI was calculated for each ROI to assess the level of disorder of the local fiber orientation map, which also corresponds to the degree of homogenization of the filaments in PLA phantoms. For each phantom, the HI value was an average of the three random ROIs. The averaged HI of the PLA phantom at 500rpm (HI = 0.93±0.008) is larger than the PLA phantom at 3000rpm (HI = 0.36±0.003). Figure 6 displays a strong negative correlation between the collecting spinning speeds of the PLA phantom and their corresponding HI.

3.2.3.   Comparison of Different Measures

In section 2.1, HI is mentioned to be advantageous over Shannon's information entropy. Here we adopted PLA phantoms experiment to test this. Figure 7 and Figure 9 shows the dependency of discrete entropy and HI on rotation increments and at different PLA spinning speeds, where the number of bins is the total number of sampling angles of linear motors. The discrete entropy, as shown in Figure 7, grows logarithmically as sampling increment declines. CVar, shown in Figure 8, tends to fluctuate. While HI is insensitive to the change of increments (shown in Figure 9). Even when only six sampling angles, which corresponds to 30$^\circ$ as increment step, are collected, HI is still close to the limiting value. Moreover, CVar and HI can both indicate the level of disorder of color map. Taking into account the fact that HI is always smaller than 1, the PLA phantom obtained at 500rmp is very close to completely random. Finally, by definition, the value of all measures only depends on the probability at each angle interval rather than coding, i.e. the angle value. Therefore, sample placement will not affect HI. This is an important property since many traditional indexes, such as variance, change as the sample is rotated, shown in Figure 2 (c) and (d).

3.3.   Application of our polarized light imaging system combined with HI in human vulvar tissues

Histologically, the homogenization of collagen fiber in the dermis has been used as a diagnostic criterion of VLS ^[27,28,29]. In VLS, significant biochemical alterations occur around the zone of homogenization with changes in the organization and expression of collagen fibers. These changes, especially the alignment quality of fibers, may form the basis for the differentiation of VLS from normal vulvar biopsies without H & E staining. Therefore, after validating the system as described in the preceding subsection, we measured the local fiber orientation distribution in unstained human vulva tissue biopsies.

14 normal vulvar tissue biopsies and 20 biopsies from patients with VLS, presented in Table 1, were obtained under a protocol approved by the Institutional Review Board (IRB) of the Second Affiliated Hospital of Chongqing Medical University (IRB No: 2013KLS002). All the involved women acknowledged and signed written consent for participation in the study. The location and number of biopsies collected were determined based on clinical necessity. As part of the standard clinical protocol, all vulvar specimens were formalin-fixed and paraffin-embedded, and then were sliced to a thickness of 4µm. In order to reduce the possible polarimetric signal artifacts, the unstained samples were chemically dewaxed before polarized imaging. Although formalin fixation has been reported to affect linear retardance and depolarization of myocardium and liver samples, use of this substance for ex-vivo polarimetry studies is acceptable as such influence in the optical polarization properties is negligible ^[64]. According to a recent study using a polarization microscopy to characterize fiber orientation of ocular tissues, the fiber orientation measurements were robust to formalin fixation, with few differences between fresh and 1 hour fixing ^[65]. After the experiments, the sections were stained with hematoxylin and eosin (H & E) to generate a pathological report.

3.4.   Results of Human Vulvar Tissues

The local fiber orientation maps of two random ROIs from one normal tissue sample and one VLS sample are shown in Figure 10(a) and Figure 10(c), respectively. The corresponding histograms are shown in Figure 10(b) and Figure 10(d). The maps from the diseased tissue showed a much larger variation in local fiber orientation angle, suggesting low alignment quality of collagen fibers. Similarly, this is reflected in the broad histograms for the diseased tissue with the more uniform distribution than those in normal tissue. The calculated HI of the ROI from the diseased tissue (HI = 0.98±0.006) is larger than the normal tissue (HI = 0.88±0.002). The box plot of averaged HI from three random ROIs for each vulvar specimen was shown in Figure 11(a). Although the value difference is not large, the distribution difference between diseased tissues and normal tissues is obvious (P < 0.01), which clearly demonstrates the differences between normal and diseased tissues. To quantify how accurately this polarized imaging diagnostic method can discriminate VLS from normal vulva, we created linear discriminant models for HI. The sensitivity and specificity of the models was computed using leave-one-patient-out cross validation, with the receiver operating characteristic (ROC) curves shown in Figure 11(b). The area under the curve (AUC) was 0.98 and the sensitivity and the specificity were 93% and 95%, respectively. In addition, it is worthwhile to mention that CVar also achieves the same accuracy in classifying VLS and normal tissues, illustrated in Figure 11(d). Therefore, in terms of characterizing disorder of angular distribution, the performance of HI is comparable with the circular statistics.

4.   Discussion

In this paper, we presented a rotating linear polarized light imaging system combined with a novel parameter HI based on entropy to quantify the alignment quality (homogeneity) of collagen fibers in biological tissues. First, the local fiber orientation was characterized in individual pixels in an image stack by finding the location of the maximum value of the fitted intensity plot over 180°. Then, the resulting map of local fiber orientations was further processed to plot a histogram of the orientations within the ROI area. Finally, based on this histogram, one can readily calculate the alignment quality of the orientation map. In the previous work, spread (standard deviation) has already been implemented for quantification of collagen fiber alignment ^[33]. However, this simple parameter is not an appropriate quantitative measure of fiber orientation distributions in our case since it is highly dependent on the placement of sample, illustrated in Figure 2. The change of sample's placement under microscope can alter the spread make it inapplicable for quantitative measurement of fiber orientation distributions. Therefore, it is critical to develop a more appropriate statistical method to distill the complex information about the degree of fiber homogenization down into a single number. To enable the clinical application, it is also essential that such index should not be influenced by system implementation and precision of device. One relevant study is to apply normalized entropy in orientation maps ^[53]. However, such index is only effective for large motor step size and may lose its generality for small case. Specifically, when the step size decreases, the index converges to one for any samples. The other solution is to employ circular variance (CVar). Nonetheless, the application of CVar has a high requirement for the circular distribution and it may fluctuate significantly in clinical practice, shown by Figure 8. Inspired by the concept of limiting entropy, HI was proposed in this manuscript to rate the isotropic level of fibers (ranging from 0 to 1). The index is constant under the rotation of sample orientation for different tests and proved to be insensitive to the change of sampling increment (angular step size). Those properties guarantee the stable HI value for the same sample regardless of variations in motors and sampling increment, also reducing the requirement of system resolution. In order to better illustrate HI, one application was presented in our manuscript to quantify the disorder of the collagen fiber orientation maps of VLS tissues and normal tissues for ex-vivo diagnosis. A polarized light imaging method, which works best for thin sample slices, was used to obtain the collagen fiber orientation maps. Besides the polarized light imaging system in this work, HI can be applied to many other methods that acquire the collagen fiber orientation maps, such as SHG. In addition, HI can also be found useful in many other scenarios and is not limited to quantifying orientation.

Our system's ability to assess the alignment quality (homogeneity) of fibers was validated using a set of synthetic phantoms, presented in section 3.1. In a recent study, we found that as the rotational speed of the collection cylinder increased ^[57], the alignment quality of the PLA phantom increased, and the corresponding homogeneity decreased. Since Figure 6 reveals the negative relationship between the HI and the spinning speed of the collection cylinder, HI can reflect the alignment quality (homogeneity) of the PLA phantom. Though the results verified the good performance of our method, limitations of this study still exist. In the PLA phantom experiments, we only practiced the influence of fiber organization to polarization state, while in a real tissue section, much more complex scatters of both cells and organelles also exist. In the future study, we will improve our phantom with both spherical and cylindrical scatters, which would be beneficial to fully simulate structural and optical properties of skin tissues.

The ability of our system to differentiate VLS was verified by the clinical application. Sensitivity of 93% and specificity of 95% are obtained to characterize VLS in vulvar biopsies without H&E staining. Compared with other optical systems for investigating VLS, our system is advantageous in many aspects. For example, Raman spectroscopy (RS) had been applied to the ex-vivo study of LS ^[66]. According to the study by Jonathan et al., RS was able to correctly differentiate VLS from other vulvar conditions in biopsies with a sensitivity of 91% and specificity of 80% ^[66]. While the RS system has a better resolution and is able to evaluate the biomolecular composition of tissue by assessing the vibrational energy levels of the constituent molecules, our system provides a simpler and lowcost method for VLS characterization with higher sensitivity and specificity. Moreover, our method provides a potential quantitative diagnosis tool for VLS and can promote understanding the roles of collagen fibers alignment during malignant transformation of VLS. One may notice that in this work, the fiber orientation of the sample is assumed to be 2D single layer rather than 3D multi-layer. This is a reasonable assumption because when preparing the biopsy, it is required that the sample should be sliced along the major orientation plane. Additionally, the thickness of sample slices (4µm) is very small compared to its length scale (1mm). Such thin sample condition guarantees that the fiber with higher inclination has smaller contribution to the imaging result. With a thickness of 4µm, the slices only contains less than 3 layers of fibers. Consequently, 2D single layer simplification is reasonable in this study. There are several limitations of the current work. First, our transmission co-polarized system works best for thin samples. If the sample is thick, reflectance-based imaging will be needed. In the future works, HI will be applied with different imaging methods to quantify the disorder of the images, which depends on the sample. Secondly, the sample size is not very large, and the control group only includes normal vulvar tissue samples. As a preliminary study, the relatively high sensitivity and specificity results from 34 samples can reveal the feasibility of our method to distinguish diseased tissues from normal tissues. However, a lager sample size with additional biopsies from patient with other vulvar conditions, i.e., vulvar intraepithelial neoplasia (VIN) and vulvar squamous cell carcinoma (VSCC) will be needed for further validation of the clinical performance in the future. Thirdly, we only consider the application of HI in 2D problems. Future study involves exploring the application of HI in 3D problems.

5.   Conclusions

This study demonstrates the feasibility of a novel parameter HI based on limiting entropy to quantify the disorder of an image. Compared with the performance of other parameters, such as cyclic variance, discrete entropy and normalized entropy, our proposed HI is more robust to the placement of the sample slice and the systematic sampling increment (angular step size). Our attempt of using this method to differentiate VLS and normal vulvar tissue in biopsies is very promising. Our method possesses several advantages with significant clinical impact: 1) the diagnosis results are objective and quantitative; 2) no need to stain the samples so that the diagnosis is faster and less expensive. Existing limitations of the study include a small number of test samples and a relatively low resolution of our system, so further study is needed. Future studies will focus on increasing the number of test subjects, exploring the effect of additional pathologies on the fiber alignment, and acquiring images at a higher resolution. Implementing HI in a clinical setting may lead to quantitative characterization and accurate diagnosis of many diseases in the future.

Acknowledgments

We gratefully acknowledge research support from the Natural Science Foundation of China (Nos.81271527 and 81327803), the Fundamental Research Funds for the Central Universities (WK2090090013), and Chongqing Health and Family Planning Committee (2015ZDXM010). We also acknowledge Dr. Bin Yang for his helpful technical suggestions.

Conflict of interest

Conflict of interest