Note. The following notebook contains code in addition to text and figures. By default, the code has been hidden. You can click the icon that looks like an eye in the toolbar above to show the code. To run the code, click the cell menu, then "run all".

Amplitude modulation

This notebook shows how the model responds to an amplitude modulated signal. See Basic Model for details of the model.


In [ ]:
###### IMPORT AND UTILITY FUNCTIONS

%matplotlib inline
from brian2 import *
from model_explorer_jupyter import *
import joblib
import ipywidgets as ipw
from collections import OrderedDict
from scipy.ndimage.interpolation import zoom
from scipy.ndimage.filters import gaussian_filter
from scipy import stats
from matplotlib import cm
import matplotlib.patches as patches
import warnings
warnings.filterwarnings("ignore")
BrianLogger.log_level_error()

# Used for smoothing the plots, applies a Gaussian filter but works properly with nan values
def nan_gaussian_filter(x, sigma, num_passes):
    z = full_like(x, nan)
    for cursigma in linspace(sigma, 0, num_passes+1)[:-1]:
        y = gaussian_filter(x, cursigma, mode='nearest')
        z[isnan(z)] = y[isnan(z)]
    return z

# Take a series of x, y points and plot a density map using kernel density estimation
# N is the grid size for the density image
def density_map(x, y, N, xmin=None, xmax=None, ymin=None, ymax=None):
    # Peform the kernel density estimate
    if xmin is None:
        xmin = amin(x)
    if xmax is None:
        xmax = amax(x)
    if ymin is None:
        ymin = amin(y)
    if ymax is None:
        ymax = amax(y)
    xx, yy = mgrid[xmin:xmax:N*1j, ymin:ymax:N*1j]
    positions = vstack([xx.ravel(), yy.ravel()])
    values = vstack([x, y])
    kernel = stats.gaussian_kde(values)
    f = np.reshape(kernel(positions).T, xx.shape)
    extent = (xmin, xmax, ymin, ymax)
    return f.T, extent

def plot_density_map(x, y, N, xmin=None, xmax=None, ymin=None, ymax=None, cmap=cm.afmhot_r, **args):
    img, extent = density_map(x, y, N, xmin=xmin, ymin=ymin, xmax=xmax, ymax=ymax)
    imshow(img, origin='lower left', aspect='auto', interpolation='nearest',
           extent=extent, cmap=cmap,
           vmin=0, vmax=amax(img)/0.7,
           **args
           )

progress_slider, update_progress = brian2_progress_reporter()

mem = joblib.Memory(cachedir="joblib", verbose=0)

In [ ]:
%%html
<!-- hack to improve styling of ipywidgets sliders -->
<style type="text/css">
.widget-label {
    min-width: 35ex;
    max-width: 35ex;
}
.widget-hslider {
    width: 100%;
}
.widget-hprogress {
    width: 100%;
}

</style>

Interactive with linear ANF rate-level function

Here, we assume that the firing rate variable $\rho$ is replaced by an amplitude modulated firing rate:

$$\rho(t)=\rho_{mean}\cdot (1+m\sin(2\pi f_m t))$$

$\rho(t)$ varies with a modulation frequency of $f_m$ with a mean value of $\rho_{mean}$ with a modulation depth of $m$. Note that compression and adaptation are not present in this model.


In [ ]:
def mtf(log_fm_range_Hz, num_fm=10,
        anf_rate_mean_Hz=200,
        modulation_depth=0.25,
        mu_mean=1.75, num_anf=50, tau_ms=6, inh=0.0,
        refractory_ms=0.6,
        repeats=50, duration_ms=1000,
        plotresults='all', plotkwds={'c': 'b'},
        show_mode_locking_diagrams=False,
        mode_locking_extent=None,
        ):
    # Set parameters
    log_fm_min_Hz, log_fm_max_Hz = log_fm_range_Hz
    fm_min_Hz = 2**log_fm_min_Hz
    fm_max_Hz = 2**log_fm_max_Hz
    duration = duration_ms*ms
    fm = fm_max_Hz*Hz
    tau = tau_ms*ms
    refractory = refractory_ms*ms
    anf_rate_mean = anf_rate_mean_Hz*Hz
    weight = mu_mean/(num_anf*tau*anf_rate_mean*(1-inh))
    lfm_min = log2(fm_min_Hz)
    lfm_max = log2(fm_max_Hz)
    # Define and run the model
    eqs = '''
    fmi = int(i/repeats) : integer
    lfm = fmi/(num_fm-1.0)*(lfm_max-lfm_min)+lfm_min : 1
    fm = 2**lfm*Hz : Hz
    am = 1+modulation_depth*sin(2*pi*fm*t) : 1
    anf_rate = anf_rate_mean*am : Hz
    anf_rate_exc = anf_rate : Hz
    anf_rate_inh = anf_rate*inh : Hz
    mu_exc = weight*num_anf*tau*anf_rate_exc : 1
    mu_inh = weight*num_anf*tau*anf_rate_inh : 1
    sigma2_exc = weight*mu_exc : 1
    sigma2_inh = weight*mu_inh : 1
    mu = mu_exc-mu_inh : 1
    sigma = sqrt(sigma2_exc+sigma2_inh) : 1
    dv/dt = (mu-v)/tau+sigma*xi*tau**-0.5 : 1 (unless refractory)
    '''
    G = NeuronGroup(repeats*num_fm, eqs, threshold='v>1', reset='v=0',
                    refractory=refractory, method='heun')
    spikemon = SpikeMonitor(G)
    statemon = StateMonitor(G, 'v', record=[0])
    run(duration)
    # Compute rMTF and tMTF
    rMTF = zeros(num_fm)
    dMTF = zeros(num_fm)
    tMTF = zeros(num_fm, dtype=complex64)
    ntMTF = zeros(num_fm, dtype=int)
    fmi = arange(num_fm)
    lfm = fmi/(num_fm-1.0)*(lfm_max-lfm_min)+lfm_min
    fm = 2**lfm
    if plotresults=='all':
        clf()
    for i, train in enumerate(spikemon.spike_trains().values()):
        f = fm[i/repeats]*Hz
        n = int(floor(duration*f))
        tmax = n/f
        train.sort()
        if show_mode_locking_diagrams:
            if i%repeats==0:
                allx = []
                ally = []
            difftrain = diff(train)
            allx.append(difftrain[:-1]/ms)
            ally.append(difftrain[1:]/ms)
            if (i+1)%repeats==0:
                if plotresults=='all':
                    figure(figsize=(2.5, 2.5))
                    title("%d Hz" % f)
                if plotresults=='all' or plotresults=='modelocking-%d' % int(f):
                    allx = hstack(allx)
                    ally = hstack(ally)
                    cx = array(allx)
                    cy = array(ally)
                    if mode_locking_extent is not None:
                        xmin, xmax, ymin, ymax = mode_locking_extent
                    else:
                        xmin = amin(cx)
                        ymin = amin(cy)
                        xmax = amax(cx)
                        ymax = amax(cy)
                    plot_density_map(cx, cy, 150,
                                     xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
                    xlabel("ISI $n$ (ms)")
                    ylabel("ISI $n+1$ (ms)")
                    tight_layout()
                    savefig('modelocking-%dHz.pdf' % int(f))
        train = train[train<tmax]
        rMTF[i/repeats] += len(train)
        dMTF[i/repeats] += float(tmax)
        tMTF[i/repeats] += sum(exp(1j*2*pi*fm[i/repeats]*train/second))
        ntMTF[i/repeats] += len(train)
    #I = ntMTF>10 # only plot tMTF if we have >10 spikes
    tMTF = abs(tMTF)/ntMTF # rMTF is the number of spikes at this point
    # Rayleigh statistic significance from Sayles et al 2013
    RS = 2*(tMTF**2)*ntMTF # Rayleigh statistic
    I = RS>13.8 # p=0.001
    rMTF /= dMTF
    
    # Plot results
    if plotresults=='all':
        figure(figsize=(10, 5))
        subplot(211)
    if plotresults=='all' or plotresults=='raster':
        plot(spikemon.t/ms, spikemon.i, ',k')
        xlabel('Time (ms)')
        ylabel('Repeat number')
        title('Spike raster plot')
        xlim(0, duration/ms)
        ylim(0, repeats*num_fm)
    
    if plotresults=='all':
        subplot(234)
    if plotresults=='all' or plotresults=='rmtf':
        semilogx(fm, rMTF, '-', basex=2, **plotkwds)
        xlabel('Modulation frequency (Hz)')
        ylabel('Firing rate (sp/s)')
        xlim(fm_min_Hz, fm_max_Hz)
        xt, _ = xticks()
        xticks(xt, map(str, map(int, xt)))
        xlim(fm_min_Hz, fm_max_Hz)
        axhline(200, ls='--', c='k')
        ylim(0, max(200, amax(rMTF)))
        title('Rate MTF')
    
    if plotresults=='all':
        subplot(235)
    if plotresults=='all' or plotresults=='tmtf':
        #semilogx(fm, tMTF, '--', basex=2, **plotkwds)
        semilogx(fm[I], tMTF[I], '-', basex=2, **plotkwds)
        xlabel('Modulation frequency (Hz)')
        ylabel('Vector strength')
        xlim(fm_min_Hz, fm_max_Hz)
        xt, _ = xticks()
        xticks(xt, map(str, map(int, xt)))
        xlim(fm_min_Hz, fm_max_Hz)
        ylim(0, 1)
        title('Temporal MTF')

    if plotresults=='all':
        subplot(236)
    if plotresults=='all' or plotresults=='modulationgain':
        #semilogx(fm, 20*log10(2*tMTF/modulation_depth), '--', basex=2, **plotkwds)
        semilogx(fm[I], 20*log10(2*tMTF[I]/modulation_depth), '-', basex=2, **plotkwds)
        xlabel('Modulation frequency (Hz)')
        ylabel('Modulation gain (dB)')
        xlim(fm_min_Hz, fm_max_Hz)
        xt, _ = xticks()
        xticks(xt, map(str, map(int, xt)))
        xlim(fm_min_Hz, fm_max_Hz)
        ylim(-25, 25)
        axhline(0, ls='--', c='k')
        title('Modulation gain')

    tight_layout();

display(ipw.interact(mtf,
             plotresults=ipw.fixed('all'),
             plotkwds=ipw.fixed({'c': 'b'}),
             mode_locking_extent=ipw.fixed(None),
             show_mode_locking_diagrams=ipw.Checkbox(
                 value=False,
                 description="Show mode locking diagrams"),
             tau_ms=ipw.FloatSlider(
                 min=0.1, max=20.0, step=0.1, value=10.0,
                 continuous_update=False,
                 description=r"Membrane time constant $\tau$ (ms)"),
             refractory_ms=ipw.FloatSlider(
                 min=0, max=5, step=0.05, value=1,
                 continuous_update=False,
                 description=r"Refractory period $t_\mathrm{ref}$ (ms)"),
             mu_mean=ipw.FloatSlider(
                 min=0, max=5, step=0.05, value=1.25,
                 continuous_update=False,
                 description=r"Mean current at mean firing rate $\mu$"),
             log_fm_range_Hz=ipw.IntRangeSlider(
                 min=0, max=13, step=1, value=(0, 10),
                 continuous_update=False,
                 description=r"Log modulation frequency range $\log_2 f_m/\mathrm{Hz}$"),
             num_fm=ipw.IntSlider(
                 min=5, max=250, step=5, value=10,
                 continuous_update=False,
                 description=r"Modulation frequency number of points"),
             num_anf=ipw.IntSlider(
                 min=1, max=100, step=1, value=50,
                 continuous_update=False,
                 description=r"Number of input synapses"),
             inh=ipw.FloatSlider(
                 min=0, max=1, step=0.05, value=0,
                 continuous_update=False,
                 description=r"Fraction of inhibition $\alpha$"),
             modulation_depth=ipw.FloatSlider(
                 min=0, max=1, step=0.05, value=0.25,
                 continuous_update=False,
                 description=r"ANF modulation depth $m$"),
             anf_rate_mean_Hz=ipw.IntSlider(
                 min=10, max=500, step=10, value=200,
                 continuous_update=False,
                 description=r"Mean ANF firing rate (sp/s)"),
             duration_ms=ipw.IntSlider(
                 min=100, max=10000, step=100, value=1000,
                 continuous_update=False,
                 description=r"Duration (ms)"),
             repeats=ipw.IntSlider(
                 min=5, max=500, step=5, value=50,
                 continuous_update=False,
                 description=r"Repeats"),
             ));

Deafferentation

The following figure shows the effect of deafferentation. This is modelled by keeping the mean weight the same, but reducing the number of input cells from 50 to 10.


In [ ]:
# Quick parameters
num_fm = 10
repeats = 50
duration_ms = 1000
# High quality parameters
# num_fm = 40
# repeats = 250
# duration_ms = 5000

figure(figsize=(7, 7))
# Sustained chopper
subplot(221)
mtf(log_fm_range_Hz=(2, 9), num_fm=num_fm,
    anf_rate_mean_Hz=200,
    modulation_depth=0.25,
    mu_mean=1.25, num_anf=50, tau_ms=10, inh=0.0,
    refractory_ms=1.0,
    repeats=repeats, duration_ms=duration_ms,
    plotresults='tmtf', plotkwds={'c': 'b', 'label': 'Normal', 'lw': 2},
    )
mtf(log_fm_range_Hz=(2, 9), num_fm=num_fm,
    anf_rate_mean_Hz=200,
    modulation_depth=0.25,
    mu_mean=1.25, num_anf=10, tau_ms=10, inh=0.0,
    refractory_ms=1.0,
    repeats=repeats, duration_ms=duration_ms,
    plotresults='tmtf', plotkwds={'c': 'r', 'label': 'Impaired', 'ls': '--', 'lw': 2},
    )
title('Sustained chopper')
#legend(loc='upper left')
# Transient chopper
subplot(222)
mtf(log_fm_range_Hz=(2, 9), num_fm=num_fm,
    anf_rate_mean_Hz=200,
    modulation_depth=0.25,
    mu_mean=1.25, num_anf=50, tau_ms=10, inh=0.5,
    refractory_ms=1.0,
    repeats=repeats, duration_ms=duration_ms,
    plotresults='tmtf', plotkwds={'c': 'b', 'label': 'Normal', 'lw': 2},
    )
mtf(log_fm_range_Hz=(2, 9), num_fm=num_fm,
    anf_rate_mean_Hz=200,
    modulation_depth=0.25,
    mu_mean=1.25, num_anf=10, tau_ms=10, inh=0.5,
    refractory_ms=1.0,
    repeats=repeats, duration_ms=duration_ms,
    plotresults='tmtf', plotkwds={'c': 'r', 'label': 'Impaired', 'ls': '--', 'lw': 2},
    )
ylabel('')
title('Transient chopper')
legend(loc='upper right')
subplot(223)
mtf(log_fm_range_Hz=(0, 10), num_fm=10,
    anf_rate_mean_Hz=200,
    modulation_depth=0.25,
    mu_mean=1.25, num_anf=50, tau_ms=10, inh=0.0,
    refractory_ms=1.0,
    repeats=repeats, duration_ms=duration_ms,
    plotresults='modelocking-47', show_mode_locking_diagrams=True,
    mode_locking_extent=(0, 30, 0, 30),
    )
title('Normal $N=50$')
subplot(224)
mtf(log_fm_range_Hz=(0, 10), num_fm=10,
    anf_rate_mean_Hz=200,
    modulation_depth=0.25,
    mu_mean=1.25, num_anf=10, tau_ms=10, inh=0.0,
    refractory_ms=1.0,
    repeats=repeats, duration_ms=duration_ms,
    plotresults='modelocking-47', show_mode_locking_diagrams=True,
    mode_locking_extent=(0, 30, 0, 30),
    )
title('Impaired $N=10$')
ax = subplot(221)
ax.text(-0.2, 1.1, 'A', transform=ax.transAxes, size=20)
ax = subplot(223)
ax.text(-0.2, 1.1, 'B', transform=ax.transAxes, size=20)
tight_layout()
savefig('fig_am_deafferentation.pdf')
savefig('fig_am_deafferentation.png');

Interactive with nonlinear ANF rate-level function

Here, we assume that the tone is amplitude modulated according to the same function, and we assume a sigmoidal auditory nerve fibre rate-level function. Specifically, we use the rate being a logistic function of the level, as suggested by Sachs and Abbas (1974). We use a function

$$\rho(\theta)=\rho_{spont}+(\rho_{sat}-\rho_{spont})S(\theta)$$

and

$$S(\theta)=\frac{1}{1+a\cdot 10^{-b\cdot\theta}}$$

Where $\theta$ is the sound level, $\rho_{spont}$ is the spontaneous rate and $\rho_{sat}$ is the saturated rate. Assuming that $\theta=0$ is the experimentally determined threshold level where $S(0)=\epsilon$, and $\theta_+$ is the experimentally determined saturation level (also the dynamic range since the threshold is 0) where $S(\theta_+)=1-\epsilon$, for some small $\epsilon$ (we use 0.05). This gives us a solution for the constants $a=1/\epsilon-1$ and $b=(2/\theta_+)\log_{10}a$. Now, write the following function for the amplitude modulation of the pressure:

$$AM(t)=1+m\sin(2\pi f_m t)$$

where $AM(t)$ varies with a modulation frequency of $f_m$ with a mean value of 1 and a modulation depth of $m$. Assume that this function is modulating a tone at level $\overline\theta$. This gives us a firing rate function (converting pressure into level and simplifying) of:

$$\rho(t)=\frac{AM(t)^{20 b}}{AM(t)^{20b}+a10^{-b\overline\theta}}$$

Note that adaptation is not present in this model.


In [ ]:
def mtf_nonlinear(log_fm_range_Hz, num_fm=10,
        anf_rate_range_Hz=(50, 300),
        anf_dynamic_range_dB=60,
        level_mean_dB=30,
        modulation_depth=0.25,
        mu_mean=1.75, num_anf=50, tau_ms=6, inh=0.0,
        refractory_ms=0.6,
        repeats=50, duration_ms=1000,
        ):
    #### Set parameters
    log_fm_min_Hz, log_fm_max_Hz = log_fm_range_Hz
    fm_min_Hz = 2**log_fm_min_Hz
    fm_max_Hz = 2**log_fm_max_Hz
    duration = duration_ms*ms
    fm = fm_max_Hz*Hz
    tau = tau_ms*ms
    refractory = refractory_ms*ms
    # rate-level parameters
    anf_rate_spont_Hz, anf_rate_sat_Hz = anf_rate_range_Hz
    anf_rate_spont = anf_rate_spont_Hz*Hz
    anf_rate_sat = anf_rate_sat_Hz*Hz
    ratelevel_epsilon = 0.05
    ratelevel_a = 1/ratelevel_epsilon-1
    ratelevel_b = 2./anf_dynamic_range_dB*log10(ratelevel_a)
    anf_rate_mean = anf_rate_spont+(anf_rate_sat-anf_rate_spont)/(1.+ratelevel_a*10**(-ratelevel_b*level_mean_dB))
    anf_rate_mid = anf_rate_spont+(anf_rate_sat-anf_rate_spont)/(1.+ratelevel_a*10**(-ratelevel_b*anf_dynamic_range_dB/2.))
    # other computed parameters
    weight = mu_mean/(num_anf*tau*anf_rate_mid*(1-inh))
    lfm_min = log2(fm_min_Hz)
    lfm_max = log2(fm_max_Hz)
    #### Define and run the model
    eqs = '''
    fmi = int(i/repeats) : integer
    lfm = fmi/(num_fm-1.0)*(lfm_max-lfm_min)+lfm_min : 1
    fm = 2**lfm*Hz : Hz
    am = 1+modulation_depth*sin(2*pi*fm*t) : 1
    ratelevel_sigmoid = am**(20*ratelevel_b)/(am**(20*ratelevel_b)+ratelevel_a*10**(-ratelevel_b*level_mean_dB)) : 1
    anf_rate = anf_rate_spont+(anf_rate_sat-anf_rate_spont)*ratelevel_sigmoid : Hz
    anf_rate_exc = anf_rate : Hz
    anf_rate_inh = anf_rate*inh : Hz
    mu_exc = weight*num_anf*tau*anf_rate_exc : 1
    mu_inh = weight*num_anf*tau*anf_rate_inh : 1
    sigma2_exc = weight*mu_exc : 1
    sigma2_inh = weight*mu_inh : 1
    mu = mu_exc-mu_inh : 1
    sigma = sqrt(sigma2_exc+sigma2_inh) : 1
    dv/dt = (mu-v)/tau+sigma*xi*tau**-0.5 : 1 (unless refractory)
    '''
    G = NeuronGroup(repeats*num_fm, eqs, threshold='v>1', reset='v=0',
                    refractory=refractory, method='heun')
    spikemon = SpikeMonitor(G)
    statemon = StateMonitor(G, 'v', record=[0])
    run(duration)
    # Compute rMTF and tMTF
    rMTF = zeros(num_fm)
    dMTF = zeros(num_fm)
    tMTF = zeros(num_fm, dtype=complex64)
    ntMTF = zeros(num_fm, dtype=int)
    fmi = arange(num_fm)
    lfm = fmi/(num_fm-1.0)*(lfm_max-lfm_min)+lfm_min
    fm = 2**lfm
    for i, train in enumerate(spikemon.spike_trains().values()):
        train.sort()
        f = fm[i/repeats]*Hz
        n = int(floor(duration*f))
        tmax = n/f
        train = train[train<tmax]
        rMTF[i/repeats] += len(train)
        dMTF[i/repeats] += float(tmax)
        tMTF[i/repeats] += sum(exp(1j*2*pi*fm[i/repeats]*train/second))
        ntMTF[i/repeats] += len(train)
    #I = ntMTF>10 # only plot tMTF if we have >10 spikes
    tMTF = abs(tMTF)/ntMTF # rMTF is the number of spikes at this point
    # Rayleigh statistic significance from Sayles et al 2013
    RS = 2*(tMTF**2)*ntMTF # Rayleigh statistic
    I = RS>13.8 # p=0.001
    rMTF /= dMTF
    # Plot results
    figure(figsize=(10, 5))
    subplot(211)
    plot(spikemon.t/ms, spikemon.i, ',k')
    xlabel('Time (ms)')
    ylabel('Repeat number')
    title('Spike raster plot')
    xlim(0, duration/ms)
    ylim(0, repeats*num_fm)
    
    subplot(234)
    semilogx(fm, rMTF, '-b', basex=2)
    xlabel('Modulation frequency (Hz)')
    ylabel('Firing rate (sp/s)')
    xlim(fm_min_Hz, fm_max_Hz)
    xt, _ = xticks()
    xticks(xt, map(str, map(int, xt)))
    xlim(fm_min_Hz, fm_max_Hz)
    axhline(200, ls='--', c='k')
    ylim(0, max(200, amax(rMTF)))
    title('Rate MTF')
    
    subplot(235)
    semilogx(fm, tMTF, '--b', basex=2)
    if sum(I):
        semilogx(fm[I], tMTF[I], '-b', basex=2)
    xlabel('Modulation frequency (Hz)')
    ylabel('Vector strength')
    xlim(fm_min_Hz, fm_max_Hz)
    xt, _ = xticks()
    xticks(xt, map(str, map(int, xt)))
    xlim(fm_min_Hz, fm_max_Hz)
    ylim(0, 1)
    title('Temporal MTF')

    subplot(236)
    if sum(tMTF>0):
        semilogx(fm, 20*log10(2*tMTF/modulation_depth), '--b', basex=2)
    if sum(I):
        semilogx(fm[I], 20*log10(2*tMTF[I]/modulation_depth), '-b', basex=2)
    xlabel('Modulation frequency (Hz)')
    ylabel('Modulation gain (dB)')
    xlim(fm_min_Hz, fm_max_Hz)
    xt, _ = xticks()
    xticks(xt, map(str, map(int, xt)))
    xlim(fm_min_Hz, fm_max_Hz)
    ylim(-25, 25)
    axhline(0, ls='--', c='k')
    title('Modulation gain')

    tight_layout()

display(ipw.interact(mtf_nonlinear,
             tau_ms=ipw.FloatSlider(
                 min=0.1, max=20.0, step=0.1, value=10.0,
                 continuous_update=False,
                 description=r"Membrane time constant $\tau$ (ms)"),
             refractory_ms=ipw.FloatSlider(
                 min=0, max=5, step=0.05, value=1,
                 continuous_update=False,
                 description=r"Refractory period $t_\mathrm{ref}$ (ms)"),
             mu_mean=ipw.FloatSlider(
                 min=0, max=5, step=0.05, value=2.0,
                 continuous_update=False,
                 description=r"Mean current at level $\theta_+/2$: $\mu$"),
             log_fm_range_Hz=ipw.IntRangeSlider(
                 min=0, max=13, step=1, value=(0, 11),
                 continuous_update=False,
                 description=r"Log modulation frequency range $\log_2 f_m/\mathrm{Hz}$"),
             num_fm=ipw.IntSlider(
                 min=5, max=250, step=5, value=40,
                 continuous_update=False,
                 description=r"Modulation frequency number of points"),
             num_anf=ipw.IntSlider(
                 min=1, max=100, step=1, value=50,
                 continuous_update=False,
                 description=r"Number of input synapses $N$"),
             inh=ipw.FloatSlider(
                 min=0, max=1, step=0.05, value=0,
                 continuous_update=False,
                 description=r"Fraction of inhibition $\alpha$"),
             modulation_depth=ipw.FloatSlider(
                 min=0, max=1, step=0.05, value=1,
                 continuous_update=False,
                 description=r"Signal modulation depth $m$"),
             anf_rate_range_Hz=ipw.IntRangeSlider(
                 min=0, max=500, step=10, value=(50, 300),
                 continuous_update=False,
                 description=r"ANF firing rate range ($\rho_{spont}$ to $\rho_{sat}$) (sp/s)"),
             anf_dynamic_range_dB=ipw.IntSlider(
                 min=5, max=100, step=5, value=60,
                 continuous_update=False,
                 description=r"ANF dynamic range $\theta_+$ (dB)"),
             level_mean_dB=ipw.IntSlider(
                 min=0, max=100, step=5, value=30,
                 continuous_update=False,
                 description=r"Tone signal level re threshold $\overline\theta$(dB)"),
             duration_ms=ipw.IntSlider(
                 min=100, max=10000, step=100, value=1000,
                 continuous_update=False,
                 description=r"Duration (ms)"),
             repeats=ipw.IntSlider(
                 min=5, max=500, step=5, value=50,
                 continuous_update=False,
                 description=r"Repeats"),
             ));