Rapid multispectral endoscopic imaging system for near real-time mapping of the mucosa blood supply in the lung

: We have developed a fast multispectral endoscopic imaging system that is capable of acquiring images in 18 optimized spectral bands spanning 400-760 nm by combining a customized light source with six triple-band filters and a standard color CCD camera. A method is developed to calibrate the spectral response of the CCD camera. Imaging speed of 15 spectral image cubes/second is achieved. A spectral analysis algorithm based on a linear matrix inversion approach is developed and implemented in a graphics processing unit (GPU) to map the mucosa blood supply in the lung in vivo . Clinical measurements on human lung patients are demonstrated. camera, but the disadvantage of limited imaging speed. The “snapshot” approach speeds up the image acquisition, but with reduced spatial resolution. In this paper, we report the development of a rapid multi-spectral endoscopic imaging system for near real-time quantitative mapping of the mucosa blood supply in the bronchial tree. Our method combines “sequential” scanning of triple-band illumination filters with standard color CCD camera based three-channel parallel image detection to speed up the imaging cube acquisition without sacrifice of the spatial resolution.


Introduction
Spectral measurement and analysis provide accurate quantification of tissue micro-vascular and morphological properties which are useful for monitoring and detecting early progression of various pathologies [1][2][3][4][5]. In particular, previous researches showed that quantifying the blood supply and morphology from reflectance spectroscopy measurements could provide useful parameters for cancer detection and treatment assessment in many organs such as pancreatic cancer [1], colon polyps [3], and lung cancer [4,5]. However, point spectroscopy does not account for tissue spatial heterogeneity, as it detects light from a single point area. This makes it ineffective for mapping a lesion.
On the other hand, spectral imaging, which involves the acquisition and analysis of a series of 2D images of reflected light sampled at different wavelengths -"spectral image cube," can provide spatial mapping of tissue morphology and physiology and account for tissue spatial heterogeneity. One potential application is developing endoscopic spectral imaging [6][7][8] which could provide critical information about the mucosa blood supply namely total hemoglobin and tissue oxygenation distribution. For example, quantifying mucosa blood supply and hypoxia related parameters from the spectral measurements could be valuable for prognosis prediction and treatment assessment of the bronchial tumors [9]. However, for practical clinical uses of spectral imaging, the acquisition and mapping of these tissue properties should be conducted rapidly and if possible in real-time to fight the movement of both tissue and the endoscope tip during the imaging procedure.
A number of systems and methods have been proposed for multispectral imaging. For example, Roblyer et al. [10] used "sequential" scanning of optical filters in the illumination and detection paths for improved detection of oral neoplasia. Kester et al. [11] built a "snapshot" multispectral camera with broadband illumination to demonstrate a 350 × 350 spatial resolution with 48 spectral bands at 15 image cubes per second. Other implementations include an acousto-optic tunable filter (AOTF), or Fabry-Perot cavities to acquire multispectral images [12,13]. The "sequential" scanning approach has the advantage of preserving the high spatial resolution provided by the camera, but the disadvantage of limited imaging speed. The "snapshot" approach speeds up the image acquisition, but with reduced spatial resolution. In this paper, we report the development of a rapid multi-spectral endoscopic imaging system for near real-time quantitative mapping of the mucosa blood supply in the bronchial tree. Our method combines "sequential" scanning of triple-band illumination filters with standard color CCD camera based three-channel parallel image detection to speed up the imaging cube acquisition without sacrifice of the spatial resolution.

Instrumentation
A spectral endoscopic imaging system (shown in Fig. 1) was developed to acquire in a noncontact fashion a series of reflectance images in the wavelength range of 400-760nm. The system is composed of a spectral light source module and a camera module connected with a flexible fiber optic endoscope. The spectral light module includes a 150W Xenon arc lamp (model #PE150AFM, PerkinElmer) for generating broad band light illumination. A filter wheel driven by a micro-stepper motor was mounted in front of the Xenon lamp to provide illumination of narrow spectral bands. Each filter mounted on the rotating wheel has dedicated special coating such that it can produce three narrow (~15 nm bandwidth) spectral bands (triple-band filter) simultaneously. A control unit in the light source module is used to synchronize the filter wheel rotation with the spectral image acquisition by the CCD camera. The filter wheel has a position encoder to communicate with the CCD camera. The customized CCD camera is based on a standard color CCD chip (Sony ICX618AQA) with pixel resolution of 659(H) × 494(V), pixel size of 5.6µm × 5.6µm, and imaging speed of 90 frames/second. As shown in Fig. 2, the three consecutive spectral bands in each filter are in the blue, green, and red/near infrared (NIR) wavelength ranges respectively. We used six such tripleband filters, producing 18 narrow spectral bands that cover the required spectral range. The camera module is made of a standard RGB color CCD camera to produce three simultaneous spectral images for each filter wheel rotation position. The eighteen spectral images (one image cube) were acquired at half video rate (~15 imaging cubes/s). The acquired spectral images are transmitted to the PC computer through USB port for data processing. The tissue reflectance spectral images were extracted from the acquired spectral images (as described in section 2.3). The video color image was generated from the integration of the reflectance spectral images convoluted with a CIE standard illuminant D65 [14] along the following three color bands (B: 415nm-490nm, G: 500nm-595nm, and R: 615nm-760nm). This is labeled as "Spectral Color" image in Fig. 4.

Spectral calibration method
The commercially available color CCD imaging detectors, including the one used in our system, typically have spectral responses for the R, G, B channels that overlap with each other due to the presence of Bayer filter mosaic. Figure 3 shows the typical spectral responses for the color CCD used in our system. Spectral calibration is thus required in order to correct for this effect before performing quantitative spectral imaging measurements and analysis using such CCD detector. The spectral calibration procedure is based on deducing calibration coefficients (matrix) that can compensate or calibrate the spectral overlapping response of the three detection channels per each triple-band filter. For example if we assume the light intensities that simultaneously reach the CCD are I b1 , I g1 , and I r1 as determined by filter F1 for the blue, green, and red narrow-band wavelengths respectively, the CCD blue channel reading can have contributions from the blue light as determined by the sensitivity T Bb1 , from the green light as determined by the sensitivity T Bg1 , and from the red light as determined by the sensitivity T Br1 . The CCD blue channel reading C B equals: Similarly the CCD green channel reading C G and red channel reading C R are given as follows: In matrix form we can express the above equations as: To spectrally calibrate the CCD response for this filter, the coefficients of the calibration matrix T RGB1 have to be obtained experimentally. This can be performed by illuminating the CCD with light I b1, I g1, and I r1 individually in separate exposures. For each single illumination, the light intensity was measured using a power meter and normalized during the analysis. The exposure time was fixed (at t = 2 seconds) during all the measurements. By inverting the above equation the spectrally calibrated light intensity as measured using this filter can be obtained from the CCD channel readings. The above calibration procedure is repeated for all the six triple-band filters to obtain T RGB1 , T RGB2 , ...., and T RGB6 in order to spectrally calibrate the CCD camera reading along the whole measured spectral range. Fig. 3. Typical spectral responses of the color CCD used in our system (obtained from Sony ICX618AQA data sheet).

Pre-processing algorithm
To account for the instrument response and to extract the tissue reflectance spectrum (R d ), the measured reflectance spectrum, R m (λ), is divided by the reference reflectance spectrum, R std (λ), obtained from a diffuse reflectance standard (Ocean Optics, FL, model #: WS-1). Thus for each pixel i in the image, the tissue reflectance R di (λ) is calculated according to [15]: Where, ki is a wavelength independent factor (referred here as the geometry correction factor) which depends on the angle of the incident light. This geometry correction factor depends on the position and orientation (distance and angle) of the endoscope tip from the tissue surface at the time of measurement. Since the position and the orientation of the endoscope tip during endoscopic navigation is not fixed and could be different between tissue measurements and standard disk measurements. This means that the geometry correction factor (ki) should be calculated in a dynamic manner at each measurement. The method we used to estimate the correction factor is based on the use of non-linear rigorous fitting of a light transport model to the measured spectra as described in details in our previous work [5], just that a different model as outlined in section 2.4 is used. However, applying this non-linear fitting inverse algorithm for each pixel in the image is very time consuming and would not be practical for real-time endoscopy imaging. For this consideration, we assumed that the effect of the endoscope tip orientation during measurement is independent of the pixel position within the image (i.e. ki = k for all pixels). The pre-processing algorithm is thus based on performing non-linear fitting using a light transport model (details in section 2.4) of the spectra we measured from one or two points only (say a point within the center of the image). At each point, the k value is based on the spatial averaging of 20 × 20 pixels at that point. The obtained value of k is then used to extract the tissue reflectance spectra R di (λ) for every pixel of the image by simple multiplication. The validity of this assumption and pre-processing method was tested in section 2.5.

Spectral analysis algorithm
The tissue reflectance spectra obtained from the above pre-processing step are then analyzed to quantify the hemoglobin volume fraction and hemoglobin oxygen saturation index. The method we used to rapidly quantify such parameters is based on the linear matrix inversion of the ratio between the absorption coefficient µ a (λ) and the reduced scattering coefficient µ' s (λ), which is given by γ(λ) = µ a (λ)/µ' s (λ), (referred here as relative absorption spectra). This relative absorption spectra γ(λ) can be obtained directly from the tissue reflectance spectra R d using a diffusion approximation model. For our analysis we have selected the δ-P1 approximation radiative transport model [16,17] in order to account for the spectral measurements obtained from the visible light illumination. Under this model assumption, the first two moments of the Henyey-Greenstein phase function (g and f respectively) are given by, f = g 2 and g* = g/(g + 1) 2 . In addition, for our broad beam illumination/point detection geometry, the 1-D δ-P1 approximation model [18] has been used. We have then reformulated the δ-P1 approximation model such that γ(λ) = µ a /µ' s can be extracted directly from the numerical inversion of the equation below:   7) and (8), and use matrix formulation we obtained the linear matrix inversion scheme below, which is used to quantify the mucosa's blood supply: Where, with (B' = B/V s ) being the relative blood volume fraction. The above equations were implemented on a graphics processing unit (GPU) (ZOTAC GeForce GTX 470) to satisfy the speed of analysis.

Experiments
We performed experimental spectral measurements and analyses to evaluate and validate the developed method and algorithms. The experimental work involved two main parts as described below. The first part of the experiment was performed to assess and evaluate the accuracy of the spectral analysis algorithm when used with the selected 18 spectral bands. For that respect we have collected about 300 in vivo spectra (400-800 nm, 10 nm resolution) using point spectrometer measurements. The device used for acquiring these in vivo point spectral data is described in our previous work [5,19] and has similar measurement geometry as the multispectral imaging device in this study. We then compared the spectral analysis results obtained from using the 18 selected bands extracted from the full spectra and the matrixinversion method versus that obtained from using the full spectra and non-linear Marquardt fitting method. The speed of analysis (computation time) was also calculated for the linear matrix-inversion method versus the non-linear fitting method with and without the using of the GPU.
The second part of the experiment and measurements were performed using our spectral imaging device. We first acquired spectral images from normal skin tissue of a healthy volunteer and from diffuse color standard (LabShpere Inc. model# BF07C) with the endoscope tip placed at different distances and angles from the targets. The objective of these measurements was to validate the calibration method we proposed for calculating the geometry correction factor (as in Eq. (5). For that respect we have calculated and compared the geometry factor obtained from using different points/pixels within the image view. Secondly, we performed preliminary in vivo endoscopic lung measurements on 10 patients to verify the system's safety and efficacy during endoscopy procedures. The power output from the endoscope tip is less than 80 mW for all the 6 illumination filters. This is much smaller than the optical power of standard white light bronchoscopy of 210 mW (Olympus CLX-F light source and BF-40D endoscope). We also calculated and quantified the mucosa blood properties and displayed the results as image-maps. The study was approved by the University of British Columbia -British Columbia Cancer Agency Research Ethics Board (Certificate #: H09-00364). Figure 4 shows 6 example images at 415 nm, 445 nm, 490 nm, 525 nm, 600 nm, and 690 nm from a spectral image cube from the lung mucosa obtained using our real-time multi-spectral imaging system. As referred above the system acquired 18 spectral images using six tripleband filters designed such that each triple-band filter has the first spectral band in the Blue spectrum range, the second spectral band in the Green spectrum range, and the third spectral band in the Red/NIR spectrum range. Figure 4 also shows the white light color image generated from the integration of the reflectance spectral images convoluted with a standard CIE illuminant (2nd image in Fig. 4, labeled as Spectral Color). For comparison, the image generated from a simple overlay of the narrow-band R, G, B images obtained with a single triple-band filter is shown as well (1st image in Fig. 4, labeled as RGB). The Spectral Color image provides better quality, objective color. Individual images at different wavebands show different intensity distributions from each other and also different from the white light color image. This reflects the wavelength dependence of tissue optical absorption and scattering properties. Thus multi-spectral imaging provides more information about the mucosa than the broad band RGB white light color imaging.  Figure 5 shows the correlation for the tissue parameters (blood oxygenation index or the tissue blood volume fraction) between the values determined by the linear matrix-inversion method using 18 spectral bands versus that determined by the non-linear least-square fitting method using continuous spectral measurements by point spectrometer. As shown in the figure, high correlations of 0.96 and 0.93 were obtained for the tissue blood oxygen saturation index and the tissue blood volume fraction respectively. The results demonstrated that the selected 18 bands analyzed using our new method and algorithms have comparable accuracy to that obtained using the non-linear least-square fitting of the continuous spectra. It should be noted that the choice of the delta-P1 approximation model was based on its improved accuracy, compared to diffusion approximation models, for predicting the light fluence generated from the propagation of the visible light within the tissue (which typically characterized by relatively high absorption coefficients). Nevertheless, the proposed method of spectral analysis is applicable for any light transport model. We have also compared the computation time from the two methods. Table 1 shows the computation speed obtained using the matrix-inversion scheme method versus that obtained using the non-linear least-square fitting method. As shown in the table, for analyzing one spectral measurement, the matrix-inversion method needs 0.6 ms using regular PC and 0.06 ms using GPU, while the least-square fitting method requires ~100 ms. This represents an increase of computation speed for the matrix inversion method of ~170-fold compared to that of the non-linear fitting method. For an image with 200x200 pixels/points, our method implemented on GPU can generate a quantitative map/image of the mucosa's blood supply in about 2.4 second. A more powerful GPU (such as NVIDIA GeForce GTX 980) can shorten the time to 450 ms which is near real-time.   Figure 6, summarizes the results obtained from the experiments used to assess and validate the accuracy of the method proposed for calibrating the spectra to account for the variation in the endoscope tip orientation during measurements. For that respect reflectance spectra were obtained from multiple points/pixels (1mm x 1mm square) along the color standard. The orientation of the endoscope tip is defined in terms of the angle (from the vertical line) and the vertical distance from the target, M (angle, distance). The figure shows the results obtained at three different orientations M1 (20°, 5 mm), M2 (20°, 2 mm), and M3 (40°, 5 mm). The mean and standard deviation of the geometry correction factors (k) obtained from these multiple measurements were calculated. This procedure was repeated for different endoscope tip orientations (as shown in Fig. 6(A)) and also for in-vivo skin target. The summary of the mean and standard deviation results obtained for different endoscopy orientations are shown in Fig. 6(B). The obtained results demonstrated that there was small variation in the calculated geometry factors (k) measured at different points/pixels within one image view. This suggests that the geometry factors are mainly dependent on the orientation of the endoscope tip during the measurement. Thus, using an arbitrary spatial point (pixel) within the image for calibration should provide acceptable sample of the geometry correction factor for all pixels within the image. This can also be understood from the observation that the correction factor is physically related to the distance and angle between the endoscope tip and the tissue and is independent of the tissue site or properties.

Results and discussions
We have also compared the measured spectra obtained from the same point on a red color standard but with the orientation of the endoscope tip changed for each measurement. The calculated geometry factor obtained for each different endoscope tip orientation was used to correct the measured spectra and to compensate for the variable measurement geometry between measurements. Example of the measured raw spectra is shown in Fig. 7(A) for the same point, but at 3 different orientations of the endoscope tip. These three spectra are quite different. After applying our calibration method, the corrected spectra are very close to each other as shown in Fig. 7(B).   Figure 8 shows tissue reflectance spectra measured at two locations within the endoscopic image view. The estimated values of the blood volume and tissue oxygenation index from each location are shown in Fig. 8(C). Both locations are normal looking mucosa. The results suggest that their blood volumes seem quite different, but the tissue oxygenation index seems quite close to each other. It should be noted that the absolute values of blood volume fraction obtained using our method could be affected by the changes in the scattering amplitude for different tissues. We have investigated a method for compensating this effect by quantifying the effective scattering volume fraction V s during the analysis and extracted the absolute blood volume fraction from the measurements. For that respect we have followed similar approach as suggested by Jacques et al [15] and assumed that the true water content of the mucosa tissue is fixed and known a priori (V w ). The scattering volume fraction V S is obtained in an iterative manner from updating the current value of V S by the multiplication factor (V w / V S ) until convergence occurred. Nevertheless, previous studies [20] have quantified the effect of varying the scattering coefficient (40 cm −1 -230 cm −1 ) on the accuracy of the blood volume absolute values and were found to be limited to +/−20%.  Figure 9 shows the blood volume and tissue oxygenation index maps for the endoscopic view showing in color image (9C). Again for this normal looking mucosa area, the blood volume seems varying quite a lot from location to location, while the tissue oxygenation index shows relatively less variations. This may be due to that the number of blood vessels could be different for different locations, while the oxygen supply in all vessels is constant within the image view with normal physiology.

Conclusion
We have presented in this paper a new multi-spectral endoscopic imaging system for near real-time mapping of mucosa's blood supply. The system can rapidly acquire spectral imaging cubes at half video rate (15 imaging cubes/s) using triple-band filters and color CCD camera. A rapid method based on linear matrix inversion scheme is used to analyze the reflectance measured from each pixel along the spectral images. With the advance of GPU processors, the proposed device and method would provide real-time mapping of mucosa's blood supply with accuracy comparable to that obtained from spectrometer point measurements. A clinical trial to test the utility of the system for improving early lung cancer detection is warranted.