Automated quantification of three-dimensional organization of fiber-like structures in biological tissues.

2.1. 3D directional variance formalism to quantify fibrillar organization The development of 3D directional variance relies on concepts from directional statistics [ 27 ]. Directional variance for circular data (for which orientations of 0° and 360° are equivalent), V 3 D , is defined as [ 27 ]: (1) V 3 D = 1 − R ¯ 3 D where R ¯ 3 D = ( C ¯ 3 D c 2 + S ¯ 3 D c 2 + Z ¯ 3 D c 2 ) 1 / 2 , and C ¯ 3 D c = ( 1 / n ) ∑ j = 1 n sin φ j cos θ j , S ¯ 3 D c = ( 1 / n ) ∑ j = 1 n sin φ j sin θ j , Z ¯ 3 D c = ( 1 / n ) ∑ j = 1 n cos φ j . The superscript c refers to the circular data, and n is the total number of voxels that contribute to the determination of the directional variance. θ and φ are the azimuthal and polar angles used to depict an orientation in 3D space, respectively ( Fig. 1a ). The azimuthal angle θ is determined by projecting fibers to a fixed 2D plane, i.e., the xy plane ( Fig. 1a ). However, according to the definition of the polar angle, φ , there is no fixed 2D plane onto which fibers could be projected for the determination of all possible φ orientations. To address this problem, we consider additional angles β and γ ( Fig. 1a ), which are azimuthal angles, like θ , and related to φ by [ 26 ]: (2) tan 2 φ = 1 / tan 2 β + 1 / tan 2 γ The determination of β and γ can be achieved using the same approach as that for θ. After projecting fibers to a fixed 2D plane, the orientation is determined from the weighted vector summation algorithm we described recently [ 18 , 26 ]. Briefly, to acquire the θ orientation of the center voxel of a n × n × n voxel cube window, we first average the cube along the z direction to project the cube window to a square window ( n × n pixels), and then calculate the orientation of the center pixel of the square, based on the basic premise shown in Figs. 1b and c . First, we define the vectors passing through the center pixel (marked in purple) of the window size of choice (11×11 pixels in this case). Then these vectors are weighted by their length and intensity fluctuations along their direction [ 26 ], as shown in Fig. 1b . The orientation of the center pixel is defined as the summation of all these weighted vectors and shown in Fig. 1c , which corresponds well to the direction of fiber alignment. The determination of angles β and γ is achieved by averaging the cube along the other two directions, followed by the identical weighted vector summation algorithm. Once we acquire the values of θ and φ for each voxel of a 3D stack containing fiber-like structures, whose orientations are equivalent for 0° and 180°, we need to modify the definition of 3D directional variance. A general approach to transforming axial (fiber like) data to circular data is to multiply angular values by 2. This strategy is sufficient for transforming the azimuthal angle θ , but not for the polar angle φ . Again, we use azimuthal angles β and γ to represent the polar angle φ , and define b j = 1 / tan 2 ( 2 β j ) + 1 / tan 2 ( 2 γ j ) to derive the modified components as: (3) C ¯ 3 D a = ( 1 / n ) ∑ j = 1 n ( b j / 1 + b j 2 ) cos ( 2 θ j ) (4) S ¯ 3 D a = ( 1 / n ) ∑ j = 1 n ( b j / 1 + b j 2 ) sin ( 2 θ j ) (5) Z ¯ 3 D a = ( 1 / n ) ∑ j = 1 n S I / 1 + b j 2 Where SI = (−1) · ( φ − 90)/| φ − 90| when φ ≠ 90°, and SI = 1 when φ = 90°. The superscript a refers to the axial data. Using these modified components, we acquire the modified R ¯ 3 D , and finally the V 3 D that is suitable to quantify the organization of 3D fiber-like structures. The custom code for the assessment of the 3D directional variance is written in MATLAB. Orientation maps of two typical planes, as marked by green and blue ( Fig. 1a ), are plotted in Figs. 1d and e respectively, which serve as a reference for assessing the orientation values. Throughout this study, optical sections are all along the xy plane, and 3D image stacks are then reconstructed based on these sections at different z depths. Note that a φ of 90° corresponds to fibers that are parallel to the optical sections. In this study, the 3D directional variance is assessed typically based on two types of analysis windows, a localized and a full image stack window. The size of the localized window is dependent on characteristics of the fibers, such as their diameter, waviness and tissue level organization. Every voxel within an image stack is assigned a directional variance value assessed based on this localized window. When applying the full image stack window, we characterize the 3D variance of all the fibers within the stack with a single value. The latter calculation is faster, but in some cases more detailed information may be needed to characterize important functional changes that may be present at smaller scales. 2.2.

Generation of simulated fiber stacks

To validate the 3D directional variance algorithm, 3 simulated fiber stacks with a size of 200×200×200 voxels are generated with different levels of fiber alignment using MATLAB. There are ~105 fibers in each stack. The first stack has completely parallel fibers with a θ orientation of 70°, and a φ orientation of 120°. An intermediate fiber organization is generated in the second stack, where half of the fibers have the θ angle ranging randomly between 40° and 70°, and the φ angle ranging randomly between 100° and 120°, whereas the other half of the fibers have the θ angle ranging between 100° and 165°, and the φ angle ranging between 40° and 80°. The third stack has randomly defined θ and φ angles ranging between 0° and 180°. 2.3. Testing the robustness of the 3D directional variance analysis with respect to changes in fiber density We generate fiber stacks with the same organization level (with a 3D directional variance of ~0.57), but with five different levels of fiber density, and assess the error in determining the fiber orientation angle and the percent error for assessing the 3D directional variance. We examine stacks with fiber densities of approximately 3.6, 4.3, 5.1, 5.6 and 6%, with the 6% stack having a density that is 67% higher than that of the 3.6% one, which is a significantly higher variation range than what we typically observe in the samples we examine ( Supplementary Table 1 ). For each fiber density, we generate 10 fiber stacks. The percent error for the 3D directional variance estimates is calculated by: ( | anticipated variance - calculated variance | ) / anticipated variance. We determine the mean and standard deviation for the error in our estimations of the orientation angle and the percent error of the extracted 3D directional variance at different fiber density levels. 2.4.

Multi-photon microscopy

SHG and TPEF images are obtained using a Leica TCS SP2 confocal microscope equipped with a tunable (710–920 nm) titanium-sapphire laser (Mai Tai; Spectra Physics; Mountain View, CA). Three non-descanned photomultiplier tube (PMT) detectors detect light in the 460 ± 20 nm, 525 ± 25 nm or 400 ± 10 nm regions. Objectives used in this study include a water-immersion 63× objective (NA 1.2; 220 μm working distance) for articular cartilage and mammary glands, a water-immersion 25× objective (NA 0.95; 2.4 mm working distance) for pancreatic cancer and collagen hydrogels, and a dry 20× objective (NA 0.70; 590 μm working distance) for the brain-like cortical tissue. Reconstruction of 3D images is performed by ImageJ (W. Rasband, National Institute of Health, USA), and image analysis is performed by the 3D directional variance algorithm in MATLAB. 2.5. Collagen hydrogel-based tissue preparation Normal Human Lung Fibroblasts (NHLF) are obtained from Lonza (Lonza, Walkersville, MD) and cultivated using FGM-2 BulletKit (Lonza, Walkersville, MD), containing basal media, serum, hFGF-β, insulin and GA-1000. Fibroblasts between passages 4–6 are used for all experiments. To study fibroblast alignment in attached collagen hydrogels, 6 well Flexcell Tissue Train culture plates (Flexcell, Burlington, NC) are employed. Pulmonary fibroblasts are detached from tissue culture flasks using 0.05% Trypsin-EDTA (ThermoFisher Scientific, Cambridge, MA) and trypsin neutralizing solution (Lonza, Walkersville, MD) is added following cellular detachment. Cells are then centrifuged and counted prior to collagen hydrogel preparation. To prepare 1 ml of cell laden collagen hydrogels, high concentration rat tail type-I collagen (Corning, NY) (final concentration: 1 mg/ml) is mixed with cold DMEM (final concentration: 0.5×), 1 M NaOH, 100 μl of pulmonary fibroblast cell suspension (10 million cells/ml) and cold sterile water. Tissue Train culture plates are removed from their sterile wrappers and placed on a Flexcell baseplate containing Trough Loaders. Application of pressure using a vacuum pump connected to reservoir results in deflection of the rubber membrane into a trough format. Within the trough, 200 μl of cell laden acid neutralized collagen is added and incubated inside a tissue culture incubator (37 °C) to form the hydrogel construct. Following polymerization, cell laden collagen hydrogels are afloat FGM-2 media, but remain attached to anchor stems on either side of the wells in the Tissue Train culture plate. To study dynamic changes in cellular and collagen alignment, cell laden collagen constructs are sacrificed at 4, 10, 18 and 30 hours post seeding and fixed using 4% paraformaldehyde for 1 hour at room temperature. Fixed constructs are stored in 1× PBS at 4 °C for further analysis. All hydrogel constructs are brought to room temperature (RT) on the day of immunostaining. Briefly, all constructs are washed three times for 20 min in PBST (1× PBS + 0.1% Tween-20). Blocking is performed for 1–2 hours at room temperature (RT) using PBST containing Normal Goat Serum (NGS; Vector Laboratories, Burlingame, CA), Bovine Serum Albumin (BSA) and 0.1% TritionX-100. Following blocking, all constructs are incubated in 1:1000 rabbit monoclonal anti-vimentin antibody (Abcam, Cambridge, MA) in PBST containing NGS and BSA for 24–48 hours. Secondary DyLight 488 conjugated goat anti-rabbit IgG antibody (Vector Laboratories, Burlingame, CA) is added at 1:500 dilution to constructs to visualize vimentin binding. Simultaneous TPEF images of fibroblasts and SHG images of collagen fibers are both acquired for 3D directional variance analysis. TPEF images are acquired with an excitation wavelength of 760 nm and recorded by the 525 nm detector, and SHG images are excited using 800 nm and recorded by the 400 nm detector. There are three samples corresponding to each time point, and six 3D image stacks of either cell or collagen per group (2 per sample) are collected for analysis. 2.6.

Articular cartilage preparation

Knee joints from 18 week old Balb/C mice are harvested by severing the femur and tibia at mid bone length. The samples are equilibrated in a 30% sucrose solution for 3 days at 4 °C before being embedded in optimal cutting temperature (OCT) compound. Embedded joints are cryosectioned sagittally in a Leica CM 1950 cryostat at 40 μm thickness using Cryofilm type IIC (10) (University of Connecticut Medical Center, Rowe Laboratory). A total of 3 mice (6 knees) are harvested as samples. These samples are rehydrated for 15 min in PBS before imaging. Twelve simultaneous TPEF and SHG 3D image stacks (2 per sample) are collected for analysis. Using 800 nm as the excitation wavelength, endogenous TPEF fluorescence images are obtained with the 525 nm detector, and SHG images are detected using the 400 nm detector. All animal procedures were approved by the Tufts University Institutional Animal Care and Use Committee (IACUC). 2.7.

Mammary glands preparation

Normal mammary tissues are collected from the fourth mammary gland of two non-parous 12 week old FVB/N female mice, while mammary tumors are obtained from two mice of the human in mouse model [ 28 ]. These tumors are composed of cancer cells of human origin that have stromalized with mouse ECM and cancer associated fibroblasts as well as immune cells, which are transplanted into the glands of NOD/SCID mice. These fresh tissues are embedded in OCT compound and snap frozen in a methanol bath. They are immersed in phosphate buffered saline (PBS) to remove the OCT compound before use. All animal procedures were approved by the Tufts University IACUC. With 800 nm as the excitation wavelength, we acquire TPEF images of endogenous fluorescence using the 525 nm detector, and SHG images using the 400 nm detector. A total of 6 simultaneous TPEF and SHG 3D image stacks per group (3 stacks per tissue sample) are collected for analysis. 2.8.

Parietal peritoneum and primary pancreatic neoplastic tissue preparation

Freshly excised biopsies from healthy parietal peritoneum and primary pancreatic neoplastic tissue are acquired during abdominal surgery from 3 different patients that have a suspected or confirmed diagnosis of pancreatic malignancy and undergo open operative resection or biopsy of the malignancy as part of their treatment plan. These tissues are imaged within 3 hours after excision. Samples are excited with 900 nm to acquire depth-resolved TPEF and SHG 3D images. TPEF signal is detected using the 525 nm detector, whereas SHG signal is detected using the 460 nm detector. Although there is some crosstalk of fluorescence in the 460 nm channel, the SHG signal dominates in this region. Twelve simultaneous TPEF and SHG 3D images per group are acquired. Biopsy acquisition and imaging are conducted according to approved institutional review board protocols from Lahey Clinic and Tufts University. 2.9.

Engineered brain-like cortical tissue preparation

Salt leached silk scaffolds are prepared and used for 3D neuronal culture as previously established [ 29 ]. Briefly, silk scaffolds are made by combining 30 ml of 6% silk solution with 60 g of 500–600 μm NaCl particles for a period of 48 hours at room temperature to allow the crosslinking of silk solution. Then, the salt is leached out for 48 hours by placing the crosslinked silk in distilled water, which leaves behind a porous silk scaffold. The scaffold is then punched using disposable biopsy punches (6 mm diameter). The designed constructs are autoclaved and coated for 1 hour at 37 °C with 0.1 mg/ml poly-D-lysine (Sigma) solution to prepare for cell seeding. Primary rat neurons dissociated from embryonic day 18 (E18) cortices are seeded at a concentration of 1 million cells/scaffold and allowed to adhere to the pores of the scaffold overnight at 37 °C. Next, collagen type I (Corning) gel is prepared as 3 mg/ml solution on ice and added at 100 μl volume per scaffold, with the constructs placed in a 96-well plate. After 30 min of gelation at 37 °C, the constructs containing primary neurons and collagen gel ECM are flooded with Neurobasal media (supplemented by 2% B27, 1% GlutaMAX and 1% penicillin streptomycin) and cultured at 37 °C for two weeks. Controlled cortical impact (CCI) of the cell seeded constructs is performed at the two week time point by modifications to a custom built set up at Massachusetts General Hospital [ 30 ]. The injury is done using a 3 mm flat tip pneumatic piston that hits each construct at a velocity of 6 m/sec. The impact lasts for a duration of 100 ms, during which the piston penetrates a depth of 0.5 mm into the construct. The injured samples (n = 5) are fixed 10 min post-injury and an equal number of uninjured samples are fixed as well to serve as controls. The uninjured and injured samples fixed with 4% paraformaldehyde (Fischer Scientific) for 30 min at room temperature, are subsequently washed for 30 min with PBS three times. Immunostaining of the fixed samples is done using a previously described protocol [ 29 ]. Primary antibody against β-III tubulin (rabbit, Sigma), a neuron-specific marker is diluted by 1:500 in a blocking solution (4% goat serum, 0.2% triton-X100, 0.05% bovine serum albumin). The samples are incubated overnight at 4 °C in the primary antibody solution, followed by 30 min washes with PBS for three times. Next, the secondary antibody Alexa 488 goat-anti-rabbit (Invitrogen) is added to the samples at a dilution of 1:250 in blocking solution. The samples are left to incubate for 2 hours at room temperature and washed three times with PBS. A total of 20 TPEF 3D image stacks of neuron axons per group (4 stacks per sample) are acquired at an excitation wavelength of 760 nm and the 525 nm detector. 2.10.

Statistical analysis

For samples with multiple groups (articular cartilage and collagen hydrogels), a one-way ANOVA with post-hoc Tukey HSD test is used to assess significant differences using JMP 12. Otherwise a two-tailed t-test is used. Results are considered significant at p < 0.05.

Supplementary Material 1 2 3