2.1.

Phantom Design

The purpose of a phantom study is to validate the segmentation method on phantoms made by well-defined layers of materials with different birefringence. Two phantoms are developed. Phantom 1 is made by tightly taping three layers of removable tape (Scotch^® Magic™, 3M, 811) around a piece of cover glass (Fisherbrand™, 12-542-C). Removable tape has moderate birefringence, whereas cover glass has no birefringence. The thickness of a single layer of the tape is measured to be $53\pm 1\mu \mathrm{m}$ by an electronic digital caliper (Marathon, CO 030200). Phantom 2 is made by covering two layers of removable tape on top of chicken breast muscle. The birefringence of chicken breast muscle is typically two to three times higher than that of tape.

2.2.

Cartilage Tissue Preparation

Articular bones from the knee joint of adult swine are obtained from a local farm and imaged within 24 h. Cartilage tissue blocks are extracted from the articular bone and cut into cubic blocks with edge length of around 3 to 5 mm. The top of the cubic tissue is the intact natural surface. In this study, PS-OCT imaging on five samples (S1 to S5) extracted from four animals is reported. The sample locations are illustrated in Fig. 1(a). All the samples are extracted from the distal femoral condyle. Samples S1, S2, S4, and S5 are extracted from the medial central location (position 1) of four different animals. Sample S3 is extracted from the lateral central location (position 2) of the same animal as S2. Thus, samples S2 and S3 compare two symmetric locations on the femoral condyle.

Fig. 1

(a) Photograph of the articular bone sample. The numbers “1” and “2” mark the locations where the cartilage tissue blocks are extracted. (b) Laser illumination and image scanning configuration. Laser illuminates from the top on the tissue surface, and $\theta$ is illumination angle. In PS-OCT imaging, each A-line is along the $Z$ direction, a B-scan image is in the XZ plane (large blue rectangle), and multiple B-scans are acquired along $Y$ axis. In SHG imaging, each frame is in the XY plane (small blue rectangles), and an image stack can be acquired at different depths.

2.3.

Polarization-Sensitive Optical Coherence Tomography Imaging

PS-OCT imaging is conducted on a previously reported PS-OCT system.¹⁰ The schematic of the system is shown in Fig. 2. The system is based on a swept-source laser (Axsun Technologies Inc., Massachusetts), with the center wavelength of $1.06\mu \mathrm{m}$, full width at half maximum bandwidth of 110 nm, and scanning rate of 100 kHz. The light from the source is split by a 90:10 single-mode optical fiber coupler, where the 90% port of the coupler is connected to a passive polarization delay unit, and the 10% port to the reference arm. With the passive polarization delay unit, two orthogonal polarization states of light are separated and delayed in time. The light then passes through a 75:25 fiber coupler, and the 25% port is sent to illuminate the sample. The backscattered signal from the sample passes through the fiber coupler again, and the 75% port is sent to the detection unit. Light from the sample arm and the reference arm are guided to a polarization diversity (PD) detection unit. The polarization sensitive features are measured by a Jones-matrix-based method.^8,10 The average power on sample is around 2 mW. The depth resolution is $\sim 8.1\mu \mathrm{m}$ in air, and the lateral resolution is $\sim 19.2\mu \mathrm{m}$.

Fig. 2

Schematic diagram of PS-OCT system. LP, linear polarizer; PC, polarization controller; GM, galvanometer scanner; BS, beam splitter; PBS, polarization beam splitter; BPD, balanced photodetector.

The illumination and scanning geometry of PS-OCT imaging are shown in Fig. 1(b). The laser beam illuminates on the natural tissue surface from top. The incident angle $\theta$ is the angle between the laser beam and the normal of the tissue surface. Sample is attached on a tilting mount so that the incident angle $\theta$ can be adjusted in the azimuth plane (XZ plane). At near normal incidence, $\theta$ is close to zero, and the beam is almost perpendicular to the tissue surface (XY plane). At tilted illumination, $\theta$ is rotated around $Y$ axis. PS-OCT scanning is performed by scanning the laser beam in $X$ (fast scan) and $Y$ (slow scan) directions, with which a 3-D volume can be acquired. Each B-scan image shows a cross-sectional view in the XZ plane.

In the PS-OCT system, two orthogonal input polarization states are illuminated on the sample, and two PD detectors measure the corresponding OCT output signals. Thus, the OCT output signals at each depth location $z$ can be written as a $2\times 2$ matrix:

A surface segmentation is then applied to obtain the depth position $z_0$ of the sample surface and the corresponding $E_{\text{out}}(z_0)$. A similar matrix $M$ of the round-trip Jones matrix of the sample can be constructed as

The round-trip phase retardation of the sample, $\delta (z)$, is then obtained as the phase difference between the two eigenvalues as

2.4.

Second-Harmonic Generation Imaging

Details about the SHG imaging system can be found in Ref. 21. A femtosecond Ti:sapphire laser is used to excite the nonlinear signal. A $60\times$ objective lens with $\mathrm{NA}=1.2$ (UPLSAP060XW, Olympus) is used to acquire high-resolution imaging. The lateral resolution of the SHG imaging is $0.4\mu \mathrm{m}$. Each SHG image frame is acquired in the XY plane as marked by the blue rectangles in Fig. 1(b). The field of view is $\sim 130\mu \mathrm{m}\times 130\mu \mathrm{m}$. Two image stacks are acquired along the tissue depth ($Z$ axis) by depth scanning. One stack is taken from the natural tissue surface into $\sim 200\text{-}\mu \mathrm{m}$ depth, and the other stack is taken after a $\sim 0.5\mathrm{mm}$ of the top layer of the tissue is removed.

4.1.

Data Processing of Slope-Based Segmentation

The segmentation algorithm is demonstrated on cartilage sample S1 as shown in Fig. 5. Figures 5(a) and 5(b) show the intensity and phase retardation image, respectively. The phase retardation is calculated based on the method described in Sec. 2.3. Details about the data processing are explained below.

Fig. 5

Segmentation results of articular cartilage sample S1 from the slope-based analysis in PS-OCT. (a) Intensity image and (b) phase retardation image. The 2-D segmented boundaries are marked by dashed lines and labeled as $S$, B1, B2, and $T$. (c) Averaged A-line profile of phase retardation, with fitting result shown in red dashed line. Blue dashed line indicates the reliable thresholding depth. (d) 3-D rendering of the segmentation results. Arrow in (b) indicates the location of the A-line shown in (c). Scale bar is $500\mu \mathrm{m}$.

Preprocessing: In Fig. 5(a), the tissue surface is clearly shown in the intensity image due to strong backscattering. The tissue surface is identified based on the intensity image and marked as $S$ in the coregistered phase retardation image, as shown in Fig. 5(b). The origin of each A-line is defined at the surface of the tissue. Each A-line in the phase retardation image is then shifted up to flatten the surface. Running averaging is performed on each A-line by averaging it with its surrounding A-lines. The averaging window is selected as 21 A-lines by 11 B-scans by trial and error. Figure 5(c) shows a representative averaged phase retardation A-line, which shows a characteristic feature of multiple linear regions along the depth. From Fig. 5(c), three regions can be observed. The first region has a relatively flat slope and is located at the depth $<200\mu \mathrm{m}$; the second region has a slow slope and is located around the depth 200 to $400\mu \mathrm{m}$; and the third region has a fast slope and is located at the depth $>450\mu \mathrm{m}$. A multisegment linear fitting will be used to fit the phase retardation A-line.

Reliability thresholding: In the intensity images, SNR drops significantly at deep depth due to light attenuation, which can affect the reliability of the analysis. Effective signal-to-noise ratio (ESNR) is specific for Jones-matrix-based PS-OCT imaging,²² which balances the SNR in the four channels acquired using two input polarization states and two PD detectors ($2\times 2=4$). To ensure reliable calculation, thresholding is applied and pixels with ESNR $<10\mathrm{dB}$ in the phase retardation image are removed.

Slope-based segmentation: Segmentation of the multiple linear regions is demonstrated on the phase retardation A-line shown in Fig. 5(c) as an example. In the segmentation process, the pixel positions $S$, B0, B1, B2, and $T$ are solved for each phase retardation A-line, where $S$ is the tissue surface, B0 and $T$ are the starting and ending point, respectively, of the reliable phase retardation accumulation, B1 is the boundary between the superficial and transitional zone, and B2 is the boundary between the transitional and deep zone. Several assumptions are applied in the algorithm based on prior knowledge about articular cartilage: articular cartilage contains three structural zones, the superficial zone has a thickness $<200\mu \mathrm{m}$, and the PS-OCT imaging depth can reach the deep zone. Here, B0 is introduced to exclude the first, a few pixels which are highly fluctuating due to surface effect and artifacts. The following steps are taken in the segmentation algorithm:

From the segmentation, three regions with different phase retardation slopes can be distinguished as R-1 ($S-\mathrm{B}1$), R-2 (B1–B2), and R-3 ($\mathrm{B}2-T$), corresponding to superficial, transitional, and deep zone, respectively.

Quantification of birefringence: The depth position and thickness of each region can be extracted and quantified from the segmentation. Birefringence is quantified by the slope of the phase retardation curve. For tissue of thickness $z$, the round-trip phase retardation accumulation can be expressed as $\delta (z)=4\pi \mathrm{\Delta }n·z/{\lambda }_0$, where ${\lambda }_0$ is the central wavelength. From the linear fitting, the slope of the phase retardation accumulation can be extracted. The slope of the phase retardation is related to the birefringence $\mathrm{\Delta }n$, as $\text{slope}=\frac{\mathrm{d}\delta (z)}{\mathrm{d}z}=\frac{4\pi \mathrm{\Delta }n}{{\lambda }_0}$.^9,14 In this study, the optic axis of the collagen fiber is unknown. Thus, the obtained birefringence should be apparent birefringence, which will be discussed in more details later. By analyzing the accumulated phase retardation and quantifying its slope over a sufficient depth range, the accumulative properties of birefringence will be obtained and used to distinguish different structural zones.

4.2.

Slope-Based Segmentation in Articular Cartilage

Figure 5 shows the PS-OCT imaging and segmentation results on sample S1 near normal illumination. The intensity and phase retardation images are shown in Figs. 5(a) and 5(b), respectively. Although the intensity image can show the tissue surface, it cannot differentiate the different structural zones. The phase retardation image shows the accumulated phase retardation along the depth. The color mapping corresponds to phase change over $0-\pi$. The band-like features indicate how the accumulated phase retardation changes, which cannot be observed in the intensity image. A representative phase retardation A-line from the B-scan image is shown in Fig. 5(c). The multisegment linear fitting results are shown as the red dashed segments. The boundaries $S$ (surface) and $T$ (maximum reliable depth) and the turning points B1 and B2 are marked. Three regions are segmented based on the multisegment linear fitting of the phase retardation A-line. From $S$ to B1, the phase retardation has not much accumulation in this short region near the surface (typically around $100\mu \mathrm{m}$). From B1 to B2, the phase retardation starts to accumulate gradually. From $\mathrm{B}2-T$, the phase retardation increases at a relatively fast speed. The segmentation algorithm is applied to a 3-D image volume (300 A-lines in the fast axis and 200 B-scans in the slow axis, around $3\mathrm{mm}\times 4\mathrm{mm}$) to segment the zones. The segmentation results are rendered in 2-D and 3-D, respectively, in Figs. 5(b) and 5(d). The layers in the 3-D rendering show the zone divisions. The global features match with the geometrical shape of the tissue, with some topographically dependent variations.

Similar segmentation results have been obtained from other samples. Figure 6 shows the results from sample S2 (left column) and S3 (right column). From top to bottom, the intensity image, phase retardation image, and phase retardation A-line are shown, respectively. The segmented boundaries are marked by the dashed lines in the phase retardation images. Comparing Figs. 5 and 6, we can see that samples S2 and S3 show similar features as S1, where the superficial zone has a relatively flat slope and the deep zone has a faster slope than the transitional zone.

Fig. 6

Segmentation results from articular cartilage sample S2 (left column) and S3 (right column). (a) Intensity image and (b) phase retardation image. Segmented boundaries are marked by dashed lines. (c) Averaged A-line profile of phase retardation, with the fitting result shown in red-dashed line. Blue-dashed line indicates the reliable thresholding depth. Arrow in (b) indicates the location of the A-line shown in (c). Scale bar is $500\mu \mathrm{m}$.

The averaged thickness and slope of each segmented region from the 3-D volume of the samples S1 to S3 are summarized in Table 1. Thickness is converted from optical path length by assuming an average refractive index $n=1.444$ in articular cartilage.²³ For sample S1, region R-1 is relatively thin with a thickness of $\sim 125\mu \mathrm{m}$. Region R-2 is twice thicker than R-1, which is consistent with the anatomical knowledge about the thickness of superficial and transitional zones.^1,2 Region R-3 is $\sim 560\text{-}\mu \mathrm{m}$ thick, and it may not represent the full thickness of the deep zone as the maximum reliable depth is limited by the light penetration. In relatively deep depth ($>1\mathrm{mm}$), the phase retardation computation becomes less reliable due to low SNR.²⁴ All the three samples have similar thickness division over the three regions. Especially, samples S2 and S3 show very close thickness in both regions R-1 and R-2. This indicates that the symmetric locations on the same cartilage sample have similar structural zone thickness.

Table 1

Averaged thickness and phase retardation slope of segmented regions in samples S1, S2, and S3.

In Table 1, the slopes of samples S1, S2, and S3 show some common characteristics. The slope increases from region R-1 to R-3. The slope in region R-1 is close to zero, with a large standard deviation that is higher than its mean value. This indicates that the phase retardation cannot accumulate sufficiently in region R-1, due to the frequent change of collagen fiber orientations in the superficial zone at different depths. The slope in region R-2 is higher than that in R-1 as collagen fibers become obliquely aligned in the transitional zone. The slope in region R-3 is about 3 to 4 times higher than that in R-2. Other cases have also been observed where the slope in the deep zone is relatively low, which are likely caused by variations in the topographical locations and illumination angles. More discussions will be provided in Sec. 5.

In our investigation, region R-2 normally shows higher fitting error than the other regions. For the phase retardation slope, the relative standard deviation compared with its average value is also higher in region R-2 than in region R-3. Those information may imply that the collagen orientation in the transitional zone is more random than that in the deep zone. The standard deviation in the superficial zone is high, likely due to a combination of several reasons, such as uneven surface, strong surface reflection, and large fiber bundles with alternating orientations at different depths.

4.3.

Validation with Histology Image

The PS-OCT segmentation is also compared with histology. In this experiment, sample S4 is imaged by PS-OCT. Afterward, the sample is sent to a histology lab to obtain the histology image. The PS-OCT phase retardation image is shown in Fig. 7(a). Three regions R-1, R-2, and R-3 can be differentiated and segmented based on the slope of the phase retardation measurement. The segmentation boundaries are marked by the dashed lines. The histology image of the same sample is shown in Fig. 7(b). Although the histology image does not have sufficient resolution to resolve collagen fibers, the shape and orientation of chondrocytes can be visualized, showing distinct features in different structural zones. Based on the averaged segmentation depth identified in the phase retardation image, the regions are marked by the brackets on the left side of Fig. 7(b). From the regions R-1, R-2, and R-3, different features of the chondrocytes can be observed. Zoomed-in images from a representative area selected from each region are shown in Figs. 7(c)–7(e), respectively. In Fig. 7(c), the chondrocytes are flattened and aligned almost parallel to the tissue surface. In Fig. 7(d), the chondrocytes are more rounded than in Fig. 7(c), and oriented obliquely. In Fig. 7(e), the chondrocytes are aligned in columns (the dashed lines indicate the approximate alignment of the chondrocytes). Those features in each selected region match with the common features of chondrocytes in the superficial, transitional, and deep zone of articular cartilage, respectively.^1,25 Thus, with the evidence from the histology results, the three segmented regions correspond to the superficial zone (R-1), transitional zone (R-2), and deep zone (R-3).

Fig. 7

Comparison of PS-OCT with histology. (a) PS-OCT phase retardation image with segmented boundaries marked by dashed lines. (b) Histology image of the same sample. Scale bars in (a) and (b) are $100\mu \mathrm{m}$; (c–e) Zoomed-in images from the ROIs indicated by the squares in (b). Red-dashed lines in (e) indicate the alignment of chondrocytes in columns. Scale bars in (c), (d), and (e) are $20\mu \mathrm{m}$.

4.4.

PS-OCT Imaging Under Various Illumination Angles

PS-OCT measurement has been reported to be sensitive to the illumination angle with respect to the collagen fiber orientation.^13,14,24,26 A multiple illumination angle study is performed on Sample S5 to further verify the multisegment slope-based segmentation in PS-OCT. Figure 8 shows the PS-OCT intensity and phase retardation images at illumination angles of $-40\mathrm{deg}$, $-35\mathrm{deg}$, $-25\mathrm{deg}$, $-10\mathrm{deg}$, 0 deg, 10 deg, 15 deg, 30 deg, and 40 deg. The banding structure in the color mapping can be used as a qualitative criterion, where tighter bands indicate higher tissue birefringence. The figures show that the tissue birefringence becomes higher as the illumination angle increases from 0 deg to larger angles.

Fig. 8

PS-OCT images measured at different illumination angles from $-40\mathrm{deg}$ to 40 deg. In each B-scan image panel, the intensity OCT image is on the top, and the phase retardation image is at the bottom. The arrows indicate the locations of the A-lines. Scale bar is $500\mu \mathrm{m}$.

To verify how the slope-based segmentation can also be applied to tilted illumination, three A-line profiles and the corresponding multisegment linear fitting are shown in Figs. 9(a)–9(c), at illumination angles of around 0 deg, $-25\mathrm{deg}$, and $-40\mathrm{deg}$, respectively. The locations of the corresponding A-lines are marked by the arrows in the phase retardation images shown in Figs. 8(a), 8(c), and 8(e), respectively. In Fig. 9, the multisegment linear fitting works very well, and three regions are successfully segmented for each case at different illumination angles. From the surface to deep location, the phase retardation slope shows a minimum slope in region R-1, a slow slope in region R-2, and a fast slope in region R-3. Furthermore, the slope in region R-3 is sensitive to the illumination angle and it increases monotonically as the illumination angle increases from 0 deg to $-25\mathrm{deg}$ and $-40\mathrm{deg}$.

Fig. 9

Phase retardation A-line profile at different illumination angles of $\theta =0\mathrm{deg}$, $-25\mathrm{deg}$, and $-40\mathrm{deg}$, respectively. Fitting result is shown in red-dashed line. Blue-dashed line indicates the reliable thresholding depth.

The segmentation algorithm is applied on 3-D volume data at different illumination angles. The averaged segmentation depth and slope for each region are shown in Figs. 10(a) and 10(b), respectively. The averaged segmentation depth for B1, B2, and $T$ are $\sim 100$, 500, and $1100\mu \mathrm{m}$, respectively. Thus, the thickness of R-1, R-2, and R-3 are $\sim 100$, 400, and $600\mu \mathrm{m}$, respectively. In tilted illumination, light travels through a slightly longer path length in a layered structure compared with normal illumination. For the maximum tilting angle of $\pm 40\mathrm{deg}$, this difference is calculated to be $\sim 10\%$ by considering refraction. In Fig. 10(a), this small difference is not observed in the segmented layer thickness at different illumination angles likely because the variation of the structural zone thickness caused by topographical variation and tissue heterogeneity is $>10\%$. Thus, in this study, the thickness is defined in the direction of the beam illumination for simplicity. At large tilting angles, the SNR is lower and the penetration depth is reduced, which can also be seen in the intensity images.

Fig. 10

(a) Averaged boundary depths and (b) averaged phase retardation slopes of the segmented regions as a function of illumination angle.

The averaged slopes of the different regions are shown in Fig. 10(b). Region R-3 has a larger slope than regions R-1 and R-2. It indicates that region R-3 has higher tissue birefringence. This is consistent with the anatomical fact that collagen fibers in the deep zone are densely packed and highly aligned.^2,25 The slope of region R-3 shows high sensitivity to the illumination angle. At 0 deg to 15 deg illumination angle, the slope is minimum at $\sim 0.15\mathrm{deg}/\mu \mathrm{m}$. At $\pm 40\mathrm{deg}$, the slope is increased to $\sim 0.4\mathrm{deg}/\mu \mathrm{m}$.

To understand the birefringence measured at different illumination angles, collagen fiber can be modeled as a uniaxial crystal. In a uniaxial crystal, the identical principal indices of refraction are called the ordinary index, $n_o$, and the distinctive index of refraction is called the extraordinary index, $n_e$. The index of refraction for the ordinary wave is simply $n_o$. The index of refraction for the extraordinary wave depends on the angle ${\theta }_c$ and is given as²⁷

In the deep zone, collagen fibers are aligned approximately in the perpendicular direction to the surface. When the illumination angle increases, a larger slope is observed, indicating higher apparent birefringence. As the collagen fibers may not be perfectly aligned and exactly parallel with the illumination light even at 0 deg tilting angle, some birefringence is still observed near normal illumination in our experiments. The slopes in regions R-1 and R-2 are found to be not sensitive to the illumination angle. This is likely because the collagen fibers in the superficial zone change orientation at different depths, and the collagen fibers in the transitional zone are bent with some randomness.

4.5.

MPM Imaging of Articular Cartilage

Although PS-OCT does not have sufficient resolution to directly visualize collagen fibers, SHG imaging is a high-resolution microscopy technique that can visualize the organization of collagen fibers. SHG imaging is carried out on sample S1 to examine the microscopic features at different depth regions. The microscopic features from the SHG imaging are then compared with the birefringence features measured by the PS-OCT for each structural zone.

Figure 11 shows the SHG images from approximately the three depths marked by 1, 2, and 3 in Fig. 5(b). Figures 11(a)–11(d) show SHG images at depth range 10 to $40\mu \mathrm{m}$, near the tissue surface in region R-1. Large fiber bundles are clearly distinguished, and the fibers are almost parallel to the tissue surface. The fiber orientation changes from frame to frame at different depths. Those features match with the collagen properties in the superficial zone. With the fiber orientation changing frequently within a short depth range, the phase retardation cannot accumulate and the measured birefringence in this layer would be low. This is consistent with the relatively flat slope in region R-1 in the PS-OCT slope-based analysis. Figures 11(i)–11(l) show SHG images at depth around 500 to $600\mu \mathrm{m}$, inside region R-3. The collagen in this region is found to form small fibrils that are densely packed. The dark circular holes correspond to chondrocytes that lack SHG contrast. Those features are consistent with the collagen properties in the deep zone, where type II collagen forms fibrils and the collagen density is high.¹ The orientation of the fibrils is hard to distinguish as limited by the resolution of the SHG imaging system. Figures 11(e)–11(h) show SHG images at depth range around 80 to $140\mu \mathrm{m}$, in the transitional area from R-1 to R-2. The fibers are much smaller than those in the superficial zone, and they are not orientated in parallel with the surface, an indication of oblique orientation. The collagen density seems to be lower than that in the region R-3. Those features are consistent with the properties of the transitional zone. Thus, with the supportive evidence from high-resolution SHG imaging, the segmentations R-1, R-2, and R-3 from the PS-OCT analysis are shown to match with the characteristics of the superficial, transitional, and deep zone, respectively, in articular cartilage.

Fig. 11

SHG images of articular cartilage. (a–d) SHG images at depths around 10, 20, 30, and $40\mu \mathrm{m}$, in the superficial zone. (e–h) SHG images at depths around 80, 100, 120, and $140\mu \mathrm{m}$, in the transitional zone. (i–l) SHG images at depths around 500 to $600\mu \mathrm{m}$ with $20\text{-}\mu \mathrm{m}$ step size. Scale bar is $30\mu \mathrm{m}$.