Achieving cellular resolution for in vivo retinal images of transgenic GFAP-GFP mice via image processing

doi:10.1364/OE.16.008250

« Show journal navigation

Achieving cellular resolution for in vivo retinal images of transgenic GFAP-GFP mice via image processing

S. Kumar, G. Ho, K. M. Woo, and L. Zhuo »View Author Affiliations

Optics Express, Vol. 16, Issue 11, pp. 8250-8262 (2008)
http://dx.doi.org/10.1364/OE.16.008250

View Full Text Article

Acrobat PDF (923 KB)

Journal Search
Article Lookup

Article Tools

Citations

Alert me when this article is cited
Export Citation/Save

Article Navigation

▸ Article Outline
– Abstract
1. Introduction
2. Materials and methods
3. Results
4. Conclusions
– References and links

Article Tools
Full Text
Article Info
References
Cited By
Figures
Metrics
Download PDF

Viewing Equations in the
HTML Article

Optimal Viewing Recommendations

Full Text
Article Info
References (28)
Cited By
Figures (6)
Metrics

Abstract

In vivo retinal images of transgenic mice, expressing GFP under the control of the GFAP (glial fibrillary acidic protein) promoter, have very poor signal-to-noise ratio (SNR) and cellular resolution such that the analysis of GFAP–GFP expressing retinal cells from these images can be a very challenging task. We report an image averaging method based on a pixel rank matching criterion which significantly enhances both these image attributes. We also show that it compares favorably against direct image averaging and a commercial averaging routine available from the Heidelberg Retinal Angiograph 2 software.

1. Introduction

The advent of the green fluorescence protein (GFP) transgenic mouse model with the glial fibrillary acidic protein (GFAP) promoter has significantly advanced in vivo retinal imaging of GFP labeled astrocytes [1

L. Zhuo, B. Sun, C. L. Zhang, A. Fine, S. Y. Chiu, and A. Messing, “Live astrocytes visualized by green fluorescent protein in transgenic mice.” Dev. Biol. 187(1), 36–42 (1997). URL http://dx.doi.org/10.1006/dbio.1997.8601. [CrossRef] [PubMed]

]. It provides valuable insight to the dynamic changes in astrocyte morphology during development and in response to physiological and pathological conditions.

The images of these fluorescently labeled cells in the optic disc of the mouse are captured via a scanning laser ophthalmoscope (SLO) where the mouse retina is scanned with a laser beam from a point source and the fluorescence emitted is then collected by a photodetector. However, the in vivo imaging of the mouse retina does not achieve diffraction limited resolution due to the wavefront aberration induced by the eye.

Despite the higher system complexity and cost incurred, the use of adaptive optics in mouse retinal imaging has reduced higher order aberration [2

D. P. Biss, D. Sumorok, S. A. Burns, R. H. Webb, Y. Zhou, T. G. Bifano, D. Ct, I. Veilleux, P. Zamiri, and C. P. Lin, “In vivo fluorescent imaging of the mouse retina using adaptive optics.” Opt. Lett. 32(6), 659–661 (2007). URL http://dx.doi.org/10.1364/OL.32.000659. [CrossRef] [PubMed]

, 3

V. Nourrit, B. Vohnsen, and P. Artal, “Blind deconvolution for high-resolution confocal scanning laser ophthalmoscopy,” J. Opt. A 7, 585–592 (2005). [CrossRef]

]. However, the poor SNR level of these images, brought about by the thermal and shot noise from the detector and its amplifier, is another important factor which adversely affects image resolution. This is due to the application of image averaging on the sequence of successively acquired retinal images where the noise level may be reduced but at the expense of losing image resolution. The respiratory motion of the mouse results in significant changes to the morphology and location of the astrocytes between successive retinal images such that the averaged composite image appears diffused and occasionally contains undesirable image artifacts.

Various denoising methods are available [4

K. W. Leszczynski, S. Shalev, and N. S. Cosby, “An adaptive technique for digital noise suppression in on-line portal imaging.” Phys. Med. Biol. 35(3), 429–439 (1990). URL http://dx.doi.org/10.1088/0031-9155/35/3/011. [CrossRef] [PubMed]

, 5

Y. Yu and S. T. Acton, “Speckle reducing anisotropic diffusion,” IEEE Trans. Image Process. 11(11), 1260–1270 (2002). URL http://dx.doi.org/10.1109/TIP.2002.804276. [CrossRef]

, 6

L. Yin, R. Yang, M. Gabbouj, and Y. Neuvo, “Weighted median filters: a tutorial,” IEEE Trans. Circuits Syst. 43(3), 157–192 (1996). URL http://dx.doi.org/10.1109/82.486465. [CrossRef]

, 7

Y. Xu, J. B Weaver, D. M. Healy, and J. Lu, “Wavelet transform domain filters: a spatially selective noise filtration technique,” IEEE Trans. Image Process. 3(6), 747–758 (1994). URL http://dx.doi.org/10.1109/83.336245. [CrossRef] [PubMed]

] but, for lack of a better method and also due to its effective noise reducing capability [8

W. Swindell and M. A. Mosleh-Shirazi, “Noise reduction by frame averaging: a numerical simulation for portal imaging systems.” Med. Phys. 22(9), 1405–1411 (1995). URL http://dx.doi.org/10.1118/1.597618. [CrossRef] [PubMed]

], image averaging has become the de facto technique for filtering noise from image sequences. Several modifications to the image averaging method have however been suggested such as the alignment of images prior to averaging [9

D. C. Ercole, C. Giuseppe, M. Alessandro, M. Carlo, and V. Marco, “Compensation of random eye motion in television ophthalmoscopy: preliminary results,” IEEE Trans. Med. Img. 6(1), 74–81 (1987). URL http://dx.doi.org/10.1109/TMI.1987.4307800. [CrossRef]

, 10

D. U. Bartsch, M. H. El-Bradey, A. El-Musharaf, and W. R. Freeman, “Improved visualisation of choroidal neovascularisation by scanning laser ophthalmoscope using image averaging.” Br. J. Ophthalmol. 89(8), 1026–1030 (2005). URL http://dx.doi.org/10.1136/bjo.2004.057364. [CrossRef] [PubMed]

, 11

K. A. Goatman, A. Manivannan, J. H. Hipwell, P. F. Sharp, N. Lois, and J. V. Forrester, “Automatic registration and averaging of ophthalmic autofluorescence images,” in Medical Image Understanding and Analysis (2001).

, 12

B. M. Ege, T. Dahl, T. Søndergaard, and O. V. Larsen, “Automatic registration of ocular fundus images,” in 1st International Workshop on Computer Assisted Fundus Image Analysis (2000).

, 13

S. B. Stevenson and A. Roorda, “Correcting for miniature eye movements in high resolution scanning laser ophthalmoscopy,” Proc. SPIE 5688A, 145–151 (2005). URL http://dx.doi.org/10.1117/12.591190. [CrossRef]

, 14

J. B. Mulligan, “Recovery of motion parameters from distortions in scanned images,” in Proceedings of the NASA Image Registration Workshop (1997).

, 15

A. V. Cideciyan, “Registration of ocular fundus images: an algorithm using cross-correlation of triple invariant image descriptors,” IEEE Eng. Med. Biol. Mag. 14(1), 52–58 (1995). URL http://dx.doi.org/10.1109/51.340749. [CrossRef]

, 16

C. R. Vogel, D. W. Arathorn, A. Roorda, and A. Parker, “Retinal motion estimation in adaptive optics scanning laser ophthalmoscopy,” Opt. Express 14, 487–497 (2006). URL http://dx.doi.org/10.1364/OPEX.14.000487. [CrossRef] [PubMed]

, 17

A. Wade and F. Fitzke, “A fast, robust pattern recognition system for low light level image registration and its application to retinal imaging,” Opt. Express 3, 190–197 (1998). URL http://www.opticsexpress.org/abstract.cfm?URI=oe-3-5-190. [CrossRef] [PubMed]

] or alternatively the rejection of some images which are significantly different from that of a predefined reference image [18

I. Graesslin, M. Mittlebach, H. Eggers, T. Schaeffter, P. Bornert, and O. Lange, “MR real-time coronary imaging using the normalized cross-correlation,” in Proceedings of the International Society for Magnetic Resonance in Medicine (2002).

, 19

I. Graesslin, M. Mittelbach, T. Schaeffter, P. Bornert, H. Eggers, and O. Lange, “Adaptive weighted averaging for real-time MR imaging,” in Proceedings of the International Society for Magnetic Resonance in Medicine (2002).

, 20

C. J. Hardy, M. Saranathan, Y. Zhu, and R. D. Darrow, “Coronary angiography by real-time MRI with adaptive averaging.” Magn. Reson. Med. 44(6), 940–946 (2000). URL http://dx.doi.org/10.1002/1522-2594(200012)44:6¡940::AID-MRM16¿3.0.CO;2-F. [CrossRef] [PubMed]

, 21

M. S. Sussman, N. Robert, and G. A. Wright, “Adaptive averaging for improved SNR in real-time coronary artery MRI.” IEEE Trans. Med. Img. 23(8), 1034–1045 (2004). URL http://dx.doi.org/10.1109/TMI.2004.828677. [CrossRef]

, 22

M. S. Sussman, A. B. Kerr, J. M. Pauly, N. Merchant, and G. A. Wright, “Tracking the motion of the coronary arteries with the correlation coefficient,” in Proceedings of the International Society for Magnetic Resonance in Medicine (1999).

Image alignment requires the identification of landmark points such as the bifurcation of vessel branches and the establishment of point to point correspondence between landmark points in successive images via automated techniques [15

, 17

, 23

J. P. Lewis, “Fast normalized cross-correlation,” in Vision Interface (1995). URL http://citeseer.ist.psu.edu/lewis95fast.html.

]. However, unlike fluorescence angiography images, the high noise levels and the poor image resolution in our endogenous GFAP–GFP fluorescence images makes it a challenging task to even define these landmark points manually let alone in an automated manner. Even the alignment of these images is error prone since the astrocytes in our images exhibit a far greater degree of freedom compared to the vessels in the fluorescence angiography images which can be easily aligned using a rigid body transformation. More powerful alignment methods are therefore required but ironically this requires a larger number of landmark points which we cannot possibly identify from our images [12

B. M. Ege, T. Dahl, T. Søndergaard, and O. V. Larsen, “Automatic registration of ocular fundus images,” in 1st International Workshop on Computer Assisted Fundus Image Analysis (2000).

]. Conversely, alignment methods which do not rely on landmark points have not been convincing [19

]. Here, we propose a novel frame averaging method based on a pixel rank matching criterion which we will henceforth refer to as “averaging via rank matching” (ARM). The two main advantages of ARM are that it achieves significantly higher SNR and cellular resolution compared to conventional frame averaging methods. This is attributed to the noise robustness of the rank matching criterion and the fact that it leverages on the a-priori estimates of the underlying cellular profile obtained by deconvolving the individual retinal images [24

H.-U. Dodt, U. Leischner, A. Schierloh, N. Jhrling, C. P. Mauch, K. Deininger, J. M. Deussing, M. Eder, W. Zieglgnsberger, and K. Becker, “Ultramicroscopy: three-dimensional visualization of neuronal networks in the whole mouse brain.” Nat. Methods 4(4), 331–336 (2007). URL http://dx.doi.org/10.1038/nmeth1036. [CrossRef] [PubMed]

] from the sequence. Section 2.4.4 describes this in detail. This is a significant finding as it enables for the first time, the in vivo analysis on the role of astrocytes in the pathogenesis of various diseases such as diabetic retinopathy, Alzheimer’s disease and Parkinson’s disease. In addition, ARM can also be easily applied in the averaging of images acquired from other modalities such as MRI and confocal microscopy.

2. Materials and methods

2.1. Transgenic GFAP–GFP mice

The generation and genotyping of the transgenic GFAP–GFP mice are as previously described by Zhuo et al. [1

]. We used adult mice (8–10 weeks old) in the FVB/N background. Animal husbandry was provided by the Biological Resource Centre in Biopolis, Singapore whereas the experimental protocol covering the current study was approved by the Institutional Animal Care and Use Committee (IACUC).

2.2. Preparation of transgenic mice for retinal imaging

Mice were anaesthetized by intra-peritoneal (i.p.) injections with 0.15ml/10g body weight of Avertin (1.5% 2,2,2-tribromoethanol; T48402) purchased from Sigma-Aldrich (St. Louis, MO, USA), and their pupils dilated with a drop of 0.5% Cyclogyl® sterile ophthalmic solution (cyclopentolate hydrochloride, Alcon®, Puurs, Belgium). Custom-made PMMA hard contact lenses (from Cantor & Nissel, Northamptonshire, UK) were used to avoid dehydration of the cornea. Careful eye examination prior to scanning laser ophthalmoscope imaging ruled out the presence of any corneal or lens opacities.

2.3. Scanning laser ophthalmoscope (SLO) imaging

We employed a commercially available scanning laser ophthalmoscope, the Heidelberg Retina Angiograph 2, HRA 2 (Heidelberg Engineering, Dossenheim, Germany), for retinal imaging of the mice. The imaging system was adapted for the optics of the mouse eye by replacing the 30° focal lens with a 55° wide angle objective. This reduces the laser beam diameter to 1.7mm and in turn allows more light to be coupled into the small mouse pupil. Maximum retinal irradiance is approximately 2.0mWcm^-2 and, therefore, lies below the limits established by the American National Standards Institute. Since our primary interest is in the endogenous GFP signal in the mouse retina, we operate the HRA 2 in the fluorescence mode where the 488nm Argon laser provides the excitation light. The emission light is collected, with the use of the 500nm barrier filter, by a photodetector. The collected signal is then fed to a frame grabber interfaced to a computer at a frame rate of 5Hz where each image sequence comprises 45 single image frames. Each image is defined by 1536×1536 pixels with an intensity resolution of 8 bits/pixel. The estimated size of a focused spot on the retina i.e. resolution limit is approximately 10µm. By Nyquist sampling criterion, the sampling density must be at least 5.0µm/pixel i.e. twice the bandwidth frequency. Since we assume that 1° of visual angle subtended at the mouse eye corresponds to about 31µm of its retina [25

A. Maass, P. L. von Leithner, V. Luong, L. Guo, T. E. Salt, F. W. Fitzke, and M. F. Cordeiro, “Assessment of rat and mouse RGC apoptosis imaging in vivo with different scanning laser ophthalmoscopes.” Curr. Eye Res. 32(10), 851–861 (2007). URL http://dx.doi.org/10.1080/02713680701585872. [CrossRef] [PubMed]

] and our field of view is about 55.9°, the actual sampling density works out to be approximately

\frac{31 \times 55.9}{1536} μ m / pixel = 1.13 μ m / pixel

2.4. Proposed frame averaging algorithm

Here, we aim to obtain a composite fluorescence image of the mouse retina, with improved SNR and cellular resolution, from a sequence of relatively noisy and poorly resolved retinal fluorescence images. This is achieved by first aligning the images such that the between-frame motion artifacts are reduced and then averaging the aligned images. To this end, the proposed ARM algorithm comprises four major sections: (i) Detecting landmark points, (ii) Establishing point-to-point correspondence, (iii) Aligning images via affine transformation and (iv) Aligning images via the novel rank matching technique.

2.4.1. Detecting landmark points

The high noise levels in the image sequence prevent us from accurately detecting the landmark points. The image sequence can be formally defined as I_i , i=1,2, …,N where the subscript i identifies a particular image in the sequence of N images and N=45. In ARM, we propose an automated method of selecting landmark points from the a-priori estimates J_i , i=1,2, …,N of the underlying fluorescence signal in each retinal image from the original sequence. The a-priori estimates are obtained by deconvolving each image I_i (Huygens Essential Software, Scientific Volume Imaging BV) [24

] in the sequence. The deconvolution software computes a theoretical point spread function (PSF), resembling a 2-D Gaussian profile in the x-y plane, based on the ophthalmoscope device settings. Also, it automatically estimates the image SNR and uses this as a regularization parameter to control the sharpness of the deconvolved image. In doing so, photon noise is effectively removed. Next, we apply a gray-scale morphological operation [26

P. Soille, Morphological Image Analysis (Springer-Verlag, 1999). URL http://web.ukonline.co.uk/soille/book1st.

] to locate the intensity maximas (peak intensities) in each J_i estimate. Naturally, the locations of the intensity maximas in the I_i images will be identical to those in the corresponding J_i images. These intensity maximas serve as landmark points where their displacement over successive images in the sequence can be used to align these images. Their localized “spot-like” neighborhoods in J_i makes them effective landmark points since displacements along any direction, in the x-y plane, can be accurately detected.

2.4.2. Establishing point-to-point correspondence between adjacent images

Classical frame averaging methods establish correspondence between the landmark points in the images I_i , from a sequence, to those in the reference image I ₁ from the same sequence. However, establishing this correspondence can be impossible in some cases since the fluorescence signal in image I_i of the sequence may have changed significantly in terms of intensity and location from the corresponding signal in the reference image I ₁. ARM establishes point-to-point correspondence between successive images I _i-1 and I_i in the sequence such that each pair of successive images in the sequence has its own set of corresponding landmark points. To this end, the normalized cross-correlation measure (NCC) [23

J. P. Lewis, “Fast normalized cross-correlation,” in Vision Interface (1995). URL http://citeseer.ist.psu.edu/lewis95fast.html.

]

{NCC}_{J_{i, u, v}, J_{i - 1, u', v'}} = \frac{〈 J_{i, u, v} - {\bar{J}}_{i, u, v}, J_{i - 1, u', v'} - {\bar{J}}_{i - 1, u', v'} 〉}{∣ J_{i, u, v} - {\bar{J}}_{i, u, v} ∣ ∣ J_{i - 1, u', v'} - {\bar{J}}_{i - 1, u', v'} ∣}

(1)

is applied as shown in Eq. 1 where J_i,u,v is the (2W_c +1)×(2W_c +1) neighborhood of the i^th deconvolved image centered around one of its landmark points (u,v) and W_c is the width of the neighborhood, J̄_i,u,v is the mean value of this neighborhood, J _{i-1,u′,v′} is the (2W_c +1)×(2W_c +1) neighborhood of the (i-1) ^th deconvolved image centered around one of its landmark points (u′,v′) and J̄ _{i-1,u′,v′} is the mean value of this neighborhood. The notation 〈·〉 denotes the dot operator whereas |·| denotes the magnitude operator. We again use the a-priori estimates J_i since they are cleaner and better resolved compared to the I_i images.

Point-to-point correspondence between landmark points (u,v) and (u′,v′) is established if (u,v) lies within the (2W_c +1)×(2W_c +1) neighborhood of (u′,v′) and their NCC measure,

{NCC}_{J_{i, u, v}, J_{i - 1, u', v'}}

, is the largest exceeding an empirically predefined threshold, NCC _th . NCC _th may vary from one sequence to another due to various factors such as SNR levels but typically its value ranges from 0.7–0.9. NCC _th should be large enough so that spurious landmark points are avoided but at the same time it should be low enough such that valid landmark points are not excluded. W_c is empirically set to 11 since it was observed that the maximum displacement in the field of view remains within this vicinity for all images in the sequence. The point-to-point correspondence in the I_i images follows exactly that determined from the J_i images.

2.4.3. Linear alignment of images via affine transformation

ARM implements the linear alignment of images in a sequential manner and takes into consideration the translation, rotation, scaling and shearing effects between successive images in the sequence. We begin with the alignment of I ₂ to I ₁ resulting in the aligned image I^a ₂ followed by the alignment of I ₃ to I ₁ and the sequential alignment process continues in this way until the last image I_n in the sequence is aligned. The alignment process at every step can be formally expressed as follows where

I_{i}^{a} = f_{i} f_{i - 1} . . . f_{2} (I_{i})

(2)

the superscript a denotes that I_i has been aligned to I ₁ and f ₂ represents the affine transformation function determined from the point-to-point correspondence of the landmark points between I ₁ and I ₂. Similarly, f_i is determined from the point-to-point correspondence of the landmark points between adjacent images I _i-1 and I_i . The motivation here is that transformation functions f_i , i=1,2, …,N between adjacent images can be more accurately determined compared to those between any given image and a fixed reference image. The reason is that the motion artifacts between adjacent images tend to be lesser as compared to images that are further apart.

But the retina is a biological tissue and as such affine transformation alone may not accurately compensate the complex respiratory motion between successive images in the sequence. In the next section, we propose a non-linear alignment approach which overcomes this limitation.

2.4.4. Non-linear alignment of images via rank matching

In the heart of ARM lies a pixel rank matching criterion which defines our non-linear alignment process. The rank of a pixel denotes the order of that pixel in an image neighborhood assuming all pixels in that neighborhood have been sorted in ascending order based on their intensity values. If image A is to be aligned to image B, the criterion requires the pixel values in image A to be updated such that the rank of the center pixel in every local neighborhood of image A is identical to that of the center pixel from a corresponding neighborhood of image B. This is a reasonable requirement since the rank of both pixels should not differ despite there being variations in the average intensity level in their neighborhood. However, this assumption is only valid in the absence of noise which could significantly alter the pixel ranks.

In order to ensure the noise robustness of our rank matching criterion, we utilize the rank information obtained from the aligned a-priori estimates J^a _i , i=1,2, …,45 computed by subjecting J_i to the same transformation functions as in Eq. 2. As in Section 2.4.3, the non-linear alignment is implemented sequentially as shown in Eq. 3 where

I_{i}^{a, r} (x, y) = I_{i}^{a} (x_{o}, y_{o}) \Leftarrow RANK I_{i}^{a} (x_{o}, y_{o}) = RANK J_{i - 1}^{ref} (x, y)

(3)

where the RANK operation is performed on a local (2W_r +1)×(2W_r +1) neighborhood centered at (x,y) and W_r is the width. The pixel (x_o ,y_o ) is located within this neighborhood and the superscript r in I^a,r _i denotes that I^a _i has been non-linearly aligned such that the ranks of its pixels match those of J^ref _i-1 and i≥2. J^ref _i is defined as

J_{i}^{ref} = \frac{(i - 1) J_{i - 1}^{ref} + J_{i}^{a, r}}{i}

(4)

where J^ref ₁=J ₁ and J^a,r _i is computed as follows:

J_{i}^{a, r} (x, y) = J_{i}^{a} (x_{o}, y_{o}) \Leftarrow RANK J_{i}^{a} (x_{o}, y_{o}) = RANK J_{i - 1}^{ref} (x, y)

(5)

Finally, the result of averaging i number of aligned images in the sequence is given by

I_{i}^{ref} = \frac{(i - 1) I_{i - 1}^{ref} + I_{i}^{a, r}}{i}

(6)

where I^ref ₁=I ₁. Eqs. 3–6 are implemented in a cyclic fashion where for every cycle, i is incremented before proceeding to the next cycle where Eqs. 3–6 are recomputed and I^ref _N denotes the final averaged image where N is the total number of images in the sequence.

3. Results

3.1. Estimating the underlying fluorescence signal

Figure 1(a) is a typical image of the optic disc in the retina of a transgenic mouse where the fluorescence signal is from the GFP protein which, in this case, labels the glial cells. As observed, the image has a very low SNR and this makes it extremely difficult to detect any suitable landmark points and it is clearly inappropriate to use such images as reference images in rank matching. An a-priori estimate is therefore required to accurately estimate the underlying fluorescence signal. We have explored several methods of determining this estimate and found that the Huygens deconvolution method (Huygens Essential Software, Scientific Volume Imaging BV) [27] gives the best results. The deconvolved result of Fig. 1(a) is shown in Fig. 1(b) and as observed, the speckle noise is removed and the fluorescence intensity profile is characterized by various intensity maximas. Conversely, a 9×9 median filtering method of estimating the a-priori estimate fails to remove the speckle noise and instead distorts the image as observed in Fig. 1(c). Figure 1(d) shows the results of applying anisotropic diffusion [5

Y. Yu and S. T. Acton, “Speckle reducing anisotropic diffusion,” IEEE Trans. Image Process. 11(11), 1260–1270 (2002). URL http://dx.doi.org/10.1109/TIP.2002.804276. [CrossRef]

] where the noise appears to be removed but at the expense of diffusing the underlying fluorescence signal resulting in the image appearing blurred overall.

Fig. 1. (a) Input image, I_i , of the GFP signal from the optic disc of the transgenic mouse (b)–(d) J_i obtained after subjecting I_i to (b) Huygens deconvolution (c) 9×9 median filtering and (d) speckle reducing anisotropic diffusion.

3.1.1. Tuning of the deconvolution parameters

The Huygens deconvolution method requires certain device information such as the pinhole diameter, numerical aperture, excitation wavelength etc. but assumes that the device in question is a wide-field or confocal microscope and not a scanning laser ophthalmoscope. We circumvent this problem by defining parameters for a hypothetical confocal microscope system such that the ratio of its Nyquist critical sampling density over the actual sampling density is comparable to that of the HRA 2 system specified in Section 2.3 (5.00µm/pixel and 1.13µm/pixel respectively). To this end, we define the confocal microscope parameters as shown in Table 1 so that the Nyquist critical and actual sampling density are 70nm and 15.8nm respectively. This ensures that both systems have the same degree of oversampling i.e.

\frac{5.00}{1.13} = \frac{70.0}{15.8}

. The theoretical point spread function for deconvolution was subsequently computed based on the parameters specified in Table 1.

Table 1. Parameters of imaginary confocal microscope for the Huygens deconvolution method

Property	Description
Microscope type	Confocal
Numerical aperture	0.9
Lens refractive index	1.515
Medium refractive index	1.515
Photon count	1
Emission wavelength, E_m (nm)	500≤E_m ≤700
Excitation wavelength, E_x (nm)	488
Backprojected pinhole radius (nm)	340
Pinhole spacing (µm)	3.4
Lateral resolution (nm)	15.8

3.2. Fast and accurate tracking of landmark points

We apply the NCC measure in Section 2.4.1 to automatically select intensity maximas as landmark points for every a-priori deconvolved image J_i in a sequence. Although a large number of landmark points may have been detected in an image J_i , only a subset of these points correspond to those in the adjacent image J _i-1 since we require that

{NCC}_{J_{i, u, v}, J_{i - 1, u', v'}} \geq 0.9

and that the corresponding points lie within a 23×23 square neighborhood of each other (W_c =11 as in Section 2.4.2).

Fig. 2. Point-to-point correspondence between: First row-(a) J ₁ and (b) J ₂, Second row-(c) J ₂ and (d) J ₃, Corresponding points in the adjacent images carry the same number.

Figure 2 illustrates this where Figs. 2(a) and (b) are a-priori estimates of the first and second images in a sequence whereas Figs. 2(c) and (d) are a-priori estimates of the second and third images in the same sequence. As observed, each pair has their own unique set of corresponding landmark points indicated by the labeled green and red dots.

3.3. Quantitative analysis of ARM

Here, we quantify the influence of the window size W_r and image SNR on the alignment accuracy of the rank-matching criterion. This is achieved by omitting the affine transformation step of Section 2.4.3 and directly applying the rank-matching criterion of Section 2.4.4 on a sequence of 45 noisy synthetic images where each 256×256 image comprises a 2-D “ball” shaped intensity profile with a diameter of approximately 36 pixels. The “ball” is vertically displaced by two pixels between adjacent images. The Huygens deconvolution method was applied to obtain the a-priori deconvolved images J_i from the sequence.

3.3.1. Influence of W_r on alignment accuracy

In each case, the final averaged image I^ref ₄₅ is compared against the ideal (noiseless) synthetic image I′ and the mean absolute error (MAE) is computed as

MAE = \frac{1}{M} \sum_{i = 1}^{M} ∣ I_{45, i}^{ref} - I_{i}^{'} ∣

and the results shown in Fig. 3(g) where i and M denote the pixel entry and total number of image pixels respectively. It is evident here that the distortion is smaller for larger window sizes and that it is minimal for W_r =19. Fig. 3(a) shows a sample image from a noisy sequence with an SNR of 0dB whereas the final averaged results in Figs. 3(b)–(e) illustrates the influence of W_r on the alignment accuracy. We observe that the averaged results in Figs. 3(b)–(d) appear truncated. This is due to the insufficient window size of (2W_r +1)×(2W_r +1) which limits the validity of the rank matching operation to the peak intensity region of the “ball” profile as observed in Figs. 3(b)–(d) although the correctly aligned region increases for larger window sizes. Accurate alignment of the entire “ball” profile can only be obtained if the window is large enough to cover the entire “ball” and its displaced counterpart in the adjacent image i.e. 2W_r +1 ≃ Ball diameter+Between-frame pixel displacement. The window size requirement is met for Fig. 3(e) and as observed, the distortion is minimal. In contrast, the 2-D “ball’ profile appears totally distorted in the direct frame averaged result of Fig. 3(f) due to the vertical displacements between adjacent images in the sequence.

3.3.2. Influence of SNR on alignment accuracy

In Figs. 4(d)–(f), we observe the efficacy of the ARM method in reducing the noise levels while minimizing distortion to image sequences with different SNR levels as shown in Figs. 4(a)–(c). A window size of W_r =19 is used in all three cases. Significantly small error values of about 1.9% and 0.5% were achieved for SNR levels of -10dB and 0dB respectively. Even for a poor SNR level of -20dB, a reasonably low MAE of about 6.5% was achieved. The efficacy of ARM is attributed to the noise robustness of the deconvolution method and the proposed rank matching criterion.

3.4. Qualitative analysis of ARM on noisy optic disc images of transgenic mice

We have also evaluated the performance of ARM on noisy optic disc images of transgenic GFAP–GFP mice and found that the results corroborate the findings made in Section 3.3 on the synthetic images.

Fig. 3. Influence of W_r on the alignment accuracy: (a) Reference image, (b)–(f) Final averaged images corresponding to (b) W_r =4, (c) W_r =9, (d) W_r =14, (e) W_r =19, (f) Direct averaging and (g) Corresponding drop in mean absolute error (MAE) error values from (b) to (e).

3.4.1. Influence of Wr on alignment accuracy

As W_r increases, the resultant improvement in the alignment accuracy can be observed from the final averaged images I^ref _N of Fig. 5 where N=45. Image artifacts are visible in I^ref _N of Fig. 5(a) due to the small window size of W_r =2. Although the underlying fluorescence profile is taking shape in Fig. 5(b), it appears that only the maxima and minima intensity regions of the underlying fluorescence have been aligned accurately resulting in the appearance of tiny fluorescence dots and narrow branches of blood vessels in the optic disc. This is consistent with the observation made in Section 3.3.1. Accurate results can only be achieved if W_r ≃17 since the sum of the “spot” diameter in the J_i images and the displacement between adjacent retinal images is approximately 35 pixels (2W_r +1=35). This window size corresponds to the result shown in Fig. 5(d). A slight deterioration in signal resolution is observed when the larger window size of W_r =22 is applied in Fig. 5(e) whereas minor image artifacts are observable for the smaller window size of W_r =12.

Fig. 4. Influence of W_r on the alignment accuracy: (a) Reference image, (b)–(f) Final averaged images corresponding to (b) W_r =4, (c) W_r =9, (d) W_r =14, (e) W_r =19, (f) Direct averaging and (g) Corresponding drop in mean absolute error (MAE) error values from (b) to (e).

Fig. 5. Final averaged image I^ref _N of a transgenic mouse retina from the ARM method (a) W_r =2, (b) W_r =7, (c) W_r =12 (d) W_r =17 and (e) W_r =22.

3.4.2. Comparison against conventional averaging methods

We also find that the final averaged image I^ref _N of Fig. 6(a), processed using ARM, has a significantly better resolution compared to those obtained using direct frame averaging in Fig. 6(b) as well as the averaging module in the commercial Heidelberg Retinal Angiograph 2 (HRA 2) software in Fig. 6(c). In addition, the noise level in Fig. 6(a) appears to be notably lower than those in Figs. 6(b) and (c). The poor resolution of the results in Fig. 6(b) and (c) makes it extremely difficult if not impossible to perform morphometric and intensity analysis of the fluorescence spots in the optic disc region.

4. Conclusions

We have presented a fully automated frame averaging method of in vivo GFAP–GFP fluorescent images of the transgenic mouse retina. It is based on a non linear alignment method we call averaging via rank matching (ARM) which ensures a significantly higher resolution and SNR of the final averaged image compared to conventional frame averaging and motion compensated frame averaging methods. Rank matching provides a significant improvement to the alignment accuracy of standard motion compensation operations such as affine transformation as it is robust against noise and variations in the fluorescence intensity of astrocytes between successive images. The method also enables an automated, accurate and robust detection of landmark points from the peak intensity points of the fluorescence profiles. These landmark points are subsequently used to establish point-to-point correspondence between successive images rather than between a fixed reference frame and all other images in a sequence since successive images are more closely aligned. Other potential applications of ARM include the averaging of images acquired from other modalities such as confocal microscopy, MRI and PET.

Fig. 6. Qualitative comparison of the composite GFAP–GFP transgenic mouse retina images, from the same sequence of noisy images, processed using (a) the proposed ARM method, (b) Direct frame averaging and (c) HRA 2’s averaging module.

Acknowledgements

We would like to thank Dr. Xiao-Shan Min for the confirmation of a proper pupil dilation and lack of the corneal or lens opacities in the anesthetized mice prior to retinal imaging. The research was supported and funded by the Institute of Bioengineering and Nanotechnology (IBN), the Agency of Science, Technology and Research (A*STAR), the Republic of Singapore.

References and links

1.	L. Zhuo, B. Sun, C. L. Zhang, A. Fine, S. Y. Chiu, and A. Messing, “Live astrocytes visualized by green fluorescent protein in transgenic mice.” Dev. Biol. 187(1), 36–42 (1997). URL http://dx.doi.org/10.1006/dbio.1997.8601. [CrossRef] [PubMed]
2.	D. P. Biss, D. Sumorok, S. A. Burns, R. H. Webb, Y. Zhou, T. G. Bifano, D. Ct, I. Veilleux, P. Zamiri, and C. P. Lin, “In vivo fluorescent imaging of the mouse retina using adaptive optics.” Opt. Lett. 32(6), 659–661 (2007). URL http://dx.doi.org/10.1364/OL.32.000659. [CrossRef] [PubMed]
3.	V. Nourrit, B. Vohnsen, and P. Artal, “Blind deconvolution for high-resolution confocal scanning laser ophthalmoscopy,” J. Opt. A 7, 585–592 (2005). [CrossRef]
4.	K. W. Leszczynski, S. Shalev, and N. S. Cosby, “An adaptive technique for digital noise suppression in on-line portal imaging.” Phys. Med. Biol. 35(3), 429–439 (1990). URL http://dx.doi.org/10.1088/0031-9155/35/3/011. [CrossRef] [PubMed]
5.	Y. Yu and S. T. Acton, “Speckle reducing anisotropic diffusion,” IEEE Trans. Image Process. 11(11), 1260–1270 (2002). URL http://dx.doi.org/10.1109/TIP.2002.804276. [CrossRef]
6.	L. Yin, R. Yang, M. Gabbouj, and Y. Neuvo, “Weighted median filters: a tutorial,” IEEE Trans. Circuits Syst. 43(3), 157–192 (1996). URL http://dx.doi.org/10.1109/82.486465. [CrossRef]
7.	Y. Xu, J. B Weaver, D. M. Healy, and J. Lu, “Wavelet transform domain filters: a spatially selective noise filtration technique,” IEEE Trans. Image Process. 3(6), 747–758 (1994). URL http://dx.doi.org/10.1109/83.336245. [CrossRef] [PubMed]
8.	W. Swindell and M. A. Mosleh-Shirazi, “Noise reduction by frame averaging: a numerical simulation for portal imaging systems.” Med. Phys. 22(9), 1405–1411 (1995). URL http://dx.doi.org/10.1118/1.597618. [CrossRef] [PubMed]
9.	D. C. Ercole, C. Giuseppe, M. Alessandro, M. Carlo, and V. Marco, “Compensation of random eye motion in television ophthalmoscopy: preliminary results,” IEEE Trans. Med. Img. 6(1), 74–81 (1987). URL http://dx.doi.org/10.1109/TMI.1987.4307800. [CrossRef]
10.	D. U. Bartsch, M. H. El-Bradey, A. El-Musharaf, and W. R. Freeman, “Improved visualisation of choroidal neovascularisation by scanning laser ophthalmoscope using image averaging.” Br. J. Ophthalmol. 89(8), 1026–1030 (2005). URL http://dx.doi.org/10.1136/bjo.2004.057364. [CrossRef] [PubMed]
11.	K. A. Goatman, A. Manivannan, J. H. Hipwell, P. F. Sharp, N. Lois, and J. V. Forrester, “Automatic registration and averaging of ophthalmic autofluorescence images,” in Medical Image Understanding and Analysis (2001).
12.	B. M. Ege, T. Dahl, T. Søndergaard, and O. V. Larsen, “Automatic registration of ocular fundus images,” in 1st International Workshop on Computer Assisted Fundus Image Analysis (2000).
13.	S. B. Stevenson and A. Roorda, “Correcting for miniature eye movements in high resolution scanning laser ophthalmoscopy,” Proc. SPIE 5688A, 145–151 (2005). URL http://dx.doi.org/10.1117/12.591190. [CrossRef]
14.	J. B. Mulligan, “Recovery of motion parameters from distortions in scanned images,” in Proceedings of the NASA Image Registration Workshop (1997).
15.	A. V. Cideciyan, “Registration of ocular fundus images: an algorithm using cross-correlation of triple invariant image descriptors,” IEEE Eng. Med. Biol. Mag. 14(1), 52–58 (1995). URL http://dx.doi.org/10.1109/51.340749. [CrossRef]
16.	C. R. Vogel, D. W. Arathorn, A. Roorda, and A. Parker, “Retinal motion estimation in adaptive optics scanning laser ophthalmoscopy,” Opt. Express 14, 487–497 (2006). URL http://dx.doi.org/10.1364/OPEX.14.000487. [CrossRef] [PubMed]
17.	A. Wade and F. Fitzke, “A fast, robust pattern recognition system for low light level image registration and its application to retinal imaging,” Opt. Express 3, 190–197 (1998). URL http://www.opticsexpress.org/abstract.cfm?URI=oe-3-5-190. [CrossRef] [PubMed]
18.	I. Graesslin, M. Mittlebach, H. Eggers, T. Schaeffter, P. Bornert, and O. Lange, “MR real-time coronary imaging using the normalized cross-correlation,” in Proceedings of the International Society for Magnetic Resonance in Medicine (2002).
19.	I. Graesslin, M. Mittelbach, T. Schaeffter, P. Bornert, H. Eggers, and O. Lange, “Adaptive weighted averaging for real-time MR imaging,” in Proceedings of the International Society for Magnetic Resonance in Medicine (2002).
20.	C. J. Hardy, M. Saranathan, Y. Zhu, and R. D. Darrow, “Coronary angiography by real-time MRI with adaptive averaging.” Magn. Reson. Med. 44(6), 940–946 (2000). URL http://dx.doi.org/10.1002/1522-2594(200012)44:6¡940::AID-MRM16¿3.0.CO;2-F. [CrossRef] [PubMed]
21.	M. S. Sussman, N. Robert, and G. A. Wright, “Adaptive averaging for improved SNR in real-time coronary artery MRI.” IEEE Trans. Med. Img. 23(8), 1034–1045 (2004). URL http://dx.doi.org/10.1109/TMI.2004.828677. [CrossRef]
22.	M. S. Sussman, A. B. Kerr, J. M. Pauly, N. Merchant, and G. A. Wright, “Tracking the motion of the coronary arteries with the correlation coefficient,” in Proceedings of the International Society for Magnetic Resonance in Medicine (1999).
23.	J. P. Lewis, “Fast normalized cross-correlation,” in Vision Interface (1995). URL http://citeseer.ist.psu.edu/lewis95fast.html.
24.	H.-U. Dodt, U. Leischner, A. Schierloh, N. Jhrling, C. P. Mauch, K. Deininger, J. M. Deussing, M. Eder, W. Zieglgnsberger, and K. Becker, “Ultramicroscopy: three-dimensional visualization of neuronal networks in the whole mouse brain.” Nat. Methods 4(4), 331–336 (2007). URL http://dx.doi.org/10.1038/nmeth1036. [CrossRef] [PubMed]
25.	A. Maass, P. L. von Leithner, V. Luong, L. Guo, T. E. Salt, F. W. Fitzke, and M. F. Cordeiro, “Assessment of rat and mouse RGC apoptosis imaging in vivo with different scanning laser ophthalmoscopes.” Curr. Eye Res. 32(10), 851–861 (2007). URL http://dx.doi.org/10.1080/02713680701585872. [CrossRef] [PubMed]
26.	P. Soille, Morphological Image Analysis (Springer-Verlag, 1999). URL http://web.ukonline.co.uk/soille/book1st.

OCIS Codes
(100.1830) Image processing : Deconvolution
(100.2000) Image processing : Digital image processing
(110.4280) Imaging systems : Noise in imaging systems
(330.2210) Vision, color, and visual optics : Vision - eye movements
(330.4300) Vision, color, and visual optics : Vision system - noninvasive assessment
(330.4460) Vision, color, and visual optics : Ophthalmic optics and devices

ToC Category:
Image Processing

History
Original Manuscript: January 15, 2008
Revised Manuscript: May 11, 2008
Manuscript Accepted: May 14, 2008
Published: May 22, 2008

Virtual Issues
Vol. 3, Iss. 6 Virtual Journal for Biomedical Optics

Citation
S. Kumar, G. Ho, K. M. Woo, and L. Zhuo, "Achieving cellular resolution for in vivo retinal images of transgenic GFAP-GFP mice via image processing," Opt. Express 16, 8250-8262 (2008)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-16-11-8250

Sort: Author | Year | Journal | Reset

References

R. A. Fields, M. Birnbaum, and C. L. Fincher, "Highly efficient Nd:YVO4 diode-laser end-pumped laser," Appl. Phys. Lett. 51, 1885-1886 (1987). [CrossRef] [PubMed]
A. Sennaroglu, "Efficient continuous-wave operation of a diode-pumped Nd:YVO4 laser at 1342 nm," Opt. Commun. 164, 191-197 (1999). [CrossRef] [PubMed]
Y. F. Chen, "Design criteria for concentration optimization in scaling diode end-pumped lasers to high powers: influence of thermal fracture," IEEE J. Quantum Electron. 35, 234-239 (1999). [CrossRef]
A. Y. Yao, W. Hou, Y. P. Kong, L. Guo, L. A. Wu, R. N. Li, D. F. Cui, Z. Y. Xu, Y. Bi, and Y. Zhou, " Double-end-pumped 11-W Nd:YVO4 cw laser at 1342 nm," J. Opt. Soc. Am. B 22, 2129-2133 (2005). [CrossRef] [PubMed]
T. Jensen, V. G. Ostroumov, J. P. Meyn, G. Huber, A. I. Zagumennyi, and I. A. Shcherbarkov, "Spectroscopic characterization and laser performance of diode-laser-pumped Nd:GdVO4," Appl. Phys. B 58, 373-379 (1994). [CrossRef]
H. J. Zhang, X. L. Meng, L. Zhu, H. Z. Zhang, P. Wang, J. Dawes, C. Q. Wang, and Y. T. Chow, "Investigation on the growth and laser properties of Nd:GdVO4 single crystal," Cryst. Res. Technol. 33, 801-806 (1998). [CrossRef]
T. Ogawa, Y. Urata, S. Wada, K. Onodera, H. Machida, H. Sagae, M. Higuchi, and K. Kodaira, "Efficient laser performance of Nd:GdVO4 crystals grown by the floating zone method," Opt. Lett. 28, 2333-2335 (2003). [CrossRef] [PubMed]
V. Lupei, N. Pavel, Y. Sato, and T. Taira, "Highly efficient 1063-nm continuous-wave laser emission in Nd:GdVO4," Opt. Lett. 28, 2366-2368 (2003). [CrossRef] [PubMed]
L. J. Qin, X. L. Meng, L. Zhu, J. H. Liu, B. C. Xu, H. Z. Xu, F. Y. Jiang, C. L. Du, X. Q. Wang, and Z. S. Shao, "Influence of the different Gd/Y ratio on the properties of Nd:YxGd1-xVO4 mixed crystals," Chem. Phys. Lett. 380, 273-278 (2003). [CrossRef]
J. Liu, X. Meng, Z. Shao, M. Jiang, B. Ozygus, A. Ding, and H. Weber, "Pulse energy enhancement in passive Q-switching operation with a class of Nd:GdxY1-xVO4 crystals," Appl. Phys. Lett. 83, 1289-1291 (2003). [CrossRef] [PubMed]
Y. F. Chen, M. L. Ku, L. Y. Tsai, and Y. C. Chen, "Diode-pumped passively Q-switched picosecond Nd:GdxY1-xVO4 self-stimulated Raman laser," Opt. Lett. 29, 2279-2281 (2004).
J. L He, Y. X. Fan, J. Du, Y. G. Wang, S. Liu, H. T. Wang, L. H. Zhang, and Y. Hang, "4-ps passively mode-locked Nd:Gd0.5Y0.5VO4 laser with a semiconductor saturable-absorber mirror," Opt. Lett. 29, 2803-2805 (2004).
S. P. Ng, D. Y. Tang, A. Q. Liu, L. J. Qin, and X. L. Meng, "Short pulse passively Q-switched NdGdYVO4 laser using a GaAs mirror, " Opt. Commun. 259, 256-260 (2006). [CrossRef]
H. H. Yu, H. J. Zhang, Z. P. Wang, J. Y. Wang, Y. G. Yu, Z. S. Shao, M. H. Jiang, and X. Y. Zhang, "Continuous wave and passively Q-switched laser performance of Nd-doped mixed crystal Nd:Lu0.5Gd0.5VO4," Appl. Phys. Lett. 90, 231110 (2007).
A. A. Kaminskii, "Modern developments in the physics of crystalline laser materials," Phys. Status Solidi (A) 200, 215-296 (2003). [CrossRef]
T. T. Basiev and R. C. Powell, "Special issue on solid state Raman lasers-Introduction," Opt. Mater. 11, 301-306 (1999). [CrossRef] [PubMed]
H. M. Pask, "The design and operation of solid-state Raman lasers," Prog. Quantum Electron. 27, 3-56 (2003). [CrossRef] [PubMed]
Y. F. Chen, "Stimulated Raman scattering in a potassium titanyl phosphate crystal: simultaneous self-sum frequency mixing and self-frequency doubling," Opt. Lett. 30, 400-402 (2005).
S. Pearce, C. L. M. Ireland, and P. E. Dyer, "Solid-state Raman laser generating <1ns, multi-kilohertz pulses at 1096 nm," Opt. Commun. 260, 680-686 (2006).
V. Pasiskevicius, C. Canalias, and F. Laurell, "Highly efficient stimulated Raman scattering of picosecond pulses in KTiOPO4," Appl. Phys. Lett. 88, 041110 (2006). [CrossRef] [PubMed]
G. A. Massey, T. M. Loehr, L. J. Willis, and J. C. Johnson, ""Raman and Electro-optic Properties of Potassium Titanate Phosphate," Appl. Opt. 19, 4136-4137 (1980). [CrossRef]
K. Kawase, M. Mizuno, S. Sohma, H. Takahashi, T. Taniuchi, Y. Urata, S. Wada, H. Tashiro, and H. Ito, "Difference-frequency terahertz-wave generation from 4-dimethylamino-N-methyl-4-stilbazolium-tosylate by use of an electronically tuned Ti:sapphire laser," Opt. Lett. 24, 1065-1067 (1999).
K. Kawase, T. Hatanaka, H. Takahashi, K. Nakamura, T. Taniuchi, and H. Ito, "Tunable terahertz-wave generation from DAST crystal by dual signal-wave parametric oscillation of periodically poled lithium niobate," Opt. Lett. 25, 1714-1716 (2000).
K. Suizu, K. Miyamoto, T. Yamashita, and H. Ito, "High-power terahertz-wave generation using DAST crystal and detection using mid-infrared powermeter," Opt. Lett. 32, 2885-2887 (2007). [CrossRef] [PubMed]
J. J. Zayhowski and A. Mooradian, "Single-frequency microchip Nd lasers," Opt. Lett. 14, 24-26 (1989). [CrossRef] [PubMed]
G. J. Dixon, L. S. Lingvay, and R. H. Jarman, "Properties of close coupled monolithic, lithium neodymium, tetraphosphate lasers," Proc. SPIE 1104, 107 (1989).
Y. F. Chen, T. M. Huang, C. F. Kao, C. L. Wang, and S. C. Huang, "Optimization in scaling fiber-coupled laser-diode end-pumped lasers to higher power: influence of thermal effect," IEEE J. Quantum Electron. 33, 1424-1429 (1997).
Y. F. Chen, "High-power diode-pumped Q-switched intracavity frequency-doubled Nd:YVO4 laser with a sandwich-type resonator," Opt. Lett. 24, 1032-1034 (1999).

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Click here to see a list of articles that cite this paper