2.1.

Cell Culture

BAEC were cultured on glass coverslips as previously described.¹⁴ When comparing the dark-field data with fluorescence imaging, cells were incubated in Hepes buffered balanced salt solution¹¹ supplemented with 100-nM MitoTracker Green (Invitrogen) for 45 min following the manufacturer’s protocol. The cover-slips containing the cells were mounted on a metal-slide for imaging as described in Ref. 14.

2.2.

Calcium Treatment

To induce subcellular morphological change, the cells were overloaded with calcium via treatment with ionomycin (Sigma-Aldrich), a calcium ionophore. Pretreatment imaging was conducted in a salt-solution¹¹ containing 1.5 mM ${\mathrm{CaCl}}_2$. The treatment solution was prepared by adding $20\text{-}\mu \mathrm{M}$ ionomycin to the same saline medium. To induce the calcium overload, the chamber containing the cells was flushed three times with this treatment solution. Taking this point as $t=0$, time-lapse post-treatment images were collected every 45 s for a maximum of 10 min. For control, a 1.5 mM ${\mathrm{CaCl}}_2$ with no ionomycin was flushed three times instead of the treatment solution.

2.3.

Optical Setup and Image Acquisition

Details of the optical setup are described in Ref. 13 and shown in Fig. 1. Briefly, an inverted microscope fitted with a $63\times$, $\mathrm{N.A.}=1.4$ oil immersion objective is used to image the sample illuminated with a laser at $\lambda =532\mathrm{nm}$. The incident laser light is passed through a spinning diffuser to average speckle. Scattered light is collected from the sample in transmission mode and filtered by displaying images of Gabor filters on a liquid crystal display (LCD) placed in a conjugate Fourier plane (${\mathrm{F}}_2$). To implement dark-field, the zeroth-order of diffraction is blocked by setting the corresponding pixels to 0 at the center of the LCD. The resulting filtered images are captured by a CCD camera placed in a conjugate imaging plane (${\mathrm{I}}_2$). The OSI image magnification on the CCD at ${\mathrm{I}}_2$ is $0.25\mu \mathrm{m}/\text{pixel}$. Through an additional port on the trinocular (not shown in the figure), the microscope also permits collection of epifluorescence images from the same sample on a different CCD camera.

Fig. 1

(a) OSI setup. The OSI technique involves the acquisition of transmission dark-field, filtered images on the charge-coupled device (CCD) camera. Two-dimensional (2-D) Gabor filters are displayed on a liquid crystal device (LCD) in the Fourier plane. ${\mathrm{L}}_1$ and ${\mathrm{L}}_2$: relay lenses, PBS: polarizing beam splitter, LP = linear polarizer, F1 and F2: conjugate Fourier planes, ${\mathrm{I}}_1$ and ${\mathrm{I}}_2$: conjugate imaging planes. The magnification on the CCD at ${\mathrm{I}}_2$ is $0.25\mu \mathrm{m}/\text{pixel}$. (b) Illustration depicting a two-dimensional Gabor filter displayed on the LCD and positioned at coordinates (U,V) in Fourier space. Pixels at the center of the LCD are turned off and act as a zeroth order block to implement dark-field.

2.4.

Gabor Filter Design

The Gabor filter images displayed on the LCD are created in MATLAB (The Mathworks, Natick, MA). In the object or image space, the filters are characterized by three parameters: period ($S$), orientation ($\phi$), and standard deviation (${\sigma }_{\text{space}}$). While deviation determines the size of the object area probed by the Gabor filter, the period and standard deviation are typically chosen not to be independent¹⁶ so that the region probed by the filter scales with the filter’s period. For this work, we choose ${\sigma }_{\text{space}}=\frac{S}{2}$, and we choose six periods, $S(\text{in}\mu \mathrm{m})=2.0$, 1.43, 1.11, 0.91, 0.77, 0.66, and four orientations $90\mathrm{deg}\le \phi <225\mathrm{deg}$ with 45-deg increment. To implement the filters in a conjugate Fourier plane, we convert these parameters to their frequency-domain equivalents within the actual optical setup. The period transforms to the spatial frequency, $F=\frac{1}{S}$ and ${\sigma }_{\text{frequency}}=\frac{1}{2\pi {\sigma }_{\text{space}}}$.¹⁷ The filters are created in MATLAB by defining a Gaussian function as shown below. Derivation of the equation is shown in Table 1 (see Appendix):

To calibrate for spatial frequency, a diffraction pattern with known spatial frequency spacing is used as described in Ref. 13. The calibration factor for the current setup is $0.0075\text{cycles}/\mu \mathrm{m}/\mathrm{LCD}\text{-}\text{pixel}$ and the distance, $F$, of the Gaussian-peaks from the center of the LCD can be calculated using this factor. We also correct for the LCD’s reflectance by correcting the Gaussian’s grayscale values using the LCD’s measured look-up-table, which gives the LCD reflectance as a function of gray scale input. As such, the grayscale values inputted into the LCD will result in the expected Gaussian-shaped reflectance. Since the Fourier transform of a Gabor function is a center-shifted Gaussian, the filters are realized by placing the peaks of the Gaussians at a distance $F$ from the center of the LCD. The coordinates of these pixels correspond to the spatial frequencies whose inverse corresponds to the prechosen periods.

Once the filtered-images are captured in the CCD, they are further processed, as shown in Fig. 2. To generate the orientedness images, local energy is calculated at each pixel first. Orientedness images are then masked using the dynamic masks generated using the maximum moment of phase-congruency (PC) images. Finally, features are extracted from both the binary masks and the masked-orientedness images. Calculation of local energy, orientedness, and PC is described in the following sections.

Fig. 2

Flowchart of the algorithm utilized for image analysis.

2.5.

Generating Local Energy and Phase Congruency (Maximum Moment) Images

In one dimension, PC is defined as $\mathrm{PC}(x)=\frac{|E(x)|}{{\sum }_n{A}_n(x)}$, where $A_n(x)$ are the amplitudes of $n$ Fourier components at a location ($x$), and local energy $E(x)$ is the vector summation of all components.¹⁵ PC acts as a highly localized operator to detect edges and corners in an image that is invariant to illumination and magnification. This technique takes advantage of the fact that Fourier components of the edge- or corner-points in an image are maximally in phase. This can be explained by a simple example of a square pulse. If all of the Fourier components of the square pulse are superimposed on each other, it can be shown that all of them will be in phase at the edges of the pulse. Thus, for the edge, the summation of the component-amplitudes is equal to the vector summation of all components, producing a PC value of 1. A minimum value of 0 may be obtained if all of the components are out of phase. Kovesi et al.¹⁵ proposed to extract the frequency information via wavelet transform instead of taking the Fourier transform. The wavelet transform is realized by even and odd symmetric log-Gabor filters of varying scale and orientation. Using the PC information, the edge-like features in the subcellular domain are then enhanced by taking the maximum moment of the PC covariance matrix. In this work, we have used an optimized code-version developed by Kovesi¹⁸ to generate the maximum moment and local energy images from each of the Gabor filtered images acquired in the setup. Thus, for each sample, 24 maximum moment images and 24 local energy images are calculated from the 24 optically filtered images acquired by the setup. These images are obtained from the “M” and “pcSum” variables found in the abovementioned code.

In addition, the user may set several input parameters pertaining to the calculation of PC. In particular, the local frequency information is obtained via digital log-Gabor wavelets of different scales and orientations. In this work, we set the number of scales (“nscale” in the code) to 4 and number of orientations (“norient”) to 4. We set the minimum wavelength of the wavelet (“minWaveLength”) to 2 pixels and the scaling between successive wavelengths (mult) to 2.1 Hence, for four scales, we obtain wavelets with wavelengths 2, 4.2, 8.82, and 18.5 pixels. To detect spatial properties at different orientations, four angles from 0 deg to 135 deg spaced 45 deg apart are used. The width of the filter function is controlled by the angular standard deviation ${\sigma }_G$ of the function that is dependent on the filter center frequency, $f_o$ as $R=\frac{{\sigma }_G}{f_o}$. We set $R=0.65$ for this work ($R$ corresponds to “sigmaOnf” in the code). Such combination of the scaling factor and $R$ ensures an even coverage of the spectrum. Finally, we set “$k$” = 1, “cutoff” = 0.5 and “gamma” (gain) = 5. As described in Ref. 15, the code utilizes a weighting function that penalizes any PC value that is not spread beyond a certain frequency. This is because PC is significant only when it occurs over a wide range of frequencies. This function is of sigmoid form and can be controlled by two parameters: cutoff fraction, which dictates the amount of frequencies considered for the PC to be significant, and gain (gamma), which controls the steepness of the sigmoid function. In addition, “$k$” corresponds to the number of standard deviations of the noise energy beyond the mean at which the noise threshold is set. Only signal energies beyond this threshold will be considered. For noisy images, the value of “$k$” can be set up to 20. However, this will increase the threshold and result in loss of useful data. In our case, we set the value as low as 1 to avoid any such loss while discarding a minimum level of noise.

2.6.

Generating Orientedness Images

Orientedness is a scattering-based parameter we had previously defined as the maximum signal over the average signal collected as a function of orientation, $\mathrm{o}\mathrm{r}\mathrm{i}\mathrm{e}\mathrm{ntedness}={\frac{\max (\text{Signal},{\phi }_i)}{⟨\text{Signal},{\phi }_i⟩}|}_{S=0.9\mu \mathrm{m}}$, and taken for data filtered with a Gabor period of $0.9\mu \mathrm{m}$.¹⁴ ${\phi }_i$ is the orientation of the Gabor filters used in the experimental setup. For a round object, there will be an equal amount of scattering in all $\phi$, resulting a minimum orientedness value of 1, whereas for an elongated object, the scatter profile will be oriented in a certain direction, producing an orientedness value greater than 1. As shown previously,¹⁴ orientedness decreases as long mitochondria fragment and round upon injury. In this work, we calculated orientedness using the local energy images, defined in the previous section, instead of the raw intensity values $\{\text{signal},{\phi }_i\}$. Here, $90\mathrm{deg}\le {\phi }_i<225\mathrm{deg}$, and the optically filtered images with a Gabor filter period of $2\mu \mathrm{m}$ were used in this calculation. Hence, the modified equation is given by $\text{orientedness}={\frac{\max (\text{local energy},{\phi }_i)}{⟨\text{local energy},{\phi }_i⟩}|}_{S=2\mu \mathrm{m}}$.

2.7.

Image Segmentation

First, the nucleus regions as well as the background outside the cells (no-cell region) are segmented out manually from the cell images so as to only analyze the scatterers within the cytoplasm. Then, following the main frame-work proposed by Jain and Farrokhnia,¹⁶ each Gabor-filtered image is preprocessed and fed to a clustering algorithm. However, we applied PC as a preprocessing step instead of the “blob-detection” approach applied in Ref. 16. Then, pixels in the maximum moment images obtained above were subjected to a K-means clustering algorithm. K-means clustering is an unsupervised learning approach that is used to assign data to a number of groups, $K$, based on the distances of data-points from $K$ initially assigned centroids. The iterative algorithm assigns each data-point to one of the $K$ groups until all of the points are assigned.¹⁹ After each iteration, each of the $K$ centroids is moved to the median intensity value of the pixels in its group. The algorithm is repeated until the sum of the point-to-centroid distances over all clusters is minimized. For this, we let the algorithm assign the initial centroids. We used the “cityblock” distance algorithm in MATLAB for the K-means distance calculation. We took the median value instead of the mean to avoid the effect of any outlier. Pixels in the 24 optically filtered images are grouped into four clusters here. We found 4 to be the highest number of clusters to be used for these images as the algorithm does not converge for a greater number of clusters. Pixels categorized in each cluster are color-labeled with one of the four colors: blue, turquoise, yellow, and red. In this analysis, all the experimental timepoints for a given cell were analyzed at once, yielding four classes evolving as a function of time.

3.1.

Cell Images Produced by the Algorithm

Images obtained at different steps are shown in Fig. 3. The Gabor-filtered image is obtained with $S=2.0\mu \mathrm{m}$, $\phi =225\mathrm{deg}$, and the corresponding PC (maximum moment), local-energy, and orientedness images are shown. The dark-field image is also shown. The dark-field image contains all spatial-frequency components; hence, it contains the scattering information from objects of all sizes and orientations. On the other hand, it is mostly the objects oriented at 225 deg that are highlighted in the Gabor-filtered image. Enhanced edges of these objects can be seen in the maximum moment image. The orientedness image maps the orientedness values of the subcellular objects in the cytoplasm. As indicated by the color scale, objects with a hue toward the blue have low orientedness and hence are round compared to the objects shown in red.

Fig. 3

(a) Dark-field, (b) Gabor-filtered ($S=2.0$, $\phi =225\mathrm{deg}$), (c) PC (maximum moment), (d) local energy, and (e) orientedness images of a representative cell. Note that edges are enhanced in the PC (maximum moment) image.

3.2.

Tracking Subcellular Organelles with Dynamic Masking

To segment the subcellular objects in the cytoplasm, we applied a K-means algorithm to classify the pixels in the maximum moment images. Before classification, the data are standardized by subtracting the mean and dividing by the standard deviation to obtain zero mean and unit variance. Figure 4 shows the four colored pixel clusters resulting from the K-means algorithm for a representative cell. To verify whether the clusters contain pixels with significant morphological information as well as to discard clusters that contain pixels from the background, we analyze the four centers, or centroids, of the four clusters [Fig. 5(a)]. Each centroid has 24 centroid components, representative of the 24 filtered intensities of pixels within the cluster associated with that centroid. The 24 centroid values are arranged as a function of optical filter period, $S$, and orientation $\phi$. As can be seen in the polar-plots [Fig. 5(a)], the centroid values of the turquoise (cluster 2) and red labels (cluster 4) vary as a function of period and orientation. On the other hand, there is no or little change in the centroid values for blue (cluster 1) and yellow labels (cluster 3). The blue and yellow clusters also have lower values. This indicates that the pixels from these two clusters are either from the background outside the cell or the nuclei, which were manually segmented out (blue in Fig. 4), or subcellular regions with significantly less intensity than the red and turquoise regions (yellow in Fig. 4). We quantified the variation in cluster values by taking the magnitude (mag) of each cluster, where $\mathrm{mag}=\sqrt{\sum _{S,\phi }{c}^2}$ and $c$ is the centroid-value corresponding to a filter $(S,\phi )$. The bigger the magnitude, the greater the variation over period and orientation. We then pick two clusters with the highest magnitudes and add them to create a single binary mask for the pretreatment condition [Fig. 5(b)]. Keeping the selected clusters the same, we repeat the combining process for each time point of data acquisition to generate masks that change over the period of time. We then create “masked orientednedness” images by multiplying the orientedness images with the binary masks at each time point (Fig. 6).

Fig. 4

Four clusters of pixels labeled in blue, turquoise, yellow, and red. Top panels: pretreatment, and bottom panels: post-treatment. Each of these clusters contains objects with different scattering profiles and intensities. Note also how the shape of the labeled objects changes with treatment; particularly elongated objects in the turquoise-cluster become smaller and rounder.

Fig. 5

(a) Color-coded centroid values are arranged in a $4\times 6$ polar plot corresponding to four orientations and six periods for each label shown in Fig 4. In both pre- and post-treatment cases, the centroid components of the turquoise and red clusters’ demonstrate variation over orientations ($\phi$) as well as periods ($S$), while those of the blue and yellow clusters have low signal with little or no variation. (b) The turquoise and red clusters are added together to produce a combined final binary mask.

Fig. 6

Masked orientedness image of a representative cell. Gradual change in the subcellular scatterers and shape of the dynamic mask can be observed upon injury at $t=0$.

As can be seen in Fig. 6, dynamic masks created by the proposed algorithm change shape over time as the cell shrinks upon severe chemical insult (Fig. 7). Most of the long objects in the masks become fragmented or round along with the cell-shrinkage and exhibit a decrease in their orientedness value. Concurrently, mitochondria become shorter and rounder upon calcium-overloading, observed in the MitoTracker-labeled fluorescent images (Fig. 8).

Fig. 7

Differential interference contrast (DIC) image of a cell under calcium-overload. (a) Long objects (yellow arrow) observed in the pretreatment image fragment and (b) become round due to excessive calcium. The cell also shrinks which is indicated by trace-marks of the cell membrane (pink arrow-head). When compared with MitoTracker-labeled fluorescent images, the fragmented and rounded long objects in the DIC images correspond to mitochondria (Fig. 8).

Fig. 8

MitoTracker labeled fluorescent images show mitochondrial fragmentation and remodeling upon calcium-overload.

3.3.

Effect of Jain’s Multichannel Filtering Approach

Jain’s framework is based on a multichannel filtering approach, originally proposed in Ref. 20, which dictates that the human visual system decomposes retinal information into a number of filtered images. The “multichannels” are realized by the Gabor filters. To understand the effect of Gabor filters on the final segmented image with classified pixels, we applied PC and clustering of the pixels on the unfiltered dark-field image. Segmented images obtained using the single-channel (DF only) and multichannel (Gabor) contain almost the same information in the cellular region, except the former has noise (Fig. 9). The multichannel approach decomposes the same information of the dark-field into multiple channels and allows selecting the channels less corrupted by noise. Hence, we can discard the channels manually (i.e., Gabor filters of particular period and orientation) that are corrupted by noise. This enabled us to generate a clearer image compared to that obtained by using dark-field only.

Fig. 9

(a) Segmented images generated using the multichannel Gabor-filtered images and (b) the single unfiltered dark-field image.

3.4.

Aspect-Ratio and Orientedness Decrease upon Injury

To quantify the subcellular morphological changes, we analyzed the subcellular objects segmented by the image masks obtained above (Fig. 5). Objects are defined as clusters of pixels with connectivity with four pixels.²¹ We measured the aspect-ratio of the objects over time to quantify the subcellular morphological changes [Fig. 10(a)]. We also calculated the median orientedness from the pixels contained in each object. Figure 10(b) shows that the object’s median orientedness decreases upon cell injury. Analysis of the aspect ratio of the MitoTracker labeled objects also shows a simultaneous decrease with treatment [Fig. 10(c)].

Fig. 10

Segmented subcellular objects become rounder after calcium overload, as indicated by the (a) aspect-ratio and (b) orientedness. The decrease in the aspect ratio of objects within the fluorescence images (c) suggests that the change in the aspect ratio of mitochondria partially accounts for the decrease in the masked subcellular objects’ aspect ratio. The data show mean and standard error for $n=6$ cells for the first six timepoints (treated) and the first seven timepoints (control). Additional timepoints were tested for one control cell and four of the treated cells.

3.5.

Texture-Features Can Detect Subcellular Morphological Dynamics Induced by Injury

Calcium-injury perturbs the spatial arrangement of the subcellular organelles. Such rearrangement is a result of changes in organelle morphology, as well as cell-shrinkage. While the treatment produces objects with a decreased aspect ratio and orientedness [Figs. 10(a) and 10(b)], relative distances among organelles also change and reduce as the cell shrinks. To quantify these two-dimensional, cell-wide, structural changes associated with the injury, we analyzed the texture of the masked orientedness images as a function of time. To this end, we extracted eight texture-features^22–24 using the code available²⁵ to analyze the spatial variation in the orientedness values. As proposed in Ref. 24, we grouped the features into three categories; features that describe the smoothness [contrast, dissimilarity, and inverse difference moment (IDM) or homogeneity], uniformity (entropy, maximum probability, and energy) and correlation (autocorrelation, correlation) of the texture.

We found that multiple texture features are sensitive to the morphological changes encoded by the orientedness images. We plotted these texture-features as a heat-map (Fig. 11) for each cell and observed the change in color over time, which indicates the features’ ability to detect the underlying structural changes. Several features changed by more than 30% in the treated cells. Autocorrelation, contrast, dissimilarity, and entropy increased in the treated cells, while they decreased slightly or remained within 10% of the starting value in the control cells. Energy and maximum probability decreased in the treated cells but remained within 10% of the starting value in the control cells. Cell-to-cell variations are also observed in both the control and treated cells. We also observed individual feature-responses for the whole dataset (Fig. 12). Autocorrelation, correlation, contrast, dissimilarity, and entropy increase while energy, maximum probability, and IDM decrease over time in the treated cells compared with the control. We select autocorrelation, contrast, entropy as the representative features from each group with positive trends (i.e., increase) and energy and IDM as negative trends (i.e., decrease). We then calculate the composite feature by subtracting the decreasing features from the increasing ones. The square root of the energy feature is taken to match the unit with that of other features. The composite feature in Fig. 13 also shows a shift from the baseline after injury.

Fig. 11

Changes in the texture-features of the masked orientedness images shown as heat-maps. Each row is a feature while each column is a time-point. Cells were monitored for different time periods, but for at least 180 s, after ionomycin treatment ($t=0$, treated), or the same saline solution but with no ionomycin ($t=0$, control). Texture-features are sensitive to the change in the spatial arrangement of the morphometric parameter orientedness. The color scale represents the magnitude of each parameter.

Fig. 12

Average responses of the individual texture features over time. The smoothness of the texture decreases as the contrast and dissimilarity increase and IDM (also known as homogeneity) decreases in the treated cells compared with controls. The uniformity of the texture decreases as the entropy increases, whereas the energy and max. probability decrease. Correlation between pixels increases. The data show mean and standard error for $n=6$ cells for the first six timepoints (treated) and the first seven timepoints (control). Additional timepoints were tested for one control cell and four of the treated cells. (a) Autocorrelation, (b) correlation, (c) contrast, (d) dissimilarity, (e) entropy, (f) energy, (g) maximum probability, and (h) IDM.

Fig. 13

Composite feature calculated from representative features: autocorrelation, contrast, entropy, energy, and IDM. The data show mean and standard error for $n=6$ cells for the first six timepoints (treated) and the first seven timepoints (control). Additional timepoints were tested for one control cell and four of the treated cells.