c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
journal homepage: www.intl.elsevierhealth.com/journals/cmpb
Novel mathematical algorithm for pupillometric data analysis Matthew C. Canver a,b,1 , Adam C. Canver c,1 , Karen E. Revere b , Defne Amado b , Jean Bennett b , Daniel C. Chung b,∗ a
Department of Bioengineering, University of Pennsylvania, Philadelphia, PA, United States F.M. Kirby Center for Molecular Ophthalmology, Scheie Eye Institute, University of Pennsylvania, Philadelphia, PA, United States c School of Biomedical Engineering, Drexel University College of Medicine, Philadelphia, PA, United States b
a r t i c l e
i n f o
a b s t r a c t
Article history:
Pupillometry is used clinically to evaluate retinal and optic nerve function by measuring
Received 3 August 2012
pupillary response to light stimuli. We have developed a mathematical algorithm to auto-
Received in revised form
mate and expedite the analysis of non-filtered, non-calculated pupillometric data obtained
20 June 2013
from mouse pupillary light reflex recordings, obtained from dynamic pupillary diameter
Accepted 18 August 2013
recordings following exposure of varying light intensities. The non-filtered, non-calculated
Keywords:
olding is used to remove data caused by eye blinking, loss of pupil tracking, and/or head
Pupillometry
movement. Twelve physiologically relevant parameters were extracted from the collected
pupillometric data is filtered through a low pass finite impulse response (FIR) filter. Thresh-
Algorithm
data: (1) baseline diameter, (2) minimum diameter, (3) response amplitude, (4) re-dilation
Retinal function
amplitude, (5) percent of baseline diameter, (6) response time, (7) re-dilation time, (8) average constriction velocity, (9) average re-dilation velocity, (10) maximum constriction velocity, (11) maximum re-dilation velocity, and (12) onset latency. No significant differences were noted between parameters derived from algorithm calculated values and manually derived results (p ≥ 0.05). This mathematical algorithm will expedite endpoint data derivation and eliminate human error in the manual calculation of pupillometric parameters from nonfiltered, non-calculated pupillometric values. Subsequently, these values can be used as reference metrics for characterizing the natural history of retinal disease. Furthermore, it will be instrumental in the assessment of functional visual recovery in humans and preclinical models of retinal degeneration and optic nerve disease following pharmacological or gene-based therapies. © 2013 Elsevier Ireland Ltd. All rights reserved.
1.
Introduction
Pupillometry is a non-invasive imaging technique that is used to assess pupillary light reflexes [1]. The pupil is innervated
∗
solely by the autonomic nervous system, allowing for assessment of visual function with respect to the central and peripheral autonomic systems [2]. Clinically, pupillometry is used as an objective means of evaluating the integrity of the retinal-cortical vision pathways. Currently, pupillometry
Corresponding author at: F.M. Kirby Center for Molecular Ophthalmology, Scheie Eye Institute, University of Pennsylvania, 422 Curie Blvd., #310, Philadelphia, PA 19104, United States. Tel.: +1 215 898 0915; fax: +1 215 573 7155. E-mail address: [email protected] (D.C. Chung). 1 These authors contributed equal work. 0169-2607/$ – see front matter © 2013 Elsevier Ireland Ltd. All rights reserved. http://dx.doi.org/10.1016/j.cmpb.2013.08.008
222
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
is used as an outcome measurement in clinical trials following gene replacement therapy for retinal degeneration, indicating the restoration of pupillary function as a direct indication of improved retinal function [3–5]. However, only recently has commercially made recording equipment been available for pupillometric evaluation involving laboratory animals, such as mice (Neuroptics Inc., San Clemente, CA). Pupillometry can be an important evaluative parameter for various retinal and neurodegenerative disorders [6,7]. This method for pupillometric data collection utilizes a dual-camera system and infrared illumination to record pupillary responses to light stimuli. The Neuroptics pupillometer tracks pupillary diameter through the differences in threshold detection, which refers to the differences in hue of the dark pupil and lighter iris tissue. When light stimuli are applied, pupillary constriction results, manifesting as a decrease in the pupillary diameter followed by a period of pupillary recovery prior to a subsequent stimulus. Constriction of the pupil is initiated after light-sensing photoreceptor cells in the retina are stimulated and an electrical impulse is relayed through bipolar, horizontal, amacrine, and ganglion cells, and ultimately conducted through the optic nerve. The signal then travels via the optic tract to synapse in the ipsilateral lateral geniculate nucleus of the thalamus and the ipsilateral pretectal region of the rostral midbrain. The latter set of afferent neurons synapse bilaterally at the Edinger-Westphal nuclei of the midbrain. These neurons contribute to pre-ganglionic parasympathetic efferent fibers in the oculomotor nerve (cranial nerve III) which projects back towards the brain and synapse in the ciliary ganglion. When stimulated, the post-ganglionic short ciliary nerves act to constrict the pupil. The pupillary light reflex involves multiple synaptic connections and is maintained only if all components of the retina and optic nerve are intact [8]. Here we demonstrate how a mathematical algorithm for pupillometric data analysis can be used to expedite, automate and standardize pupillometry data analysis and establish normative pupillometry parameters that are highly accurate for mouse models and humans. Normative values can ultimately be used to assess efficacy of pharmacologic and gene therapy treatments of retinal or neurodegenerative diseases.
2.
Methods and materials
2.1.
Experimental design
Mouse pupillometric measurements were used to develop the algorithm. The Neuroptics Pupillometer (San Clemente, CA, USA) utilizes a dual-camera system and infrared wavelengths to record pupillary responses to light stimuli. The pupillometer tracks pupil diameter through threshold detection. The Neuroptics software sampled at 15 Hz frequency. The data used to develop the program involved a pulse duration of 100 ms and a relaxation time of 9900 ms. Wild-type C57Bl6 mice were dark-adapted for a minimum of 6 h prior to testing and were not anesthetized during the procedure. Mice were manipulated multiple times prior to measurement to acclimate to hand held restraining methods. The mice were held in a prone position on the pupillometer platform and cameras were aligned optically to each eye. Three light pulses were
Fig. 1 – Representative tracings of pupillary light response during the three recording periods of pupillary constriction, as mapped: time vs. pupillary diameter.
administered at 10 s intervals for each of three intensities: 4.6, 37.5, 300 W/cm2 . Each eye was flashed independently, but consensual responses were obtained for each eye during unilateral light flashes. Each intensity was performed sequentially on a given animal for a total experiment duration of 180 s (3 trials per eye and 30 s per trial). Five mice were utilized to obtain pupillometric data points, and separate investigators were blinded to specific data sets and responsible for the calculation using either a standard manual mathematical calculations or our novel mathematical algorithm. The mice were maintained in clean micro-isolator cages and given standard laboratory feed and water ad libitum. All experimental protocols were reviewed and approved by the Institutional Animal Care and Use Committee (IACUC) at the University of Pennsylvania and followed the guidelines set forth by the National Institutes of Health “Guide for Care and Use of Laboratory Animals.”
2.2.
Data importation and digital filtration
The calculation algorithm was developed as a custom program created for MATLAB MathWorks Software (MathWorks, Natick, MA, USA). The initial non-calculated pupillometric values were uploaded into MATLAB and separated by left and right pupils. Data between stimuli were separated into individual stimulus-response segments for analysis (Fig. 1). The non-filtered, non-calculated values obtained from the pupillometer, also referred to as raw data, from the left and right pupils were plotted against time. Noise was removed from the signal using a 20-term low-pass digital finite impulse response (FIR) filter with a cutoff frequency (fc ) of 1 Hz. The cutoff frequency was established based on the calculation of the Fast Fourier Transform (FFT) and a plot of the magnitude response (data not shown).
2.3. Removal of data caused by eye blinking, loss of pupil tracking, and/or head movement During data collection, eye blinking, loss of pupil tracking, and/or head movement can cause non-physiological changes
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
223
in pupil diameters. This data was removed based on establishment of a threshold value. The data was numerically differentiated using the central difference approximation: f (x) =
f (xi+1 ) − f (xi−1 ) 2h
and the forward difference equation: f (xi ) =
f (xi+1 ) − f (x) h
and the backward difference equation: f (xi ) =
f (x) − f (xi−1 ) h
where x represents the first point in a 2 point estimation and h represents the change in x. The first derivative (velocity) and the second derivative (acceleration) were calculated from the raw data. The average acceleration and its standard deviation were calculated. The threshold was set to be a multiple of the standard deviation according to the formula: threshold = k(standard deviation), where k is a constant and does not need to be an integer. If any point within a given stimulus response was above the threshold, then the entire stimulus-response segment was excluded from the parameter calculations. This thresholding methodology efficiently distinguished between physiologically relevant accelerations and acceleration resulting from procedural artifacts such as from eye blinking, loss of pupil tracking, and/or head movement. The data is plotted in both pre- and postfiltration form, as well as after the removal of data (referred to as the useable data plot), as represented by the algorithm in Fig. 2.
Fig. 2 – Schematic algorithm outlining the steps involved in obtaining final parameter values from raw pupillary light response data.
2.4.
(8) Average constriction velocity was calculated as:
Parameter calculation
Each parameter of interest from the data was calculated for each stimulus-response segment between light flashes. The following parameters were extracted: (1) baseline diameter (BD), (2) minimum diameter (MD), (3) response amplitude (RA), (4) re-dilation amplitude (DA), (5) percent of baseline diameter (%BD), (6) response time (RT), (7) re-dilation time (DT), (8) average constriction velocity, (9) average re-dilation velocity, (10) maximum constriction velocity, and (11) maximum re-dilation velocity, and (12) onset latency. Each parameter was tabulated along with the mean and standard deviation. (1) BD was determined to be the average of the pupil diameter 100 ms prior to the light stimulus [9]. (2) MD was the lowest pupil diameter after pupillary constriction (Fig. 3). (3) RA was calculated as:
and DT was calculated as: tstimulus(n+1) − tminimum of stimulus(n)
average constriction velocity =
RA RT
(9) Average dilation velocity was calculated as: average redilation velocity =
DA DT
RA = BD − MD (4) DA was calculated as: DA = BDn+1 − MDn (5) The percent constriction was calculated as: percent of BD =
MD BD
100%
(6) and (7) RT was calculated as: tstimulus − tminimum amplitude ,
Fig. 3 – Illustrates the corresponding regions for calculated pupillary light response parameters. Minimum diameter (MD), response amplitude (RA), re-dilation amplitude (DA), response time (RT), re-dilation time (DT) and latency (L).
224
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
(10) The maximum constriction velocity was calculated by obtaining the first derivative (velocity) prior to the point of minimum amplitude. (11) The maximum re-dilation velocity was calculated by determining the maximum of the first derivative (velocity) after the point of minimum amplitude, but prior to the subsequent light stimuli. (12) Onset latency was calculated by determining the maximum acceleration following the stimulus [9].
2.5.
Evaluation of algorithm
Each of the 12 parameters was obtained from 40 stimulusresponse segments (n = 40) using the algorithm and subsequently compared to the same parameters calculated manually using a paired Students t-test (˛ = 0.05). P-values greater than 0.05 signified that the manually calculated parameters do not differ from the algorithm generated values. Manual calculations were initiated by obtaining baseline pupillary (BD) and minimum diameters (MD), and calculations for RA, DA and percent of BD were calculated. Time of amplitude change were obtained to calculate RT, DT, velocity means and latency. Formulas used were identical to those programmed into the novel algorithm.
2.6.
Sensitivity and specificity
Artifacts from eye blinking, loss of pupil tracking, and/or head movement were visually identified in forty-two data files and induced aberrations in the raw data. Stimulus-response segments from each eye were treated independently, thereby making n = 252 (6 stimulus responses per mouse, 3 from both the left and right pupil). A true positive was defined as accepting stimulus-response segments that were not contaminated by artifacts. A true negative was defined as the removal of stimulus-response segments that were contaminated by artifacts. A false positive was defined as the removal of stimulus-response segments not contaminated by artifact. A false negative was defined as the maintenance of data that was contaminated by artifacts. Sensitivity and specificity were calculated using the formulas: # of true negatives specificity = # of true negatives + # of false positives specificity =
# of true positives . # of true positives + # of false negatives
3.
Results
3.1.
Evaluation of algorithm
The calculated parameters from the novel algorithm were compared to manually-derived values from the raw data obtained from the pupillometry apparatus. The comparison of the derived parameters from the two sources was conducted using the paired Students t-test, which yielded non-significant p-values for all parameters (n = 40) (Table 1). The comparison of the data obtained by the custom algorithm and the manual method demonstrate that the novel
Table 1 – Comparison of calculated parameters from the novel algorithm compared to manually-derived values from the raw data obtained from the pupillometry apparatus. The comparison of the derived parameters from the two sources were conducted using the paired Students t-test, which yielded non-significant p-values for all parameters (n = 40). Parameter RA MD BD Percent constriction DA Average constriction velocity Average redilation velocity RT DT Maximum constriction velocity Maximum re-dilation velocity Onset latency
p-Value 0.266 0.323 0.277 0.207 0.0879 0.317 0.320 0.422 0.201 0.135 0.661 0.299
algorithm is accurate and expeditious in calculating the parameters of interest. The algorithm provides a streamlined, time- and cost-effective means in which to validate the accuracy of the gathered data. In Table 1, the data suggests that there is no statistically significant (p ≥ 0.05) difference between manually calculated and novel algorithm derived parameters.
3.2.
Sensitivity and specificity
Percent sensitivity and specificity are plotted as a function of the standard deviation multiplier k used to establish the threshold. When artifacts are visually present, it appears that the optimal threshold is between one and two standard deviation multipliers. As k increases, specificity increases and sensitivity decreases. A larger value of k should be used conservatively, as an increase in k decreases the likelihood of data removal.
4.
Discussion
The p-values (p ≥ 0.05) demonstrate that there is no significant difference between the parameter values generated by the novel algorithm as compared to the parameter values calculated manually. Thus the described algorithm can be used to reliably extract a variety of parameters from raw pupillometric data generated by the Neuroptics Pupillometer. The algorithm can be used to establish normative values of pupillary response for numerous mouse strains, both wild type and transgenic models. The parameters were chosen based on their relevance in disorders of retinal degeneration and neurodegenerative disease [6,7]. Subsequent application of these values can be applied to assessment of visual changes following pharmacological or gene-based interventions. The 12 parameters involved in the pupillary light reflex are vital for the interpretation and significance of each pupillometric response. These parameters give insight into the health of the retina as well as into the autonomic afferent and
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
efferent neurological pathways necessary for the pupillary light reflex and for vision. In addition to using this approach to characterize the retinal phenotype of animal models of retinal degeneration, this algorithm will also be useful for evaluating effects of intervention on both the afferent and efferent portions of the retinal/cortical pathways. Thus, it will serve as an outcome measure for experimental therapies for retinal degeneration, for evaluating disease and disease treatment strategies in animals with optic nerve disorders such as optic neuritis, and for evaluating the integrity of the efferent parasympathetic pathway. Furthermore, pupillometry provides a more sensitive objective measure of retinal function than full field electroretinograms [3–6]. Pupillometry can also assess optic nerve in inflammatory states (optic neuritis) [10,11] and in glaucoma [10,12–14]. Finally, pupillometry can also assess the integrity of the central nervous system efferent pathway.
5.
Conclusion
Since pupillometry is beset with background and subjectderived aberrations, this algorithm can reliably remove data altered by eye blinking, loss of pupil tracking, and/or head movement. Based on our findings, the optimal threshold for removal of aberrant data is a value of k between one and two (one and two standard deviations) when artifacts are visually present. However, a threshold system also allows the experimenter to control the sensitivity and/or specificity of the algorithm. Accurate and expeditious derivation of data is essential to any study involving pre-clinical models of retinal or optic nerve disease. This novel algorithm demonstrates both accuracy and speed in calculating 12 pupillometric parameters. Thus, it has great potential to further the use of pupillometry as an expeditious, non-invasive method for the evaluation of retinal function and neural recovery in humans and preclinical models of retinal degeneration and optic nerve disease.
Acknowledgments The authors would like to thank Dr. Kenneth Foster, Professor of Bioengineering, School of Engineering and Applied Science, University of Pennsylvania, for his advice concerning filter choice and design. The authors would also like to thank Dr. Gui-Shuang Ying for his advice concerning statistical analysis. This work was supported by the Foundation Fighting Blindness which sponsors the CHOP-Penn Pediatric Center for Retinal Degenerations, Research to Prevent Blindness (JB and departmental funds), anonymous donor, NIH K08 EY017024, Hope for Vision Foundation, Fidelity Charitable Trust, (DC), NIH 1F30AG030961-01 (DA), the Paul and Evanina Mackall Foundation Trust at the Scheie Eye Institute, and the FM Kirby Foundation.
225
references
[1] R. Hussain, S. Hopkins, E. Frohman, T. Eagar, P. Cravens, B. Greenberg, S. Vernino, O. Stuve, Direct and consensual murine pupillary reflex metrics: establishing normative values, Autonomic Neuroscience 151 (2009) 164–167. [2] P. Bitsios, R. Prettyman, E. Szabadi, Changes in autonomic function with age: a study of pupillary kinetics in healthy young and old people, Age and Ageing 25 (1996) 432–438. [3] A. Maguire, F. Simonelli, E. Pierce, E. Pugh, F. Mingozzi, J. Bennicelli, et al., Safety and efficacy of gene transfer for Leber’s congenital amaurosis, The New England Journal of Medicine 358 (2008) 2240–2248. [4] M. Maguire, K. High, A. Auricchio, J. Wright, E. Pierce, F. Testa, F. Mingozzi, J. Bennicelli, G. Ying, S. Rossi, A. Fulton, K. Marshall, S. Banfi, D. Chung, J. Morgan, B. Hauck, O. Zelenaia, X. Zhu, L. Raffini, F. Coppieters, E. De Baere, K. Shindler, N. Volpe, E. Surace, C. Acerra, A. Lyubarsky, T. Redmond, E. Stone, J. Sun, J. McDonnell, B. Leroy, F. Simonelli, J. Bennett, Age-dependent effects of RPE65 gene therapy for Leber’s congenital amaurosis: a phase 1 dose-escalation trial, Lancet 374 (2009) 1597–1605. [5] F. Simonelli, A. Maguire, F. Testa, E. Pierce, F. Mingozzi, J. Bennicelli, S. Rossi, K. Marshall, S. Banfi, E. Surace, J. Sun, T. Redmond, X. Zhu, K. Shindler, J. Morgan, G. Ying, C. Acerra, J. Wright, J. McDonnell, K. High, J. Bennett, A. Auricchio, Gene therapy for Leber’s congenital amaurosis is safe and effective through 1.5 years after vector administration, Molecular Therapy 18 (2010) 643–650. [6] S. Shindler, K. Revere, M. Dutt, S. Ying, D. Chung, In vivo detection of experimental optic neuritis by pupillometry, Experimental Eye Research 100 (July) (2012) 1–6. [7] T. Aleman, S. Jacobson, J. Chico, M. Scott, A. Cheung, et al., Impairment of the transient pupillary light reflex in Rpe65(−/−) mice and humans with leber congenital amaurosis, Investigative Ophthalmology and Visual Science 45 (2004) 1259–1271. [8] J. Surakka, J. Ruutiainen, A. Romberg, P. Puukka, E. Kronholm, H. Karanko, Pupillary function in early multiple sclerosis, Clinical Autonomic Research 18 (2008) 150–154. [9] D. Purves, G. Augustine, D. Fitzpatrick, W. Hall, A. LaMantia, J. McNamara, L. White, Neuroscience, (4th ed.), Sinauer Associates, Inc., Sunderland, 2008. [10] O. Bergamin, R. Kardon, Latency of the pupil light reflex: sample rate, stimulus intensity, and variation in normal subjects, Investigative Ophthalmology and Visual Science 44 (2003) 1546–1554. [11] C. Ellis, The afferent pupillary defect in acute optic neuritis, Journal of Neurology, Neurosurgery & Psychiatry 42 (1979) 1008–1017. [12] R. Beck, P. Cleary, M. Anderson, J. Keltner, W. Shults, D. Kaufman, E. Buckley, J. Corbett, M. Kupersmith, N. Miller, et al., A randomized, controlled trial of corticosteroids in the treatment of acute optic neuritis. The optic neuritis study group, New England Journal of Medicine 326 (1992) 581–588. [13] Y. Chen, H. Wyatt, W. Swanson, M. Dul, Rapid pupil-based assessment of glaucomatous damage, Optometry & Vision Science 85 (2008) 471–481. [14] L. Kalaboukhova, V. Fridhammar, B. Lindblom, Relative afferent pupillary defect in glaucoma: a pupillometric study, Acta Ophthalmologica Scandinavica 85 (2007) 519–525.
journal homepage: www.intl.elsevierhealth.com/journals/cmpb
Novel mathematical algorithm for pupillometric data analysis Matthew C. Canver a,b,1 , Adam C. Canver c,1 , Karen E. Revere b , Defne Amado b , Jean Bennett b , Daniel C. Chung b,∗ a
Department of Bioengineering, University of Pennsylvania, Philadelphia, PA, United States F.M. Kirby Center for Molecular Ophthalmology, Scheie Eye Institute, University of Pennsylvania, Philadelphia, PA, United States c School of Biomedical Engineering, Drexel University College of Medicine, Philadelphia, PA, United States b
a r t i c l e
i n f o
a b s t r a c t
Article history:
Pupillometry is used clinically to evaluate retinal and optic nerve function by measuring
Received 3 August 2012
pupillary response to light stimuli. We have developed a mathematical algorithm to auto-
Received in revised form
mate and expedite the analysis of non-filtered, non-calculated pupillometric data obtained
20 June 2013
from mouse pupillary light reflex recordings, obtained from dynamic pupillary diameter
Accepted 18 August 2013
recordings following exposure of varying light intensities. The non-filtered, non-calculated
Keywords:
olding is used to remove data caused by eye blinking, loss of pupil tracking, and/or head
Pupillometry
movement. Twelve physiologically relevant parameters were extracted from the collected
pupillometric data is filtered through a low pass finite impulse response (FIR) filter. Thresh-
Algorithm
data: (1) baseline diameter, (2) minimum diameter, (3) response amplitude, (4) re-dilation
Retinal function
amplitude, (5) percent of baseline diameter, (6) response time, (7) re-dilation time, (8) average constriction velocity, (9) average re-dilation velocity, (10) maximum constriction velocity, (11) maximum re-dilation velocity, and (12) onset latency. No significant differences were noted between parameters derived from algorithm calculated values and manually derived results (p ≥ 0.05). This mathematical algorithm will expedite endpoint data derivation and eliminate human error in the manual calculation of pupillometric parameters from nonfiltered, non-calculated pupillometric values. Subsequently, these values can be used as reference metrics for characterizing the natural history of retinal disease. Furthermore, it will be instrumental in the assessment of functional visual recovery in humans and preclinical models of retinal degeneration and optic nerve disease following pharmacological or gene-based therapies. © 2013 Elsevier Ireland Ltd. All rights reserved.
1.
Introduction
Pupillometry is a non-invasive imaging technique that is used to assess pupillary light reflexes [1]. The pupil is innervated
∗
solely by the autonomic nervous system, allowing for assessment of visual function with respect to the central and peripheral autonomic systems [2]. Clinically, pupillometry is used as an objective means of evaluating the integrity of the retinal-cortical vision pathways. Currently, pupillometry
Corresponding author at: F.M. Kirby Center for Molecular Ophthalmology, Scheie Eye Institute, University of Pennsylvania, 422 Curie Blvd., #310, Philadelphia, PA 19104, United States. Tel.: +1 215 898 0915; fax: +1 215 573 7155. E-mail address: [email protected] (D.C. Chung). 1 These authors contributed equal work. 0169-2607/$ – see front matter © 2013 Elsevier Ireland Ltd. All rights reserved. http://dx.doi.org/10.1016/j.cmpb.2013.08.008
222
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
is used as an outcome measurement in clinical trials following gene replacement therapy for retinal degeneration, indicating the restoration of pupillary function as a direct indication of improved retinal function [3–5]. However, only recently has commercially made recording equipment been available for pupillometric evaluation involving laboratory animals, such as mice (Neuroptics Inc., San Clemente, CA). Pupillometry can be an important evaluative parameter for various retinal and neurodegenerative disorders [6,7]. This method for pupillometric data collection utilizes a dual-camera system and infrared illumination to record pupillary responses to light stimuli. The Neuroptics pupillometer tracks pupillary diameter through the differences in threshold detection, which refers to the differences in hue of the dark pupil and lighter iris tissue. When light stimuli are applied, pupillary constriction results, manifesting as a decrease in the pupillary diameter followed by a period of pupillary recovery prior to a subsequent stimulus. Constriction of the pupil is initiated after light-sensing photoreceptor cells in the retina are stimulated and an electrical impulse is relayed through bipolar, horizontal, amacrine, and ganglion cells, and ultimately conducted through the optic nerve. The signal then travels via the optic tract to synapse in the ipsilateral lateral geniculate nucleus of the thalamus and the ipsilateral pretectal region of the rostral midbrain. The latter set of afferent neurons synapse bilaterally at the Edinger-Westphal nuclei of the midbrain. These neurons contribute to pre-ganglionic parasympathetic efferent fibers in the oculomotor nerve (cranial nerve III) which projects back towards the brain and synapse in the ciliary ganglion. When stimulated, the post-ganglionic short ciliary nerves act to constrict the pupil. The pupillary light reflex involves multiple synaptic connections and is maintained only if all components of the retina and optic nerve are intact [8]. Here we demonstrate how a mathematical algorithm for pupillometric data analysis can be used to expedite, automate and standardize pupillometry data analysis and establish normative pupillometry parameters that are highly accurate for mouse models and humans. Normative values can ultimately be used to assess efficacy of pharmacologic and gene therapy treatments of retinal or neurodegenerative diseases.
2.
Methods and materials
2.1.
Experimental design
Mouse pupillometric measurements were used to develop the algorithm. The Neuroptics Pupillometer (San Clemente, CA, USA) utilizes a dual-camera system and infrared wavelengths to record pupillary responses to light stimuli. The pupillometer tracks pupil diameter through threshold detection. The Neuroptics software sampled at 15 Hz frequency. The data used to develop the program involved a pulse duration of 100 ms and a relaxation time of 9900 ms. Wild-type C57Bl6 mice were dark-adapted for a minimum of 6 h prior to testing and were not anesthetized during the procedure. Mice were manipulated multiple times prior to measurement to acclimate to hand held restraining methods. The mice were held in a prone position on the pupillometer platform and cameras were aligned optically to each eye. Three light pulses were
Fig. 1 – Representative tracings of pupillary light response during the three recording periods of pupillary constriction, as mapped: time vs. pupillary diameter.
administered at 10 s intervals for each of three intensities: 4.6, 37.5, 300 W/cm2 . Each eye was flashed independently, but consensual responses were obtained for each eye during unilateral light flashes. Each intensity was performed sequentially on a given animal for a total experiment duration of 180 s (3 trials per eye and 30 s per trial). Five mice were utilized to obtain pupillometric data points, and separate investigators were blinded to specific data sets and responsible for the calculation using either a standard manual mathematical calculations or our novel mathematical algorithm. The mice were maintained in clean micro-isolator cages and given standard laboratory feed and water ad libitum. All experimental protocols were reviewed and approved by the Institutional Animal Care and Use Committee (IACUC) at the University of Pennsylvania and followed the guidelines set forth by the National Institutes of Health “Guide for Care and Use of Laboratory Animals.”
2.2.
Data importation and digital filtration
The calculation algorithm was developed as a custom program created for MATLAB MathWorks Software (MathWorks, Natick, MA, USA). The initial non-calculated pupillometric values were uploaded into MATLAB and separated by left and right pupils. Data between stimuli were separated into individual stimulus-response segments for analysis (Fig. 1). The non-filtered, non-calculated values obtained from the pupillometer, also referred to as raw data, from the left and right pupils were plotted against time. Noise was removed from the signal using a 20-term low-pass digital finite impulse response (FIR) filter with a cutoff frequency (fc ) of 1 Hz. The cutoff frequency was established based on the calculation of the Fast Fourier Transform (FFT) and a plot of the magnitude response (data not shown).
2.3. Removal of data caused by eye blinking, loss of pupil tracking, and/or head movement During data collection, eye blinking, loss of pupil tracking, and/or head movement can cause non-physiological changes
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
223
in pupil diameters. This data was removed based on establishment of a threshold value. The data was numerically differentiated using the central difference approximation: f (x) =
f (xi+1 ) − f (xi−1 ) 2h
and the forward difference equation: f (xi ) =
f (xi+1 ) − f (x) h
and the backward difference equation: f (xi ) =
f (x) − f (xi−1 ) h
where x represents the first point in a 2 point estimation and h represents the change in x. The first derivative (velocity) and the second derivative (acceleration) were calculated from the raw data. The average acceleration and its standard deviation were calculated. The threshold was set to be a multiple of the standard deviation according to the formula: threshold = k(standard deviation), where k is a constant and does not need to be an integer. If any point within a given stimulus response was above the threshold, then the entire stimulus-response segment was excluded from the parameter calculations. This thresholding methodology efficiently distinguished between physiologically relevant accelerations and acceleration resulting from procedural artifacts such as from eye blinking, loss of pupil tracking, and/or head movement. The data is plotted in both pre- and postfiltration form, as well as after the removal of data (referred to as the useable data plot), as represented by the algorithm in Fig. 2.
Fig. 2 – Schematic algorithm outlining the steps involved in obtaining final parameter values from raw pupillary light response data.
2.4.
(8) Average constriction velocity was calculated as:
Parameter calculation
Each parameter of interest from the data was calculated for each stimulus-response segment between light flashes. The following parameters were extracted: (1) baseline diameter (BD), (2) minimum diameter (MD), (3) response amplitude (RA), (4) re-dilation amplitude (DA), (5) percent of baseline diameter (%BD), (6) response time (RT), (7) re-dilation time (DT), (8) average constriction velocity, (9) average re-dilation velocity, (10) maximum constriction velocity, and (11) maximum re-dilation velocity, and (12) onset latency. Each parameter was tabulated along with the mean and standard deviation. (1) BD was determined to be the average of the pupil diameter 100 ms prior to the light stimulus [9]. (2) MD was the lowest pupil diameter after pupillary constriction (Fig. 3). (3) RA was calculated as:
and DT was calculated as: tstimulus(n+1) − tminimum of stimulus(n)
average constriction velocity =
RA RT
(9) Average dilation velocity was calculated as: average redilation velocity =
DA DT
RA = BD − MD (4) DA was calculated as: DA = BDn+1 − MDn (5) The percent constriction was calculated as: percent of BD =
MD BD
100%
(6) and (7) RT was calculated as: tstimulus − tminimum amplitude ,
Fig. 3 – Illustrates the corresponding regions for calculated pupillary light response parameters. Minimum diameter (MD), response amplitude (RA), re-dilation amplitude (DA), response time (RT), re-dilation time (DT) and latency (L).
224
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
(10) The maximum constriction velocity was calculated by obtaining the first derivative (velocity) prior to the point of minimum amplitude. (11) The maximum re-dilation velocity was calculated by determining the maximum of the first derivative (velocity) after the point of minimum amplitude, but prior to the subsequent light stimuli. (12) Onset latency was calculated by determining the maximum acceleration following the stimulus [9].
2.5.
Evaluation of algorithm
Each of the 12 parameters was obtained from 40 stimulusresponse segments (n = 40) using the algorithm and subsequently compared to the same parameters calculated manually using a paired Students t-test (˛ = 0.05). P-values greater than 0.05 signified that the manually calculated parameters do not differ from the algorithm generated values. Manual calculations were initiated by obtaining baseline pupillary (BD) and minimum diameters (MD), and calculations for RA, DA and percent of BD were calculated. Time of amplitude change were obtained to calculate RT, DT, velocity means and latency. Formulas used were identical to those programmed into the novel algorithm.
2.6.
Sensitivity and specificity
Artifacts from eye blinking, loss of pupil tracking, and/or head movement were visually identified in forty-two data files and induced aberrations in the raw data. Stimulus-response segments from each eye were treated independently, thereby making n = 252 (6 stimulus responses per mouse, 3 from both the left and right pupil). A true positive was defined as accepting stimulus-response segments that were not contaminated by artifacts. A true negative was defined as the removal of stimulus-response segments that were contaminated by artifacts. A false positive was defined as the removal of stimulus-response segments not contaminated by artifact. A false negative was defined as the maintenance of data that was contaminated by artifacts. Sensitivity and specificity were calculated using the formulas: # of true negatives specificity = # of true negatives + # of false positives specificity =
# of true positives . # of true positives + # of false negatives
3.
Results
3.1.
Evaluation of algorithm
The calculated parameters from the novel algorithm were compared to manually-derived values from the raw data obtained from the pupillometry apparatus. The comparison of the derived parameters from the two sources was conducted using the paired Students t-test, which yielded non-significant p-values for all parameters (n = 40) (Table 1). The comparison of the data obtained by the custom algorithm and the manual method demonstrate that the novel
Table 1 – Comparison of calculated parameters from the novel algorithm compared to manually-derived values from the raw data obtained from the pupillometry apparatus. The comparison of the derived parameters from the two sources were conducted using the paired Students t-test, which yielded non-significant p-values for all parameters (n = 40). Parameter RA MD BD Percent constriction DA Average constriction velocity Average redilation velocity RT DT Maximum constriction velocity Maximum re-dilation velocity Onset latency
p-Value 0.266 0.323 0.277 0.207 0.0879 0.317 0.320 0.422 0.201 0.135 0.661 0.299
algorithm is accurate and expeditious in calculating the parameters of interest. The algorithm provides a streamlined, time- and cost-effective means in which to validate the accuracy of the gathered data. In Table 1, the data suggests that there is no statistically significant (p ≥ 0.05) difference between manually calculated and novel algorithm derived parameters.
3.2.
Sensitivity and specificity
Percent sensitivity and specificity are plotted as a function of the standard deviation multiplier k used to establish the threshold. When artifacts are visually present, it appears that the optimal threshold is between one and two standard deviation multipliers. As k increases, specificity increases and sensitivity decreases. A larger value of k should be used conservatively, as an increase in k decreases the likelihood of data removal.
4.
Discussion
The p-values (p ≥ 0.05) demonstrate that there is no significant difference between the parameter values generated by the novel algorithm as compared to the parameter values calculated manually. Thus the described algorithm can be used to reliably extract a variety of parameters from raw pupillometric data generated by the Neuroptics Pupillometer. The algorithm can be used to establish normative values of pupillary response for numerous mouse strains, both wild type and transgenic models. The parameters were chosen based on their relevance in disorders of retinal degeneration and neurodegenerative disease [6,7]. Subsequent application of these values can be applied to assessment of visual changes following pharmacological or gene-based interventions. The 12 parameters involved in the pupillary light reflex are vital for the interpretation and significance of each pupillometric response. These parameters give insight into the health of the retina as well as into the autonomic afferent and
c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 3 ( 2 0 1 4 ) 221–225
efferent neurological pathways necessary for the pupillary light reflex and for vision. In addition to using this approach to characterize the retinal phenotype of animal models of retinal degeneration, this algorithm will also be useful for evaluating effects of intervention on both the afferent and efferent portions of the retinal/cortical pathways. Thus, it will serve as an outcome measure for experimental therapies for retinal degeneration, for evaluating disease and disease treatment strategies in animals with optic nerve disorders such as optic neuritis, and for evaluating the integrity of the efferent parasympathetic pathway. Furthermore, pupillometry provides a more sensitive objective measure of retinal function than full field electroretinograms [3–6]. Pupillometry can also assess optic nerve in inflammatory states (optic neuritis) [10,11] and in glaucoma [10,12–14]. Finally, pupillometry can also assess the integrity of the central nervous system efferent pathway.
5.
Conclusion
Since pupillometry is beset with background and subjectderived aberrations, this algorithm can reliably remove data altered by eye blinking, loss of pupil tracking, and/or head movement. Based on our findings, the optimal threshold for removal of aberrant data is a value of k between one and two (one and two standard deviations) when artifacts are visually present. However, a threshold system also allows the experimenter to control the sensitivity and/or specificity of the algorithm. Accurate and expeditious derivation of data is essential to any study involving pre-clinical models of retinal or optic nerve disease. This novel algorithm demonstrates both accuracy and speed in calculating 12 pupillometric parameters. Thus, it has great potential to further the use of pupillometry as an expeditious, non-invasive method for the evaluation of retinal function and neural recovery in humans and preclinical models of retinal degeneration and optic nerve disease.
Acknowledgments The authors would like to thank Dr. Kenneth Foster, Professor of Bioengineering, School of Engineering and Applied Science, University of Pennsylvania, for his advice concerning filter choice and design. The authors would also like to thank Dr. Gui-Shuang Ying for his advice concerning statistical analysis. This work was supported by the Foundation Fighting Blindness which sponsors the CHOP-Penn Pediatric Center for Retinal Degenerations, Research to Prevent Blindness (JB and departmental funds), anonymous donor, NIH K08 EY017024, Hope for Vision Foundation, Fidelity Charitable Trust, (DC), NIH 1F30AG030961-01 (DA), the Paul and Evanina Mackall Foundation Trust at the Scheie Eye Institute, and the FM Kirby Foundation.
225
references
[1] R. Hussain, S. Hopkins, E. Frohman, T. Eagar, P. Cravens, B. Greenberg, S. Vernino, O. Stuve, Direct and consensual murine pupillary reflex metrics: establishing normative values, Autonomic Neuroscience 151 (2009) 164–167. [2] P. Bitsios, R. Prettyman, E. Szabadi, Changes in autonomic function with age: a study of pupillary kinetics in healthy young and old people, Age and Ageing 25 (1996) 432–438. [3] A. Maguire, F. Simonelli, E. Pierce, E. Pugh, F. Mingozzi, J. Bennicelli, et al., Safety and efficacy of gene transfer for Leber’s congenital amaurosis, The New England Journal of Medicine 358 (2008) 2240–2248. [4] M. Maguire, K. High, A. Auricchio, J. Wright, E. Pierce, F. Testa, F. Mingozzi, J. Bennicelli, G. Ying, S. Rossi, A. Fulton, K. Marshall, S. Banfi, D. Chung, J. Morgan, B. Hauck, O. Zelenaia, X. Zhu, L. Raffini, F. Coppieters, E. De Baere, K. Shindler, N. Volpe, E. Surace, C. Acerra, A. Lyubarsky, T. Redmond, E. Stone, J. Sun, J. McDonnell, B. Leroy, F. Simonelli, J. Bennett, Age-dependent effects of RPE65 gene therapy for Leber’s congenital amaurosis: a phase 1 dose-escalation trial, Lancet 374 (2009) 1597–1605. [5] F. Simonelli, A. Maguire, F. Testa, E. Pierce, F. Mingozzi, J. Bennicelli, S. Rossi, K. Marshall, S. Banfi, E. Surace, J. Sun, T. Redmond, X. Zhu, K. Shindler, J. Morgan, G. Ying, C. Acerra, J. Wright, J. McDonnell, K. High, J. Bennett, A. Auricchio, Gene therapy for Leber’s congenital amaurosis is safe and effective through 1.5 years after vector administration, Molecular Therapy 18 (2010) 643–650. [6] S. Shindler, K. Revere, M. Dutt, S. Ying, D. Chung, In vivo detection of experimental optic neuritis by pupillometry, Experimental Eye Research 100 (July) (2012) 1–6. [7] T. Aleman, S. Jacobson, J. Chico, M. Scott, A. Cheung, et al., Impairment of the transient pupillary light reflex in Rpe65(−/−) mice and humans with leber congenital amaurosis, Investigative Ophthalmology and Visual Science 45 (2004) 1259–1271. [8] J. Surakka, J. Ruutiainen, A. Romberg, P. Puukka, E. Kronholm, H. Karanko, Pupillary function in early multiple sclerosis, Clinical Autonomic Research 18 (2008) 150–154. [9] D. Purves, G. Augustine, D. Fitzpatrick, W. Hall, A. LaMantia, J. McNamara, L. White, Neuroscience, (4th ed.), Sinauer Associates, Inc., Sunderland, 2008. [10] O. Bergamin, R. Kardon, Latency of the pupil light reflex: sample rate, stimulus intensity, and variation in normal subjects, Investigative Ophthalmology and Visual Science 44 (2003) 1546–1554. [11] C. Ellis, The afferent pupillary defect in acute optic neuritis, Journal of Neurology, Neurosurgery & Psychiatry 42 (1979) 1008–1017. [12] R. Beck, P. Cleary, M. Anderson, J. Keltner, W. Shults, D. Kaufman, E. Buckley, J. Corbett, M. Kupersmith, N. Miller, et al., A randomized, controlled trial of corticosteroids in the treatment of acute optic neuritis. The optic neuritis study group, New England Journal of Medicine 326 (1992) 581–588. [13] Y. Chen, H. Wyatt, W. Swanson, M. Dul, Rapid pupil-based assessment of glaucomatous damage, Optometry & Vision Science 85 (2008) 471–481. [14] L. Kalaboukhova, V. Fridhammar, B. Lindblom, Relative afferent pupillary defect in glaucoma: a pupillometric study, Acta Ophthalmologica Scandinavica 85 (2007) 519–525.
