A Simple and Accurate Model to Predict Responses to Multi-electrode Stimulation in the Retina
Implantable electrode arrays are widely used in clinical studies , clinical practice and basic neuroscience research and have advanced our understanding of the nervous system ., Implantable electronic devices can be used to record neurological signals and stimulate the nervous system to restore lost functions ., Sensing electrodes have been used in applications such as brain-machine interfaces 1 and localization of seizure foci in epilepsy 2 ., Stimulating electrodes have been used for the restoration of hearing 3 , sight 4 , 5 , bowel control 6 , and balance 7 , and in deep brain stimulation ( DBS ) to treat a range of conditions 8 ., Most neuroprostheses operate in an open-loop fashion; they require psychophysics to tune stimulation parameters ., However , devices that can combine both sensing and stimulation are desirable for the development of a new generation of neuroprostheses that are controlled by neural feedback ., Feedback in neuroprostheses is being explored in applications such as DBS for the enhancement of memory 9 , abatement of seizures 10 , control of Parkinson’s disease 11 , and the control of brain machine interfaces 12 ., Models that can accurately characterize a neural system and predict responses to electrical stimulation are beneficial to the development of improved stimulation strategies that exploit neural feedback ., Volume conductor models are typically used to describe retinal responses to electrical stimulation , however these are computationally intensive and can be difficult to fit to neural response data 13–15 ., Simpler models that can be constrained using neural recordings are necessary for real-time applications ., Linear-nonlinear models based on a spike-triggered average ( STA ) have been successfully used to characterize retinal responses to light 16–19 ., Models that incorporate higher dimensional components identified through a spike-triggered covariance ( STC ) analysis have been explored to describe higher order excitatory and suppressive features of the visual system 20–25 ., Generally , STA and STC models make use of white noise inputs and have the advantage that a wide repertoire of possible inputs patterns can be explored ., White noise models have previously been explored to describe the temporal properties of electrical stimulation in the retina 26 , 27 ., Spatial interactions between stimulating electrodes has not been previously investigated ., An example of a stimulation algorithm that could benefit from an accurate description of the spatial interactions is current steering , which attempts to improve the resolution of a device by combining stimulation across many electrodes to target neurons at a particular point 28 ., Two benefits obtained by using neural feedback algorithms are ( 1 ) the accurate prediction of the response to an arbitrary stimulus across the electrode array and ( 2 ) the ability to fit the device to individual patients from the recorded neural responses to a set of stimuli presented in a reasonable amount of time ., Here , we combine whole cell patch clamp recordings from individual retinal ganglion cells ( RGCs ) with stimulation using a multi-electrode array to demonstrate a model with the above advantages ., We find that a simple linear-nonlinear model accurately captures the effects of multi-electrode interactions and estimates the spatial relationship between stimulus and response ., The approach is scalable to a large number of electrodes , which is prohibitive to accomplish with psychophysics ., In contrast to conventional volume conductor models of electrical stimulation 13–15 , our model is straightforward to constrain using neural response data and is orders of magnitude more computationally efficient , making it suitable for use in real-time applications ., Methods conformed to the policies of the National Health and Medical Research Council of Australia and were approved by the Animal Experimentation Ethics Committee of the University of Melbourne ( Approval Number: 1112084 ) ., Data were acquired from retinae of Long-Evans rats ranging from 1 to 6 months of age ., Long-Evans rats were chosen for several reasons ., First , rat RGC morphological types have been examined in detail 29 , 30 and have similarities to RGCs found in other species , including the macaque monkey 31 and cat 32 , 33 ., Second , the size of the rat retina is larger than the mouse retina and so we were able to cover the entire stimulating electrode array with half of the retina ., The animals were initially anesthetized with a mixture of ketamine and xylazine prior to enucleation ., After enucleation , the rats were euthanized with an overdose of pentobarbital sodium ( 350 mg intracardiac ) ., Dissections were carried out in dim light conditions to avoid bleaching the photoreceptors ., After hemisecting the eyes behind the ora serrata , the vitreous body was removed and each retina was cut into two pieces ., The retinae were left in a perfusion dish with carbogenated Ames medium ( Sigma ) at room temperature until used ., Pieces of retina were mounted on a multi-electrode array ( MEA ) with ganglion cell layer up and were held in place with a perfusion chamber and stainless steel harp fitted with Lycra threads ( Warner Instruments ) ( Fig 1A ) ., Once mounted in the chamber , the retina was perfused ( 4–6 mL/min ) with carbogenated Ames medium ( Sigma-Aldrich , St . Louis , MO ) at room temperature ., The chamber was mounted on the stage of an upright microscope ( Olympus Fluoview FV1200 ) equipped with a x40 water immersion lens and visualized with infrared optics on a monitor with x4 additional magnification ., Electrical stimulation was applied subretinally through a custom-made MEA fabricated on a glass coverslip consisting of 20 platinum stimulating electrodes ( Fig 1 ) ., Each electrode had an exposed disc of 400 μm diameter , and a vertical pitch of 1 mm ., The stimulating area of the MEA spanned an area of 3 . 5 x 3 . 5 mm2 ( excluding the outer ring which was not used ) ., Glass coverslips were cleaned in an oxygen plasma chamber for 20 minutes ( Fig 1B ) ., Next , a positive ( UV-removable ) photoresist ( AZ1415H , Microchemicals ) was spin-coated onto the surface at 4000 revolutions per minute for 60 seconds ( Fig 1C ) ., A laser-printed chrome on soda glass photolithography mask was used to expose a pattern in the photoresist , then developed chemically ( MIF726 , Microchemicals ) revealing openings for electrode pads and tracks ( Fig 1D and 1E ) ., The developed cover slips were loaded into an electron beam deposition chamber ( Thermionics ) and pumped to a vacuum pressure of 1 . 5×10−6 mbar ., A 20 nm titanium adhesion layer was deposited at a rate of 0 . 2 Å/sec , followed by a platinum deposition of 130 nm at a rate of 0 . 6 Å/sec ., Residual photoresist was removed by soaking in acetone for 30 minutes , rinsing with deionized water , and finally oxygen plasma cleaning for 10 minutes ., For electrode isolation , a negative ( UV-curing ) photoresist ( SU8-2002 , Microchemicals ) was spin-coated onto the coverslip and exposed through a different photolithography mask leaving only metal exposed for stimulating electrodes and contact pads ( Fig 1F and 1G ) ., The entire device was then cured at 200°C on a temperature-controlled hotplate ., Whole cell intracellular recordings were obtained using standard procedures 34 at room temperature ., The main reason for recording at room temperature was to ensure that recordings lasted for many hours ., To obtain a whole cell recording , a small hole was made in the inner limiting membrane to expose a small number of RGC somas ., A pipette was then filled with internal solution containing ( in mM ) 115 K-gluconate , 5 KCl , 5 EGTA , 10 HEPES , 2 Na-ATP , 0 . 25 Na-GTP ( mosM = 273 , pH = 7 . 3 ) , Alexa Hydrazide 488 ( 250 mM ) , and biocytin ( 0 . 5% ) ( Fig 1A ) ., Initial pipette resistance in the bath ranged between 5–10 MΩ ., Prior to recording , the pipette voltage was nulled , pipette resistance was compensated with the bridge balancing circuit of the amplifier , and capacitance was compensated on the amplifier head stage ., Voltage recordings were collected in current clamp mode and amplified ( SEC-05X , NPI electronic ) , digitized with 16-bit precision at 25 kHz ( National Instruments BNC-2090 ) , and stored for offline analysis ., Intracellular recordings lasting up to 4 hours were obtained ., Stimulation artefacts that were present in the intracellular recording were removed offline by setting the membrane potential to a constant value for the duration of the stimulus ., Spikes in the remaining membrane potential waveform could be easily detected by finding peaks that crossed a set value ., Spike times were calculated as the time that the action potential reached its peak value ., Spike delay times were calculated by taking the difference between the spike time and the preceding stimulus offset time ., Intrinsic physiological differences , such as spike width , membrane time constant , and input resistance , among RGC types have been described 35 , 36 , which could lead to differences in response latencies to electrical stimulation ., Therefore , we performed a k-means cluster analysis on the spike latency from stimulation offset time ., The number of clusters ( k ) to fit was set manually by visual inspection of the clusters ., From the cluster analysis , we could detect if there were two or more clusters that might be attributed to direct activation or indirect activation via activation of presynaptic neurons ., Unless otherwise stated , responses to electrical stimulation were evaluated by analyzing the short-latency responses ., Short-latency responses were spikes that fell within two standard deviations of the mean of the shortest-latency cluster ., Long-latency responses fell within the cluster that occurred directly after the short-latency response ., Stimulation consisted of biphasic pulses of equal phase duration ( 500 μs ) with an interphase gap ( 50 μs ) and random amplitude ., The random amplitudes were sampled from a Gaussian distribution with variance σ2 ., Fig 2 illustrates the random amplitude pulses applied to all electrodes ., Stimulation waveform signals were generated by a custom-made MATLAB ( MathWorks version 2014a ) interface commanding a multi-channel stimulator ( Tucker Davis Technologies: RZ2 base station and IZ2 multichannel stimulator ) ., All stimulus amplitudes were bounded by the limits of the stimulator ( ±300 μA ) ., Biphasic pulses were applied to all electrodes at a frequency of 10 Hz and the numbers of short-latency responses were recorded ., To avoid overstimulation of a cell , an appropriate value of σ was chosen for each cell ., Three stimulus trains of 500 pulses were initially generated with fixed σ = 50 μA and applied to the tissue ., Next , a new set of stimulus trains were generated using a σ that varied between 50 μA and 250 μA in steps of 50 μA ., The number of spikes detected within 5 ms from the stimulus time was used to compute a response probability ., A sigmoidal curve was fit to the data of σ versus response probability to find the value of σ that resulted in half the maximum level of response ., For cells where the maximum response probability was close to 1 , σ was chosen to be a value that resulted in a response probability of 0 . 5 ., For other cells that saturated at a response probability less than 1 , σ was a lower value ., Once an appropriate value for σ was chosen for the cell , a stimulus vector ,, St→, , of length 20 ( equal to the number of electrodes ) was generated by sampling each element from a Gaussian distribution ., If the amplitude of stimulation of an electrode exceeded the stimulator limits ( ±300 μA ) , then the amplitude value was discarded and a new value was generated from the distribution ., Each stimulus was applied 3–5 times before a new, St→, was generated ., The experiment continued for as long as the cell’s responses remained stable ( usually 3–4 hours ) ., Once cells started to show signs of deterioration ( e . g . unstable high frequency spontaneous activity ) , the experiment was terminated ., After recording , the retinal tissue was removed from the chamber and mounted onto filter paper ., The tissue was then fixed for ~45 min in phosphate-buffered 4% paraformaldehyde and stored for up to 2 weeks in 0 . 1 M phosphate-buffered saline ( PBS; pH 7 . 4 ) at 4°C ., The tissue was then immersed in 0 . 5% Triton X-100 ( 20 μg/ml streptavidin conjugated Alexa 488; Invitrogen ) in PBS overnight to expose biocytin-filled cells ., Tissue was washed thoroughly in PBS , mounted onto Superfrost+ slides , and coverslipped in 60% glycerol ., Cells were then reconstructed in 3D using a confocal microscope ( FluoView FV1200 ) ., RGC types were initially identified by their focal light response at the beginning of each experiment ., ON cells showed an increase in spike rate at the onset of light; OFF cells showed an increase in spike rate at the offset of light; ON-OFF cells showed an increase in spike rate at the onset or offset of light ., Additionally , RGC classification was correlated with morphology based on dendritic stratification in the inner plexiform layer ( IPL ) 29 , 30 ., The level of stratification was defined as 0–100% from the level of the inner nuclear layer to the level of the ganglion cell layer ., The stratification depth ( s ( x ) ) was quantified as a percentage of the inner plexiform layer ( IPL ) thickness , according to, s ( x ) =100 ( Ls−xLs−Le ) ., ( 1 ), Here , x refers to the depth of a terminal dendrite and Ls and Le represent the IPL-GCL border and the GCL-INL border of the inner plexiform layer , respectively , where depth decreases from the ganglion cell layer towards the photoreceptor layer ., Cells that stratified in the inner part of the IPL ( s ( x ) ≤ 40% ) are denoted as OFF-cells ., Cells that stratified in the outer part of the IPL ( s ( x ) ≥ 60% ) are referred to as ON-cells ., For all cells in this study , the physiological and morphological classifications correlated well ., Dendritic field sizes were calculated by tracing out a region connecting the dendritic tips of each cell and fitting an ellipse to the region ., The major axis of the ellipse was used to estimate the dendritic field size ., Our objective was to find a mathematical description able to accurately capture the spiking probability of RGCs to subretinal stimulation using a MEA ., We characterized neural responses using an N-dimensional linear subspace of the stimulus space , combined with a nonlinearity describing the intrinsic nonlinear firing properties of neurons ., Using STC analysis , we derived the lower dimension stimulation subspace that led to a short-latency response in the neuron ., By projecting the raw and spike-triggered stimuli onto the lower dimension subspace , we estimated the intrinsic nonlinearity ., The probability of generating a spike , given stimulus, St→, , was estimated as, P ( R=spike|St→ ) =NN ( v→1∙St→ , v→2∙St→ , … , v→N∙St→ ) ,, ( 2 ), where, N, represents the static nonlinear function operating on arguments in μA and, v→i, ( i = 1 , 2 , … , N ) represent the N significant principal components ., To find, v→i, ( i = 1 , 2 , … , N ) , the stimulus data were first separated into a matrix containing only stimuli generating a short-latency response , SD , and a matrix containing all stimuli , ST ( Fig 3A ) ., We found the low-dimensional linear subspace that best captured the difference between the spike-triggered stimuli and the raw ensemble by performing principal component analysis ( PCA ) on the covariance matrix of the spike-triggered ensemble ,, Cs=cov ( SD ) ,, ( 3 ), and comparing it to the variance of the raw ensemble which was approximately σ2 in all stimulus directions due to the Gaussian nature of ST . PCA on Cs produce a set of eigenvectors giving a rotated set of axes in stimulus space and a corresponding set of eigenvalues giving the variance of the spike-triggered ensemble along each of the axes ., Eigenvalues that are greater than the variance of the input stimuli reveal the directions where the spike-triggered stimuli have experienced an increase in variance from the raw ensemble ., Similarly , eigenvectors that are smaller than the variance of the input stimuli reveal directions where the spike-triggered stimuli have experienced a decrease in variance from the raw ensemble ., The eigenvalues that are sufficiently different from the raw ensemble correspond to eigenvectors (, v→i , i=1 , 2 , … , N, ) pointing in directions in the stimulus space that carry information about the spiking probability of the neuron ., To test if the neural response could be accurately characterized by a one-dimensional model , we examined how many eigenvalues resulting from PCA were significantly different to chance 20 ., We compared the eigenvalues obtained by PCA on spike-triggering stimuli to a distribution of eigenvalues for a randomly chosen ensemble of stimuli ., This was done by randomly time-shifting the spike train and performing PCA on the corresponding randomized spike-triggered stimuli to recover a new set of eigenvalues ., By repeating these 1000 times , we construct a distribution of eigenvalues and set a confidence criterion outside of which we presumed the magnitude of the true eigenvalues to be significant ., The confidence criterion used was two standard deviations , or a 95% confidence interval ., If the greatest or least eigenvalue fell outside these bounds , we rejected the null hypothesis that the spike-triggered ensemble was not significantly different to the full ensemble ., We then projected out the axis corresponding to this eigenvalue from the raw data ., We repeated the test until all remaining eigenvalues fell within the bounds of the null distribution , suggesting that the remaining eigenvalues were insignificant in affecting the variance of the neuron ., Components having an eigenvalue significantly greater than the variance of the randomly time-shifted ensemble were considered to be components that generate an excitatory response on the cell ., Conversely , components that are significantly smaller than the variance of the raw ensemble were considered to be components that suppressed the cell’s response ., For the majority of cells , we found that a simplification to one-dimension (, v→1, ) accurately captured the spike-triggering information , thereby reducing the equation to one dimension ., Using this simplification , Eq ( 2 ) becomes, P ( R=spike|St→ ) =N1 ( v→1⋅St→ ) ., ( 4 ), Results in the literature indicate that response thresholds to electrical stimulation for some cell types might differ depending on pulse polarity 37 ., To explore difference in response to pulse polarity , we allowed the probability to be described by two different nonlinear functions and we found the electrical receptive fields ( ERFs ) for stimuli having a net effect that was either cathodic-first or anodic-first ., Eq ( 4 ) then becomes, P ( R=spike|St→ ) =N+ ( w→+⋅St→ ) +N− ( w→−⋅St→ ) +cRs ,, ( 5 ), where, N+, and, N−, represent static nonlinear functions and, w→+, and, w→−, represent the ERFs for stimuli with positive projections (, v→1∙St→>0, , net anodic-first ) and negative (, v→1∙St→<0, , net cathodic-first ) , respectively ., Rs represents the spontaneous firing rate in Hz and c represents a scaling factor of units Hz-1 ., To find the nonlinearities and the ERFs , the first principal component (, v→1, ) was used to divide the stimulus space into positive and negative regions by projecting all stimuli of SD and ST onto the first principal component ( Fig 3B ) ., Positive and negative regions were defined by the stimuli having either a positive or negative projection onto the first principal component ., This produced two spike-triggered stimulus matrices ,, SD+, and, SD−, ., The means of the matrices are analogous to the spike-triggered average for net anodic-first and net cathodic-first stimuli 16 , and provide an estimate of the ERFs ,, w→+, and, w→−, , respectively ., Fig 3B shows an example of the stimuli projected onto the first two principal components and the ERFs ,, w→+, and, w→−, ., After the stimuli were separated into two regions , the nonlinear functions ,, N+, and, N−, , were recovered by projecting all stimuli onto the normalized vectors, w→+, and, w→−, and segmenting the projected stimuli into 30 bins ( 15 for the net anodic-first and 15 for the net cathodic-first regions ) such that each bin contained an equal number of spikes ., By comparing the number of spike-eliciting stimuli to the total number of stimuli in each bin , an estimate of the spike probability was generated ., The nonlinear function from Eq ( 5 ) was then fit to the data , with the following equations for the sigmoidal curves:, N+ ( x+ ) =a+1+exp ( −b+ ( x+−c+ ) ), ( 6 ), N− ( x− ) =a−−a−1+exp ( −b− ( x−−c− ) ) ,, ( 7 ), where, x+=w→+∙St→, and, x−=w→−∙St→, ., Coefficients a+ and a− represent scaling factors that determine the saturation amplitudes , b+ and b− represent the gain of the sigmoidal curves , and c+ and c− represent the thresholds ( 50% of the saturation level ) for the net anodic-first and net cathodic-first stimulation , respectively ., Note that the vectors, w→+, and, w→−, might not necessarily be parallel to each other , nor parallel to, v→1, ., This may result in electrodes that differentially influence the neuron’s response to anodic-first or cathodic-first stimulation ., To test the similarity between the positive and negative ERFs , we calculated the correlation coefficient between them ., A correlation coefficient close to -1 indicated that the two ERFs are approximately equal in magnitude but opposite in sign , and therefore the cell was equally influenced by the same combination of electrodes ., A value closer to 0 indicates that the two ERFs have no correlation , and therefore the cell is not influenced by the same electrodes ., Positive correlation coefficients were not expected and did not occur ., The spatial extent over which a cell was influenced by electrical stimulation was estimated by computing a weighted mean of the distance between the cell and the electrodes ., The distance between the cell and each electrode center was weighted by the electrode’s influence on the cell as given by the ERFs ., The weighted mean for both ERFs was given by ,, D+=∑i=120wi+di∑i=120wi+, ( 8 ), D−=∑i=120wi−di∑i=120wi−, ( 9 ), where, wi+, and, wi−, are the weights given by, w→+, and, w→−, respectively , and di is the distance between the cell and electrode i ., To test which electrodes significantly affected the cells response ,, w→+, and, w→−, were recalculated 1000 times by projecting the data onto the first eigenvector of the randomly time-shifted distribution of eigenvectors from the significance test ., A distribution for, w→+, and, w→−, was constructed from which the true, w→+, and, w→−, could be compared ., Electrodes from the true, w→+, and, w→−, were compared to the root mean square ( RMS ) of the distribution and electrodes that were larger than the RMS bounds were considered significant ., For cells where more than one significant principal component was obtained from the significance test , we compared the variance explained by the first principal component to that of the next most significant component ., This was done by comparing the separation of the first eigenvalue e1 from the mean of the randomized distribution of eigenvalues (, e¯rnd, ) with the separation between the next most significant eigenvalue ( e2 ) and the same mean ., The strength was defined as, G=|e1−e¯rnd||e2−e¯rnd| ,, ( 10 ), and gives a relative measure of how much larger e1 is compared to the next most significant eigenvalue ., e¯rnd, was calculated from the first iteration of the hypothesis test ., For each cell , 80% of the data were used to fit the model parameters , while the remaining data were used to validate the model ., To obtain a quantitative estimate of the performance of the model , the probability of response given a stimulus was calculated from the validation data and compared to the model prediction ., The validation stimuli were assigned 1 if they produced a direct response and 0 otherwise ., Each stimulus was also assigned a predicted probability using the model ( Eq ( 5 ) ) recovered from the training data ., The stimuli were then binned into segments in the range of 0 to 1 depending on their predicted probability and an actual probability for each bin was calculated by the fraction of stimuli assigned a 1 ., The mean square error ( EMS ) was then calculated ,, EMS=1B∑i=1B ( P^i−Pi ) 2 ,, ( 11 ), where B is the number of bins ,, P^i, is the predicted probability , and Pi is the calculated probability from the data for a particular bin ., For all cells , B was equal to 10 ., The root mean square error ( ERMS ) of the model , given by, ERMS=EMS ,, ( 12 ), was used as a quantitative measure of the model accuracy ., We also compared the error of a one-dimensional model to that of a two-dimensional model ., The two-dimensional spike probability was estimated by, P ( R=spike|St→ ) =N2 ( v→1⋅St→ , v→2⋅St→ ) ,, ( 13 ), where, v→2, represented the next most significant component , either the second ( excitatory ) or last ( suppressive ) principal component ., To find the two-dimensional nonlinearity (, N2, ) , a surface was fit to the spike-triggered data projected onto these two most significant components ., The surface fit was obtained using a cubic spline interpolation on MATLAB’s curve fitting toolbox ., Once the surface was fit , the validation data was used to calculate the mean model error calculated using Eqs ( 11 ) and ( 12 ) ., Intracellular recordings lasting up to 4 hours were obtained from 25 cells ., This population included 7 ON , 13 OFF , 3 ON-OFF , and 2 cells where 3-D morphological reconstructions were not possible ., Our comparison of histological and physiological results were consistent with those of Huxlin and Goodchild 29: ON center cells stratify in the inner IPL ( 40–100% depth ) , while OFF center cells stratify in the outer IPL ( 0–40% depth ) ., ON-OFF types stratify in both the inner and outer layers of the IPL ., Fig 4 shows an example of an ON-OFF RGC with dendrites stratifying in both inner and outer layers of the IPL ., A summary of the stratification depths for the ON , OFF , and ON-OFF cells are given in Table 1 ., To fit the model parameters , cells were simultaneously stimulated with biphasic pulses on all electrodes , where the amplitude of the pulses were randomly chosen from a Gaussian distribution of zero mean and standard deviation σ ( here after white noise stimuli ) ., To determine an appropriate value of σ for each cell , three short stimulus trains ( approximately 3 min each ) of white noise stimuli with different σ were initially presented to the cell ( σ varied from 50–250 μA in steps of 50 μA ) ., The number of times the cell responded within 5 ms was used to obtain a response probability ., Each cell responded with a different maximum response probability when stimulated with white noise at the highest value of σ; some cells could respond with a spike probability close to one , while others only responded with a spike probability less than one ., However , cells that responded to fewer pulses tended to show an increased level of long-latency activity ( > 5 ms ) , most likely due to intensified network activation ., The value of σ used for white noise stimulation for each cell in the rest of the experiment was the value corresponding to half the saturation level ., Fig 5 shows examples of two cells with different σ values ., Cell 2 responded with a spike probability close to one even at low σ values while cell 1 responded maximally with a spike probability of around 0 . 6 ( Fig 5A ) ., The value of σ used for white noise stimulation for cell 1 was 85 μA and for cell 2 was 145 μA ., Note that we used this method to calibrate our experiments and the nonlinear curves do not show the maximum probability of firing , as each point is an average over a variety of stimulus amplitudes ., Following this calibration , longer trains of white noise stimulation ( approximately 2 minutes each ) with the corresponding value of σ for each cell were used to obtain data for recovering the model parameters ., The corresponding Gaussian distributions for cells 1 and 2 are shown in Fig 5B ., The experiment for each cell lasted approximately 3–4 hours ., Stimulation artefacts were present in the recordings that could be removed by blanking without affecting the ability to detect the cells’ spikes ., Fig 6A shows examples of some of the spiking patterns observed during experiments:, ( i ) a failed anodic-first stimulus ,, ( ii ) a successful short-latency anodic-first stimulus ,, ( iii ) a successful short-latency cathodic-first stimulus , and, ( iv ) a successful long-latency cathodic-first stimulus ., The top panel in each subplot shows the raw recording and the bottom panels show the same signals with the artefact removed by blanking ., Also shown in the bottom panels are the thresholds used to detect spikes ( horizontal lines ) ., These figures show that spikes could be easily identified without interference from the stimulus artefact ., The spike latencies after a stimulus pulse were analyzed for each cell ., Some cells produced a bimodal distribution attributed to the short- and long-latency responses ( N = 13 ) , with four cells showing overlapping distributions for the two latencies ., The remaining cells only produced short-latency spikes that were close to the timing of the stimulus pulse ( N = 8 ) ., Fig 6B depicts the spike latencies for all cells ., The average short-latency cluster mean for all cells was 1 . 75 ms from stimulus offset ( SD 1 ms ) ., The longest short-latency cluster for a cell had a mean of 4 . 35 ms ( SD 1 . 37 ms ) ., Fig 6C and 6D show the distributions of spike latencies for two sample cells , along with fitted Gaussian distributions obtained from the cluster analyses ., Fig 6C shows a cell with two distinct clusters , with a short-latency cluster mean at 1 . 95 ms . Fig 6D shows two overlapped clusters with the short-latency cluster mean at 4 . 35 ms . Our aim was to find a mathematical description that could accurately capture the response probability of neurons to concurrent stimulation using a MEA ., To do this we first performed a principal components analysis on the ensemble of stimuli that triggered a short latency spike ., For all cells we found that the neural response could be well predicted by projection onto a subspace spanned by the first principal component ,, v→1, ., The variance explained by, v→1, was significantly higher than that of next greatest component ,, v→2, , suggesting that the spiking information was well captured by, v→1, ., Fig 7A shows the spike-triggered probabilities projected onto, v→1, and, v→2, from the sample cell in Fig 3B ., The histograms show the number of stimuli ( gray ) and responses ( black ) along each axis; the ratio of the bars of the two histograms is used to determine the spike probabilities along each axis ., From the histograms , it is clear that the distribution of the spike-triggered stimuli was bimodal in the, v→1, axis; however , it remained unimodal along, v→2, , similar to Gaussian distribution of the full stimulus ensemble ., A statistical hypothesis test was used to determine how many eigenvalues recovered by PCA revealed a significant amount of the spike-eliciting information ., The test compares the eigenvalues recovered from the data , to a set of eigenvalues produced by randomly time-shifting the spike train and performing PCA on the new set of stimuli ., From the set of time-shifted eigenvalues , a 95% confidence limit was set to determine
Implantable electrode arrays are widely used in therapeutic stimulation of the nervous system ( e . g . cochlear , retinal , and cortical implants ) ., Currently , most neural prostheses use serial stimulation ( i . e . one electrode at a time ) despite this severely limiting the repertoire of stimuli that can be applied ., Methods to reliably predict the outcome of multi-electrode stimulation have not been available ., Here , we demonstrate that a linear-nonlinear model accurately predicts neural responses to arbitrary patterns of stimulation using in vitro recordings from single retinal ganglion cells ( RGCs ) stimulated with a subretinal multi-electrode array ., In the model , the stimulus is projected onto a low-dimensional subspace and then undergoes a nonlinear transformation to produce an estimate of spiking probability ., The low-dimensional subspace is estimated using principal components analysis , which gives the neuron’s electrical receptive field ( ERF ) , i . e . the electrodes to which the neuron is most sensitive ., Our model suggests that stimulation proportional to the ERF yields a higher efficacy given a fixed amount of power when compared to equal amplitude stimulation on up to three electrodes ., We find that the model captures the responses of all the cells recorded in the study , suggesting that it will generalize to most cell types in the retina ., The model is computationally efficient to evaluate and , therefore , appropriate for future real-time applications including stimulation strategies that make use of recorded neural activity to improve the stimulation strategy .
Implantable multi-electrode arrays ( MEAs ) are used to record neurological signals and stimulate the nervous system to restore lost function ( e . g . cochlear implants ) ., MEAs that can combine both sensing and stimulation will revolutionize the development of the next generation of devices ., Simple models that can accurately characterize neural responses to electrical stimulation are necessary for the development of future neuroprostheses controlled by neural feedback ., We demonstrate a model that accurately predicts neural responses to concurrent stimulation across multiple electrodes ., The model is simple to evaluate , making it an appropriate model for use with neural feedback ., The methods described are applicable to a wide range of neural prostheses , thus greatly assisting future device development .
medicine and health sciences, action potentials, engineering and technology, signal processing, membrane potential, ocular anatomy, retinal ganglion cells, electrophysiology, neuroscience, surgical and invasive medical procedures, mathematics, functional electrical stimulation, ganglion cells, algebra, white noise, membrane electrophysiology, bioassays and physiological analysis, research and analysis methods, animal cells, electrophysiological techniques, cellular neuroscience, retina, electrode recording, linear algebra, cell biology, anatomy, physiology, neurons, biology and life sciences, cellular types, physical sciences, afferent neurons, ocular system, eigenvalues, neurophysiology