A cascade model of information processing and encoding for retinal prosthesis

Abstract

Retinal prosthesis offers a potential treatment for individuals suffering from photoreceptor degeneration diseases. Establishing biological retinal models and simulating how the biological retina convert incoming light signal into spike trains that can be properly decoded by the brain is a key issue. Some retinal models have been presented, ranking from structural models inspired by the layered architecture to functional models originated from a set of specific physiological phenomena. However, Most of these focus on stimulus image compression, edge detection and reconstruction, but do not generate spike trains corresponding to visual image. In this study, based on state-of-the-art retinal physiological mechanism, including effective visual information extraction, static nonlinear rectification of biological systems and neurons Poisson coding, a cascade model of the retina including the out plexiform layer for information processing and the inner plexiform layer for information encoding was brought forward, which integrates both anatomic connections and functional computations of retina. Using MATLAB software, spike trains corresponding to stimulus image were numerically computed by four steps: linear spatiotemporal filtering, static nonlinear rectification, radial sampling and then Poisson spike generation. The simulated results suggested that such a cascade model could recreate visual information processing and encoding functionalities of the retina, which is helpful in developing artificial retina for the retinally blind.

Introduction

As neural entrance of human visual system, the retina contains five different cell types: ganglion, amacrine, horizontal, bipolar cells and photoreceptor. These visual cells are organized into a layered architecture. There are many subclasses within each major cell type (Goetz and Trimarchi, 2012). Not the same as in the brain, most neurons in the retina are much localized. They respond through gradual changes of the membrane potential. The only exception to this is ganglion cells, whose long axons converge to make up the optic nerve (Siegert et al., 2009).

Visual information from the surrounding environment is captured by photoreceptors and subsequently converted into electrical neural signals, which then are processed by different excitatory and inhibitory pathways coexisting inside the retinal network and ultimately the output of the retina to higher brain areas by ganglion cells is encoded as spike trains (Seung and Sümbül, 2014).

Nowadays, approximately 1 in 4 blind individuals worldwide suffers from selective photoreceptor loss, which often arises from blinding retinal diseases, such as retinitis pigmentosa or age-related macular degeneration. However, no effective treatment currently exists for retinitis pigmentosa, while progression of age-related macular degeneration can be slown but not prevented. In this case, bypass the lesioned neural tissue, using electrical stimulation to artificially deliver visual information to the surviving retina has become increasingly accepted as the best near-term hope for restoring useful vision in these patients, such as the use of retinal prosthesis (Weiland et al., 2011, 2014).

According to the location of their device, retinal prosthesis generally falls into two categorical types: epiretinal and subretinal prosthesis. No matter what kind of retinal prosthesis to be practical, an inevitable question that must be solved is designing a device that can perform necessary biological function as efficiently as retina. Hence, modeling and simulating the biological retina how to transform incoming light signal into spike trains that can be properly decoded by the brain is a challenge for the implementation of retinal prosthesis.

Amongst the existing models of retina, some structural models built on a full set of synaptic and cellular parameters, they concern more with a detailed reproduction of anatomic connections in successive layers (Bálya et al., 2002; Bonin et al., 2005; Passaglia et al., 2009; Lefebvre et al., 2011; Lideberg, 2013). Other functional models place special emphasis on recreating functional computations of retina, each physiologic phenomenon being modeled with a series of filters (Pelayoa et al., 2004; Wohrer and Kornprobst, 2009; Adesnik et al., 2012). Most of these only focus on stimulus image compression, edge detection and reconstruction, but few studies have given models that generate spike trains corresponding to input visual information.

In the present study, a cascade model of the retina including the out plexiform layer (OPL) for information processing and the inner plexiform layer (IPL) for information encoding was presented, which integrates both anatomic connections and functional computations of retina. In this model, spike trains corresponding to stimulus image were numerically computed by linear spatiotemporal filtering, static nonlinear rectification, radial sampling and then Poisson spike generation based on MATLAB software. Subsequently, we give a general discussion of the advantages and disadvantages of the cascade model.

Materials and Methods

The cascade model

[Figure 1] exhibits the global structure of our model comprised of two cascade synaptic layers: the OPL and the I. The OPL defines a signal processing stage through the interaction of light receptors, horizontal cells and bipolar cells. It is modeled as a linear spatiotemporal filter to mediate spatial oppositions that allow the retinal output to enhance image edges. The IPL deals with signal encoding through the interplay of bipolar, amacrine and ganglion cells, which produces the spike rate and then the sampled ganglion cells generating spike trains. Furthermore, information stream flowing between neighboring cells from different layers is depicted by arrows.{Figure 1}

Computation by the OPL

In the OPL, the visual signal from receptors to bipolar cells can take two parallel paths. The direct path from receptors to bipolar cells results in the excitatory center signal. The indirect path transits through horizontal cells, and hence through two synapses, this leads to the inhibitory surround signal. In our model, synaptic transmissions are considered linear and both paths will correspond to a linear spatiotemporal filter. The spatial filtering is produced by Gaussian kernel as follows:

here A and σ are, respectively, the peak sensitivities and standard deviations of gaussian function. Their temporal component is the convolution of several exponential kernels as follows:

[INLINE:2]

where α denotes the total number of exponential kernel, and τ is the decay constant. For t = ατ, the resulting filter peaks. When applied to a signal, the kernel induces a temporal delay of value (α + 1)τ.

Finally, always antagonistic center and surround signals are mediated by the bipolar cells layer (Tokutake and Freed, 2008; Cho and Choi, 2014). The band-pass filtering properties of bipolar cells layer are given by the following equations:

[INLINE:3]

[INLINE:4]

[INLINE:5]

where I(x, y, t) is the input luminosity visual information; ω represents the relative weight between center and surround signal. According to the previous physiological experiments recording, A C > A S , σC < σS, d C <<i> d S.

Computation by the IPL

In the IPL, ganglion cells are thought as completely independent spiking encoders. Each ganglion cell responds to the stimulus, by modulating its firing frequency, within its receptive field and then transmits that information to the next visual stage independent of other ganglion cells.

In contrast to the OPL, signal transmission in the IPL is approximately static nonlinear and is implicit in most of the experimental literature (Naka and Machuca, 1995). Biologically, static nonlinearities in the retinal neuron network arise from multiple different reasons, such as spiking threshold and saturation. A simplified synaptic transmission function is applied to approximate this translation process (Keat et al., 2001). This step is expressed by the following equations:

[INLINE:6]

[INLINE:7]

[INLINE:8]

here VBip is normalization processing of VBip-1 ; V Gan is the average spiking rate of the ganglion cells; VMax is the maximum of the spiking rate; V0 is spiking threshold; λ is static nonlinear coefficient; and µ is a constant, 0 ≤ µ ≤ 1.

In an effort to make the pixels of spike rate image matching with the ganglion cells for one-to-one correspondence, the distribution mode of ganglion cells must be considered. In a human retina, ganglion cells are not uniformly distributed, but show a particular radial structure, whose density decreases exponentially with the distance increase from its central region, called the fovea (Curcio and Allen, 1990; Harman et al., 2000; Masland, 2012). Hence, the spike rate image was sampled in the polar coordinate fashion given by the following equations:

[INLINE:9]

[INLINE:10]

here i is the index of sampling circular; j is the index of sampling point; R is sampling radius; θ is the sampling polar angle; M is the maximum of sampling radius.

Measuring Inter-Spike Intervals (ISI) has been extensively used for the computational studies of spike generation process of real neurons (Dayan and Abbott, 2005; Gollisch and Meister, 2008). For a more biological ISI distribution, the refractory period following each spike is taken into account, which ensures that two successive spikes can never get closer than a few milliseconds. We calculate the firing times of ganglion cells by:

[INLINE:11]

[INLINE:12]

If (11) again (12)

where VRand is a random number generated by MATLAB and set the bounds from zero to one; TRef is the absolute refractoriness of retinal ganglion cell, which obeys the normal distribution of N (3 ms, 1 ms) (Field and Chichilnisky, 2007).

Numeric computation was performed with MATLAB 2007b software (MATLAB, Natick, MA, USA).

Results

The visual stimulus

In this work, a static Lena grayscale image was implemented as original visual stimulus to our cascade model. The size of image was 256 × 256 pixels and gray-level was 256. The stimulus intensity was represented by pixel value.

Information processing in the OPL

When visual stimulus image was given to the model, the OPL processed it by using Gaussian kernel and exponential cascade kernel for spatiotemporal filtering.

The Gaussian function is a common way to model electrical coupling between cells and dendritic spread of retinal cells. [Figure 2]A shows the low-pass spatial filter, which is a normalized and averaging filter. The low-pass temporal filter is also normalized to have an integral of one, which only performs temporal averaging on incoming information, without a linear gain. When the parameter α is greater 0, the shape of exponential cascade filter shows high variability [Figure 2]B.{Figure 2}

In equation (5), ω defines the relative weight between surround and center signal. The best-fitting weights for center and surround single are biologically measured to display a strong cell-to-cell variability, which ranges from 0.85 to 1 (Enroth-Cugell and Robson, 1966).

The simulated result of information processing in the OPL with different relative surround weights are shown in [Figure 3]. As ω increases, the filtering characteristic of OPL is affected obviously. When ω = 0.95, spatiotemporal filtering exports clear edge and good contrast image. The intensity of each pixel in the difference image indicates the change in spiking rate for the ganglion cell, whose receptive field is centered at that pixel. The main parameters we used in the OPL are: AC = AS = 1, σC = 1.5 pix, αS = 4.5 pix, (αC + 1)τC = 30 ms, (αS + 1)τS = 60 ms.{Figure 3}

Information encoding in the IPL

Static nonlinear rectification has been demonstrated as a very distinct feature in visual image undergoes neural transformation. As shown in [Figure 4], the experimental curves reported early (Chichilnisky, 2001) can be fitted well with the simplified synaptic transmission function for static nonlinearity in the IPL. Visual neurons in the retina could effectively adjust their gain to maintain an optimal sensitivity at varying light levels. The sigmoid function was carefully adapted to assure a reasonable physiological range for the spike rate of ganglion cells. The main parameters we used in the static nonlinear rectification are: VMax = 70, λ = 0.05, V0 = 0.06.{Figure 4}

Processing the spatiotemporal filtering image resulting from the out plexiform layer with this function gives the average spike rate image, the result is shown in [Figure 5]A and what follows is sampling them according to the topography of ganglion cells in the retina. To this end, a uniform radial sampling pattern of all spike rates for one circular level and 25 samples per circle were adopted [Figure 5]B.{Figure 5}

The average firing rate of sampled ganglion cells is taken as the input parameter (VGan ) for Poisson spike generator. Once the absolute refractory is passed, a spike is emitted. Finally, as shown in [Figure 6], we get the firing times of 25 sampled ganglion cells within one second. Dark vertical lines indicate a neuron fires at certain time. This temporal sequence of spikes accords with the sparse encoding feature of neurons. When transforming a continuous firing rate into a discrete and random set of spike trains, the Poisson process is considered as a good alternative for coding, in terms of the statistical links between the spike train and the underlying firing rate (Kara et al., 2000; Hartveit and Heggelund, 1994). So we believe that the output is a sequence that indicates the spiking pattern we need to induce on the ganglion cell in order to achieve the desired visual response to our input image in brain.{Figure 6}

Discussion

The goal of this study was to construct a cascade model to simulate how the biological retina processes and encodes incoming light signal. The main constraint that has been found to be imposed by retinal connectivity is the nonseparability of filtering in the OPL. Because of its cellular implementation, the spatial surround signal is necessarily delayed temporally. In comparison with the conventional difference of Gaussian model, a supplementary temporal component is induced in this paper, which makes our non-separable spatiotemporal filter not only acting as an edge detector but also preserving more detailed information of the original image. Integration of center and surround signals at the level of bipolar cells implies a signal weighting parameter ω that concentrates the whole balance of a channel of unites between low pass and band pass properties. To indeed display a visible mixture of both characteristics, we found the best weight ω to be about 0.95, but there is no evidence that weight is a constant for different subtypes of bipolar cells, nor that it is a constant spatially. The static nonlinearity naturally provides a better fit to real cellular response and allows the model to function under the wide variety of illuminations in our environment. Even in real experiments, it is observed that retinal spike trains have more or less variability than a Poisson emission process, at least during the periods of high firing activity. In case of no other information is known about a spike train than its average firing rate, we applied the homogeneous Poisson process to minimize artificial, model-related correlations, in addition, neuronal absolute refractory period is embedded in the precise spike emissions of a Poisson process. It efficiently reflects the information processing and encoding mechanism of retina in response to a static visual stimulus. A further advantage is that the model allows for the simultaneous analysis of the responses of a large number of neurons for a given stimulus condition.

In order to reproduce the most important characteristics of retina, it was attempted to include all relevant details on the basis of the available experimental data. However, full biological realism is difficult to achieve in a model, as it is impossible to include every known fine detail of retinal anatomy and physiology. Therefore, specific assumptions and simplifications were made in order to ensure a high degree of realism. This is sometimes difficult as the available experimental data is not always conclusive or contradictory. Also, many details of the model are based on experimental evidence from different species. For example, all neurons were implemented as single-compartment structures, and effects of neighboring cells in the same layer were neglected, such as photo transduction and cellular adaption (Chen et al., 2005). Moreover, translation equation between distinct layers was simplified to allow for an easier mathematical treatment, and information loss is inevitable. Additionally, an important information processing property of the retina, contrast gain control (Victor, 1987) was not involved in our model. Empirical models have been proposed, some including divisive normalization of the signal by a surround component, but no definitive explanation about the usefulness or biological origin of this phenomenon exists yet. Temporal properties of ganglion cells as well as their global filtering gain depend on the local level of contrast. Empirical models have been proposed, some including divisive normalization of the signal by a surround component, but no definitive explanation about the usefulness or biological origin of this phenomenon exists yet. This effect could result in a nonlinear compression function that would be dynamic, depending on the local contrast of the scene (Berry et al., 1999; Jarsky et al., 2011; Scholl et al., 2012).

The overall motivation for the work is to restore partial vision to people who are blind due to loss of photoreceptor function. If we can construct a more complicated retina model including all basic retinal neuron classes and their realistic interconnection patterns based on better anatomical and physiological data, it will provide greater assistance in the realization of a clinically viable prosthesis.[32]