Unravelling the temporal properties of human eCAPs through an iterative deconvolution model

OBJECTIVE
The electrically evoked compound action potential (eCAP) has been widely studied for its clinical value in evaluating cochlear implants (CIs). However, to date, single-fiber recordings have not been recorded from the human auditory nerve, and many unknowns remain about the firing properties that underlie the eCAP in patients with CIs. In particular, the temporal properties of auditory nerve fiber firing might contain valuable information that may be used to estimate the condition of the surviving auditory nerve fibers. This study aimed to evaluate the temporal properties of neural firing underlying human eCAPs with a new deconvolution model.


DESIGN
Assuming that each auditory nerve fiber produces the same unitary response (UR), the eCAP can be seen as a convolution of a UR with a compound discharge latency distribution (CDLD). We developed an iterative deconvolution model that derived a two-component Gaussian CDLD and a UR from recorded eCAPs. The choices were based on a deconvolution fitting error minimization routine (DMR). The DMR iteratively minimized the error between the recorded human eCAPs and the eCAPs simulated by the convolution of a parameterised UR and CDLD model (instead of directly deconvolving recorded eCAPs). Our new deconvolution model included two separate steps. In step one, the underlying URs of all eCAPs were derived, and the average of these URs was called the human UR. In step two, the CDLD was obtained by using the DMR in combination with the estimated human UR. With this model, we investigated the temporal firing properties of eCAPs by analysing the CDLDs, including the amplitudes, widths, peak latencies, and areas of CDLDs. The differences of the temporal properties in eCAPs between children and adults were explored. Finally, we validated the two-Gaussian component CDLD model with a multiple-Gaussian component CDLD model.


RESULTS
The estimated human UR contained a sharper, narrower negative component and a wider positive phase, compared to the previously described guinea pig UR. Furthermore, the eCAPs from humans could be predicted by the convolution of the human UR with a two-Gaussian component CDLD. The areas under CDLD (AUCD) reflected the number of excited nerve fibers over time. Both the CDLD magnitudes and AUCDs were significantly correlated with the eCAP amplitudes. Furthermore, different eCAPs with the same amplitude could lead to greatly different AUCDs. Significant differences of the temporal properties of eCAPs between children and adults were found. At last, the two-Gaussian component CDLD model was validated as the most optimal CDLD model.


CONCLUSION
This study described an iterative method that deconvolved human eCAPs into CDLDs, under the assumption that auditory nerve fibers had the same electrically evoked UR. Based on human eCAPs, we found a human UR that was different from the guinea pig UR. Furthermore, we found that CDLD characteristics revealed age-related temporal differences between human eCAPs. This temporal information may contain valuable clinical information on the survival and function of auditory nerve fibers. In turn, the surviving nerve condition might have prognostic value for speech outcomes in patients with CIs.


Introduction
A cochlear implant (CI) is a device that restores hearing by directly applying electrical stimulation to the auditory nerve fibers inside the cochlea. A CI can also be used to record auditory nerve activity via a telemetry function; this recording yields the electrically evoked compound action potential (eCAP). The eCAP is an objective measure that can be used to assess the quality of the electrode-nerve interface (Zhu et al., 2002;Miller et al., 2008;Botros and Psarros, 2010) and the physiological status of the auditory nerve (Ramekers et al., 2014;Strahl et al., 2016). Clinically, the eCAP is generally evaluated by examining the main peaks, namely the first negative peak (N1) and the first positive peak (P1) (Stypulkowski and van den Honert, 1984;Lai and Dillier, 2000;Abbas et al., 1999;Kim et al., 2011;Alvarez et al., 2011;He et al., 2017). However, the temporal properties of the eCAP are often overlooked. It has been shown that the acoustically evoked compound action potential (CAP) amplitude was linearly correlated with the number of activated nerve fibers (Goldstein and Kiang, 1958;Versnel et al., 1992a). It is generally assumed that the neural response of each single nerve fiber, called the unitary response (UR), is constant and that all URs contribute equally to the CAP (Goldstein and Kiang, 1958;Prijs, 1986;Versnel et al., 1992a). In this study, we assume that this unitary response concept also holds for the eCAPs (e.g., van Gendt et al., 2019), since the eCAP is the superposition of many action potentials from individual auditory nerve fibers in response to an electric stimulus over time. Hence, the eCAP can be described as the convolution of a UR with a compound discharge latency distribution (CDLD), according to equation (1): where t is time, CDLD is a probability density function, and t is the variable of integration. The CDLD weights all URs of each excited nerve fiber over time, and it reflects the synchronicity (i.e., the temporal properties) of the excited nerve fibers. The area under the CDLD (AUCD) indicates the exact number of excited fibers.
The temporal information contained in the CDLD can potentially reflect additional, valuable information that the eCAP amplitude does not show directly. For instance, in patients with CIs, speech perception has been related to auditory nerve fiber survival and function and the number of spiral ganglion cells (Khan et al., 2005;Fayad and Linthicum, 2006;Ramekers et al., 2014;Seyyedi et al., 2015). Consequently, the AUCD might serve as a predictor of the survival and function of auditory nerves. Additionally, the CDLD can be used to study the mechanisms underlying the double peaks in eCAPs (Stypulkowski and van den Honert, 1984;van de Heyning et al., 2016). In these double-peaked eCAPs, the identity of the firing neuron population remains unclear; i.e., it remains unknown whether each peak represents a distinct population or both peaks are evoked by the same group of neurons.
Several studies have reported that there is a relationship between acoustically evoked CAPs and the underlying single fiber discharge patterns. This relationship was found in CAPs recorded in guinea pigs with the convolution model given in equation (1) (e.g., Wang, 2005;Dolan et al., 1983). From those recordings, the UR of guinea pigs was derived (Versnel et al., 1992a, 1992b. In some studies, the inverse problem was studied; i.e., predicting the firing properties by directly deconvolving CAPs and eCAPs with a known UR (Hoke et al., 1979;Charlet de Sauvage et al., 1987;Strahl et al., 2016). In the study by Strahl et al., this method was applied to investigate the CDLD of human eCAPs with the guinea pig UR, and a two-Gaussian component CDLD was derived. When we reproduced their method on human patient data, we found a physiologically unrealistic CDLD, with negative phases and sharp peaks (Fig. 2). Strahl et al. (2016) corrected the negative phases and the high-frequency components by filtering and shifting the CDLD. However, the collective URs that contribute to the eCAP dictate that the CDLD starts after the onset of the electric stimulus. Therefore, the CDLD should be strictly zero before the onset of the stimulus and positive after its onset. Alternatively, it would be better to improve the deconvolution model to obviate the need to postprocess the CDLDs.
Here, we present an iterative method to model the deconvolution computation. In this method, an eCAP, calculated by convoluting a UR model and a two-Gaussian component CDLD model, was optimized to match a recorded eCAP by minimizing the fitting error. The recorded eCAPs were used as input for this iterative method to obtain the UR and CDLDs. Based on this method, two steps were performed to investigate the temporal information in human eCAPs: in the first step, human UR was investigated; in the second step, the two-Gaussian component CDLDs were derived. Some studies reported that the eCAP amplitude has a proportional relationship with the number of excited nerve fibers (e.g., Versnel et al., 1992a;Miller et al., 1998Miller et al., , 1999Westen et al., 2011). However, the synchronicity of the excited nerve fibers could also affect the eCAP amplitude. In comparison to the eCAP amplitude, the AUCD, however, can more accurately reflect the number of the excited nerve fibers and give information on the synchronicity of the excited nerve fibers. Thus, we investigated whether the AUCD rises proportionally with the eCAP amplitude. The differences of the temporal information in eCAPs between children and adults were investigated. After the derivation of this CDLD model of Strahl et al. (2016), we further explored the optimal number of Gaussian components to parameterize the CDLD. To this end, we designed a multiple-Gaussian component CDLD model for predicting recorded eCAPs, and we varied the number of components.
Accordingly, in this study, we aimed to develop an iterative deconvolution model that did not depend on any CDLD postprocessing to explore the temporal information contained in human eCAPs.

Patients and recordings
The eCAPs used in this study were obtained intraoperatively from 111 patients that had undergone CI implantations at the Leiden University Medical Center (Table 1). These eCAPs were recorded as part of the clinical intraoperative routine to assess CI function. All patients received a HiRes90K device (Advanced Bionics, Sylmar, CA), either with a 1J or a Mid-Scala electrode array. These electrode arrays consisted of 16 electrode contacts (numbered from 1 to 16 in apical to basal order). The eCAPs were recorded with the forward masking paradigm provided in the Research Studies Platform Objective Measures (RSPOM) software program (Advanced Bionics, Sylmar, CA). The eCAPs were measured on eight odd electrode contacts with stimulus levels ranging from 50 to 500 CU. The eCAP signal analysis was performed automatically by the RSPOM program (for details, see Biesheuvel et al., 2018). In brief, the eCAPs were evoked using monopolar charge-balanced, biphasic pulses (32 ms/phase) and recorded with a sampling rate of 56 kHz and a gain of 300. Raw eCAP recordings were 1.7 ms in duration and were filtered with a zero-phase shift, low-pass filter, using a cut-off frequency of 8 kHz. The N1 peak was identified as the minimum over the period from 180 to 490 ms, and P1 as the maximum from 470 to 980 ms after the end of stimulation.
The eCAP amplitude was defined as the voltage difference between P1 and N1. After an automated analysis, the identified N1 and P1 peaks were visually inspected.
The noise level was defined as the average of the tail section of the eCAP, i.e., the last 30 samples of the recorded eCAP. It was assumed that no possible remaining neural response or stimulus artifact was present in this tail section (Biesheuvel et al., 2018). Similarly, at a baseline level of eCAPs, there is no neural response or remaining artifact such that this should be mathematically equal to zero (e.g., Prijs, 1986;Charlet de Sauvage et al., 1987). Thus, we used the average level of the tail section as the baseline of the recorded eCAPs. The signal-to-noise ratio (SNR) was defined as the eCAP amplitude divided by the noise amplitude. Then, the eCAP was verified using a semiautomatic method programmed using MAT-LAB (Mathworks, Natick, MA, USA) with two criteria: the eCAP amplitude was larger than 20 mV; the SNR of the eCAP exceeded þ13 dB. If eCAP recordings did not meet these criteria, they were excluded.
The eCAP waveforms were pre-processed before we analysed them with the deconvolution model. The baseline of each recorded eCAP was corrected to zero. Then, 50 additional samples were added to the start and end of the eCAP waveforms by performing a linear extrapolation to zero, to ensure that the entire eCAP waveform was included in the deconvolution analysis and to avoid introducing distortion with the deconvolution algorithm. This extrapolation only influenced the CDLD before and after the recording window (Strahl et al., 2016). We analysed a total of 4982 eCAPs.

Deconvolution model
To explore the temporal information in eCAPs, according to Eq.(1), we modelled the eCAPs as the convolution of a UR model  with a CDLD model. Because the human UR was thought to be similar to the guinea pig UR (Briaire and Frijns, 2005;Whiten, 2007), we applied a guinea pig UR model in the present study (Versnel et al., 1992a), as shown in equation (2).
The UR consisted of a negative (N) and positive (P) phase. The transition point between the negative phase and the positive phase was defined as t 0 . Thus, U ¼ U N and s ¼ s N for t < t 0 ; and U ¼ U P and s ¼ s P for t > t 0 , where s N and s P described the widths (s) of the negative and positive phases of the UR, respectively. The U N and U P described the magnitudes (V) of the two peaks.
Consistent with Strahl et al. (2016), the CDLD model consisted of two Gaussian components, as shown in equation (3).
where N represents a Gaussian distribution; the variables a 1 ; m 1 and s 1 belong to the early Gaussian component (in time), and the variables, a 2 ; m 2 and s 2 belong to the late Gaussian component. The a 1 and a 2 are the peak amplitudes; the m 1 and m 2 are the peak latencies; and the s 1 and s 2 are the peak widths.
Subsequently, the UR and CDLD were used to predict the recorded eCAP waveforms with a deconvolution fitting error minimization routine (DMR). The DMR iteratively optimized the parameters of both UR and CDLD by minimizing the fitting error with a least-squares curve fit using MATLAB. The UR had to be solved before the temporal information could be derived. To this end, we performed two steps, as shown in Fig. 3.

The derivation of the human UR
In step one, a human UR was estimated. The parameters of the UR model (Eq. (2)) and the CDLD model (Eq. (3)) were simultaneously, iteratively adjusted with the DMR to approximate the recorded eCAPs (Fig. 3). To obtain realistic CDLDs and URs, the boundaries of the variables for the UR and CDLD models were iteratively varied to restrict the DMR. The boundary limits of the deconvolution model were based on the parameters of guinea pig contributing auditory nerve fibers and that the UR was identical between human subjects, we derived a human UR by averaging all the URs estimated from eCAPs, across subjects, electrode contacts, and stimulus levels.

The derivation of CDLDs
In step two, the temporal properties of eCAPs recorded in humans were analysed with our iterative deconvolution method. With a fixed human UR, as derived in step one (Fig. 3), we could optimize the parameters of the CDLD model. Because of the fixed UR, the UR and CDLD models could not interact with each other, so that all the temporal information in eCAPs was forced into CDLDs. According to guinea pig UR, the boundaries of the variables of the eCAP waveforms consist of an unusually large P1 and a small N1, and the ratio of the P1 to the N1 is larger than 1. These deviant eCAPs cannot be predicted by our deconvolution model, because the convolution of the human UR, consisting of a large negative phase and a small positive phase (green line, Fig. 1), with a strictly positive CDLD, cannot generate such eCAP waveforms. An example of the deviant eCAPs was shown in Fig. 4C (green line). Therefore, these 322 eCAPs were excluded.

Analysis of the temporal information in eCAPs in CDLDs
As explained in section 1 (Introduction), we expected the temporal information in eCAPs to be captured in CDLDs. First, the histograms of 6 CDLD parameters (in Eq. (4)) derived from 4660 eCAPs were plotted individually in Fig. 5. Second, because a indicated the CDLD magnitude, we assumed that a was positively associated with the eCAP amplitudes. Therefore, we evaluated the association between the a of the CDLD and the eCAP amplitude. Third, more excited nerve fibers led to both a larger eCAP amplitude and a larger AUCD. However, only the AUCD, which was calculated by integrating the CDLD over time, reflected the exact number of activated nerve fibers. Hence, we explored the AUCD as the best proxy for the exact number of activated nerve fibers over time. We also examined the correlation between the AUCD and the eCAP amplitude. The correlation analysis in this section was assessed using Pearson's coefficient using MATLAB.

Differences of the temporal information in eCAPs between children and adults
To explore the differences of the temporal information of the excited auditory nerve fibers between children and adults, we compared the differences of 6 CDLD parameters between child group (<12 years) and adult group (! 12 years) in Table 1 using the Wilcoxon Mann-Whitney U test. The significance level of each comparison was adjusted to 0.0083 using the Bonferroni correction (0.05 divided by 6 comparisons).

The validation of the two-Gaussian component CDLD model
We designed a multiple-Gaussian component CDLD model to determine whether the two-Gaussian component CDLD model was optimal. The formula for the CDLD model was: ða n *Nðm n ; s n ÞÞ where N represents a Gaussian distribution, m represents the number of Gaussian components, a n represents the amplitude, m n Fig. 3. Deconvolution model flow-chart. In step one, both the parameterised unitary response (UR) model (URp) and the parameterised compound discharge latency distribution (CDLD) model (CDLDp) could be manipulated with the deconvolution fitting error minimization routine (DMR, asterisk). In this step, the URs of all eCAPs were derived, and the average of these URs was defined as the human UR (UR h , black square). In step two, the UR h was fixed, and only the CDLDp could be manipulated with the DMR. Then, the CDLDs of all eCAPs were calculated (CDLDs, black square).
represents the peak latency, and s n represents the variance of the latencies in the Gaussian component n.
The fitting errors of simulations using different multiple-Gaussian component CDLD models were assessed by calculating the mean squares error (MSE) in MATLAB.

The unitary response of human auditory nerve fibers
To determine the human UR, we averaged all the URs obtained from the available 4982 eCAPs, by performing the DMR (step one in Fig. 3). Determined with Eq.(2), the final parameters of the mean human UR with standard deviations were: U N ¼ 0.155 ± 0.003 mV, s N ¼ 0.038 ± 0.002 ms, U P ¼ 0.022 ± 0.002 mV, s P ¼ 0.155 ± 0.009 ms, t 0 ¼ À0.128 ± 0.003 ms (Fig. 2). Compared to the guinea pig UR (Versnel et al., 1992a, see Fig. 1), the s N of the negative phase of the human UR was 68% narrower, but 30% higher in magnitude, and the s p of the positive phase of the human UR was slightly broader, and 51% smaller in magnitude.

Typical cases in the deconvolution model
With the DMR method and the human UR derived in step one, we could predict the recorded eCAPs (step two, Fig. 3). The morphological properties of the waveforms, 4660 eCAPs were classified according to visual inspection into two categories (Fig. 4), as described by Lai and Dillier (2000): single-peak eCAPs (75%) and double-peaked eCAPs (19%). Subsequently, we estimated the CDLDs from these eCAPs with the deconvolution model. We found 322 deviant eCAPs (6%), with a ratio of the P1 to N1 larger than 1. As explained above, they could not be predicted with our deconvolution model, and were excluded. The remaining 4660 eCAPs were used in subsequent analyses. Examples of these three eCAP categories (Fig. 4A, B and C) and the corresponding CDLDs (Fig. 4D, E and F, respectively) as predicted with the deconvolution model were shown.

The CDLD parameters
To investigate the synchronicity of the excited nerve fibers, we evaluated the distributions of all CDLD parameters for all eCAPs, recorded at different electrode contacts and different stimulus levels (Fig. 5). We found that all the distributions were skewed; that is, all the parameters of the early and late components (a 1 and a 2 , m 1 and m 2 , s 1 and s 2 ) were not normally distributed, based on a two-sample Kolmogorov-Smirnov test (p < 0.001). The median amplitudes of the two Gaussian components were slightly, but not significantly different (p ¼ 0.15) using the Wilcoxon Mann-Whitney U test. The mean latency, m 1 was significantly different from m 2 (p < 0.05). Furthermore, m 1 displayed a smaller degree of dispersion than m 2 (standard deviations: 0.05 and 0.12 ms, respectively). The average width of the early Gaussian component of CDLDs was significantly different from the average width of the late Gaussian component using the Wilcoxon Mann-Whitney U test (p < 0.01), but the s 1 displayed a smaller standard deviation than s 2 (0.03 and 0.07 ms, respectively). Moreover, the two CDLD amplitude parameters (a 1 , a 2 ) were correlated with the eCAP amplitudes using Spearman correlation coefficient ( Fig. 5; linear regression, r 1 ¼ 0:85; p 1 < 0:001; r 2 ¼ 0:66; p 2 < 0:001). Table 2 shows the average (with standard deviation) and the median (with median deviation) of the 6 CDLD parameters.

Relationship between the eCAP amplitude and the AUCD
We investigated whether the AUCD increases proportionally with the eCAP amplitude. Fig. 6 shows the AUCD plotted against the eCAP amplitude. As anticipated, the AUCD was significantly correlated with the eCAP amplitude (r ¼ 0.83, p < 0.001). Of note, different eCAPs with the same amplitudes could lead to very different CDLDs. For instance, two different eCAPs with the same amplitude (1 mV) had corresponding AUCDs that ranged from 200 to 500. This result indicated that the electrical stimulation did not necessarily activate the same number of nerve fibers each time, even when two eCAPs displayed the same amplitude.

Temporal information in CDLDs between children and adults
We compared the differences of 6 CDLD parameters between the group of children and adults. The significance level of each comparison was corrected to 0.0083 using Bonferroni correction.
Four CDLD parameters (a 1 , m 2 , s 1 and s 2 ) showed significant differences between children and adults. We did not observe significant differences for the parameters a 2 and m 1 . The averages (with standard deviation) of CDLD parameters between the two groups are shown in Table 3.

Validation of the two-Gaussian component CDLD model
Next, we determined whether our two-Gaussian component CDLD model was the most optimal model. We tested models with 1e6 Gaussian components in the CDLD model (m in Eq. (4)), and the fitting error after performing a DMR was evaluated for all the modelled eCAPs (mean squares error, MSE). It turned out that when the m was increased above 2, the fitting errors dropped just slightly (Fig. 7). Apparently, a multi-Gaussian component CDLD model gained no substantial benefit by increasing the number of components beyond 2, meaning that a two-Gaussian component CDLD (Eq. (3)) was the best model.

Discussion
In this study, a model was developed and tested that deconvolved human eCAPs into CDLDs, based on the UR assumption of  auditory nerve fibers (Goldstein and Kiang, 1958;Strahl et al., 2016). As a part of this model, we estimated a human UR that proved to be different from the UR of guinea pigs (Versnel et al., 1992a). To the best of our knowledge, this study was the first to describe a human UR. We modelled the CDLDs underlying human eCAPs to describe the number of electrically excited auditory nerve fibers and their latency. Using the CDLD model, we were able to show differences in temporal characteristics of eCAPs between children and adults were found.

The UR of human auditory nerve fibers
To derive reliable CDLDs from human eCAPs through deconvolution, a representative human UR is critical (Kiang et al., 1976;Wang, 2005;Schoonhoven et al., 1989;Versnel et al., 1992a). However, the human UR waveform had not been previously published. To the best of our knowledge, no modelling studies or electrophysiological recordings have described the human UR in auditory nerve fibers. Previous studies assumed that the human UR was similar to that of guinea pigs (Briaire and Frijns, 2005;Whiten, 2007;Strahl et al., 2016). However, this may not hold true, given that the cochlea in guinea pigs is quite different in size and shape from the cochlea in humans (Nadol, 1988). In addition, the cell bodies of spiral ganglion cells are not myelinated in humans, but they are in guinea pigs. These differences could lead to a different UR. Therefore, we aimed to derive a human UR based on human eCAPs. As a starting point, we used the UR function of guinea pigs published by Versnel et al. (1992a), in combination with wider boundary limits for fitting. We found that the human UR differed from the guinea pig UR. Hence, our modelled human UR was the first attempt to describe the UR of human auditory nerves. Compared to the guinea pig UR, the modelled human UR had a steeper negative component, but a slightly wider and shallower positive peak (Fig. 1). A possible explanation is that the absence of myelin in the human cell body can reduce the neural conduction velocities (Susuki, 2010) and result in a delayed UR. Nadol (1988) has reported that the different cochlear morphology of cochlea in human and guinea pig, i.e., the size of the cochlea and the number of cochlear turns, may lead to different action potential waveforms. To further understand the differences of UR between human and guinea pig, more anatomical and electrophysiological studies are needed.
Of note, the assumption that URs are identical between fibers has not been fully validated. For instance, fibers have different fiber diameters, fiber-to-electrode distances, and response properties, which might trigger different URs, and these URs might contribute to eCAPs differently. However, some research has suggested that nerve fiber diameters were comparable at different locations in the cochlea (Liberman and Oliver, 1984) and that the URs contributed by different fibers along the cochlea were not significantly different (Miller et al., 1999). In the present study, we found significant correlations between the eCAP amplitudes and the CDLD parameters a 1 and a 2 , and between the eCAP amplitudes and the AUCD.
These findings indicate that the eCAP amplitude increases when more auditory nerve fibers are excited by electrical stimulation and that these nerve fibers fire with a higher level of synchronicity. These outcomes are consistent with the assumption that the (e)CAP amplitude is linearly correlated with the number of activated nerve   fibers (Goldstein and Kiang, 1958;Versnel et al., 1992a). In this study, the recorded eCAPs from children and adults were effectively predicted using the same UR. This finding supports the assumption that the UR is identical between fibers and across subjects (e.g., Goldstein and Kiang, 1958;Prijs, 1986 ;Versnel et al., 1992a). Based on these observations, we cautiously assumed that the UR of human auditory nerve fibers was constant. To address this assumption further, more modelling studies or electrophysiological recordings studies on human auditory nerve fibers are required.

The temporal information of eCAPs contained in CDLDs
In this study, we validated that the two-Gaussian component CDLD was the best model. We constructed a multi-Gaussian component CDLD model (with 1e6 components) to predict the recorded eCAPs by performing DMR. When the number of Gaussian components (n, in Eq. (4)) rose from 1 to 2, the fitting outcome showed a reduced fitting error (78%). When the n increased from 2 to 3, the fitting outcome showed that little additional benefit was gained (4%; Fig. 7). These results indicated that the two-Gaussian component CDLD was the best model. This finding is consistent with the findings of Strahl et al. (2016), who also described a two-Gaussian component CDLD. The human UR and the CDLD could interact with each other in step one as the parameters of the UR and the CDLD were both manipulated. Consequently, the temporal information in eCAPs could be demonstrated both in URs and CDLDs and locally but not globally optimal parameters of the CDLD model were derived. Therefore, the CDLDs derived in step one cannot reliably reflect the temporal information in eCAPs. To address this issue, we performed our iterative deconvolution again using a fixed UR so that only the parameter of CDLD can be optimized in step two (Fig. 3). During this iterative procedure, all the temporal information was encoded into CDLDs and these CDLDs accurately reflected the temporal information in eCAPs.
Most of the eCAP waveforms, both the single peak eCAPs and the double peak eCAPs, in our study were fit better with a double Gaussian component CDLD model (Fig. 4) than with a single Gaussian component model. Our finding suggested that eCAP waveforms that appear to have a single peak could arise from a two-Gaussian component CDLD. If true, it follows that CDLDs might consist of two independent components that originate from two separate groups of neural responses. This hypothesis was in line with the concept of double group neural responses proposed by Stypulkowski and van den Honert (1984) and with the simulations including a combination of neural responses arising from axons and peripheral processes (Lai and Dillier, 2000). In their findings, the early Gaussian component of CDLDs could be attributed to direct excitation of the axonal process in the modiolus proximal to the spiral ganglion cell, and the late CDLD component could be attributed to the activation of the axon peripheral to the cell body of the bipolar ganglion cell. This hypothesis was supported by our finding that the time interval between m 1 and m 2 (0.2 ms, see in Fig. 5) was shorter than the absolute refractory period of these nerve fibers (approximately 0.45 ms), as reported by He et al. (2017). Therefore, we could rule out the possibility that the neural responses in the two-component of CDLD might have originated from the same group of auditory nerve fibers. A recent study, by Finley et al., 2019(presented at CIAP 2019, indicated that multiple neural response sites with different waveform morphologies, latencies, magnitudes, and scalar distributions could contribute to differences in the eCAPs measured in the cochlea. A limitation of our deconvolution model was that it relied on the eCAP waveforms. In some cases, the eCAPs had deviant waveforms that could not be simulated by our deconvolution model. However, those instances were rare (approximately 6%).
We determined the distributions of the CDLD parameters shown in Eq. (3). The distribution of the s 1 showed smaller means and variations compared to the distributions of the s 2 . The average of a 2 /a 1 (0.96) in our study was quite similar to that reported by Strahl in humans (apex: 0.96, middle: 0.86, base: 0.85). However, our m 1 (0.36 ms), m 2 (0.55 ms), s 1 (0.071 ms), and s 2 (0.15 ms) values were smaller than those reported by Strahl (on average: m 1 ¼ 0.52 ms, m 2 ¼ 0.9 ms, s 1 ¼ 0.14, and s 2 ¼ 0.27). When we used the guinea pig UR to derive the CDLD from our eCAPs, the m 1 and m 2 changed slightly in the direction of the value reported by Strahl et al. (2016): shifting from 0.36 ms to 0.55 mse0.41 ms and 0.69 ms, respectively, with a larger variance. This indicated that using a guinea pig UR leads to a poorer fitting of human CDLDs. Thus, these differences might be attributable to the UR used in Strahl's study, which was derived from guinea pig eCAPs, rather than human eCAPs. The distance between the electrode contact and the nerve fibers may affect CDLDs. Previous studies reported that the perimodiolar electrode arrays can yield a lower threshold in comparison to electrodes located close to the outer wall (e.g., Frijns et al., 1995;Briaire and Frijns, 2000). When an electrode is located closer to the modiolus, less current is required to excite the auditory nerve fibers (e.g., Kang et al., 2015). Therefore, a shorter distance to modiolus can lead to activation of more nerve fibers and a larger a 1 , a 2 and AUCD can be obtained. Conversely, a larger distance could result in a smaller a 1 , a 2 and AUCD.
We found that analysing CDLDs had at least two advantages over analysing eCAP amplitudes directly. First, when studying latency effects, the CDLD could more precisely reflect the latency of eCAPs over time than the N1 and P1 of eCAPs. Second, the AUCD (i.e., the integral of CDLD over time) could provide more accurate information about the number of excited nerve fibers. Because only healthy fibers can be activated, the AUCD might also reflect the survival of nerve fibers (Khan et al., 2005;Fayad and Linthicum, 2006). On the other hand, we found that eCAPs with equal amplitudes could lead to different AUCDs in the deconvolution model (Fig. 6), indicating that different eCAPs with same amplitude could arise from very different numbers of excited nerve fibers. Contrary to the unitary response concept, this would indicate that the eCAP amplitude could not accurately indicate the number of excited auditory nerve fibers.
In this study, significant differences of the temporal information between children and adults were revealed by calculating CDLDs using iterative deconvolution (Table 3). These differences may be attributable to the observation that auditory nerve fibers of hearing-impaired adults undergo significant degeneration over the years (e.g., Abbas and Brown, 1991). For instance, in comparison to the group of children, we observed larger peak widths (s 1 and s 2 ) in adults, presumably reflecting that the excited nerve fibers showed a lower level of synchronicity in this group. Compared with adults, a larger a 1 in the group of children may indicate that more nerve fibers can be excited in their central axon. A larger m 2 observed in the adult group likely implies more severe degeneration of the peripheral process compared with the children.
The outcomes of the present study suggested some potential applications for future clinical practice. First, although many studies have investigated the relationship between speech performance and nerve fiber survival (Kawano et al., 1998;Khan et al., 2005;Fayad and Linthicum, 2006;Xu et al., 2012), the findings were inconsistent. A likely explanation might be that most investigators focused mainly on the eCAP magnitudes, which could not precisely indicate the number of activated nerve fibers and their latency. With our deconvolution model, it is possible to investigate the relationship between the temporal information in eCAPs and speech perception in patients with CIs. For instance, our deconvolution model provided the growth function of the AUCD and the threshold of the AUCD, which might associate well with speech perception in patients with CIs. Additionally, consistent with the finding by Strahl et al. (2016) that some of the CDLD parameters could indicate the degeneration of nerve fibers, we found that the AUCDs derived with our method might be more accurate than the eCAP amplitude for indicating nerve fiber survival in patients because only the AUCD indicates the number of excited nerve fibers. Furthermore, our deconvolution model has low computational complexity. With the estimated UR, the computation of each CDLD could be completed in 0.1 s from the recorded eCAP using our deconvolution model in MATLAB. Thus, our deconvolution could be potentially integrated into clinical software to derive the temporal information of eCAPs in near-real-time of CI recipients.

Conclusions
This study described an iterative deconvolution model, based on the UR hypothesis, to derive the CDLD from recorded human eCAPs. We estimated a human version of the UR, which was not available previously. Importantly, we found that the human UR differed from the guinea pig UR. With the estimated human UR, we derived the CDLDs of 4660 eCAPs. We demonstrated that CDLDs had advantages over the more commonly used eCAP amplitude because they better reflected the temporal properties of eCAPs. Therefore, CDLDs provided better estimates of the number of excited auditory nerve fibers and their firing latencies.