Modelling the response of the retinal cells to the stimuli
# of the jupyter notebook, but not needed in the source code
from theonerig.core import *
from theonerig.processing import *
from theonerig.utils import *
from theonerig.plotting import *

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import axes3d

Functions used for modelling

sigmoid[source]

sigmoid(x, sigma, amp, x0, y0)

Sigmoid function params:

- x: 1D numpy array at which to evaluate the points
- sigma: steepness of the sigmoid
- amp: amplitude of the sigmoid
- x0: shift in x of the sigmoid
- y0: shift in y of the sigmoid
x = np.linspace(-10,10,100)
plt.figure()
plt.plot(x, sigmoid(x, sigma=1, amp=1, x0=0, y0=0), label="sigma:1 amp:1 x0:0 y0:0")
plt.plot(x, sigmoid(x, sigma=2, amp=1, x0=0, y0=0), label="sigma:2 amp:1 x0:0 y0:0")
plt.plot(x, sigmoid(x, sigma=1, amp=2, x0=0, y0=0), label="sigma:1 amp:2 x0:0 y0:0")
plt.plot(x, sigmoid(x, sigma=1, amp=1, x0=2, y0=0), label="sigma:1 amp:1 x0:2 y0:0")
plt.plot(x, sigmoid(x, sigma=1, amp=1, x0=0, y0=1), label="sigma:1 amp:1 x0:0 y0:1")
_ = plt.legend()

gaussian[source]

gaussian(x, sigma, amp, x0, y0)

Gaussian function params:

- x: 1D numpy array at which to evaluate the points
- sigma: width of the gaussian
- amp: amplitude of the gaussian
- x0: shift in x of the gaussian
- y0: shift in y of the gaussian
x = np.linspace(-10,10,100)
plt.figure()
plt.plot(x, gaussian(x, sigma=1, amp=1, x0=0, y0=0), label="sigma:1 amp:1 x0:0 y0:0")
plt.plot(x, gaussian(x, sigma=2, amp=1, x0=0, y0=0), label="sigma:2 amp:1 x0:0 y0:0")
plt.plot(x, gaussian(x, sigma=1, amp=2, x0=0, y0=0), label="sigma:1 amp:2 x0:0 y0:0")
plt.plot(x, gaussian(x, sigma=1, amp=1, x0=2, y0=0), label="sigma:1 amp:1 x0:2 y0:0")
plt.plot(x, gaussian(x, sigma=1, amp=1, x0=0, y0=1), label="sigma:1 amp:1 x0:0 y0:1")
_ = plt.legend()

sum_of_gaussian[source]

sum_of_gaussian(t, sigma_1, amp_1, x0_1, sigma_2, amp_2, x0_2, y0)

Sum of gaussian, using gaussian params:

- t: 1D numpy array at which to evaluate the points
- sigma_1: width of the 1st gaussian
- amp_1: amplitude of the 1st gaussian
- x0_1: shift in x of the 1st gaussian
- sigma_2: width of the 2nd gaussian
- amp_2: amplitude of the 2nd gaussian
- x0_2: shift in x of the 2nd gaussian
- y0: shift in y of the gaussian (shared for the two gaussians)

diff_of_gaussian[source]

diff_of_gaussian(t, sigma_1, amp_1, x0_1, sigma_2, amp_2, x0_2, y0)

Difference of gaussian, using gaussian params:

- t: 1D numpy array at which to evaluate the points
- sigma_1: width of the 1st gaussian
- amp_1: amplitude of the 1st gaussian
- x0_1: shift in x of the 1st gaussian
- sigma_2: width of the 2nd gaussian
- amp_2: amplitude of the 2nd gaussian
- x0_2: shift in x of the 2nd gaussian
- y0: shift in y of the gaussian (shared for the two gaussians)
x = np.linspace(-10,10,100)
plt.figure()
plt.plot(x, sum_of_gaussian(x, sigma_1=1, amp_1=1, x0_1=0, 
                       sigma_2=2, amp_2=-0.5, x0_2=0, y0=0), label="sigma_1:1 amp_1:1 x0_1:0 sigma_2:2 amp_2:-0.5 x0_2:0")

plt.plot(x, sum_of_gaussian(x, sigma_1=1, amp_1=1, x0_1=0, 
                       sigma_2=2, amp_2=-0.5, x0_2=-2, y0=0), label="sigma_1:1 amp_1:1 x0_1:0 sigma_2:2 amp_2:-0.5 x0_2:-2")
_ = plt.legend()

gaussian_2D[source]

gaussian_2D(xz, sigma_x, sigma_z, amp, theta, x0, z0, y0)

Two dimensional Gaussian function params:

- xz: meshgrid of x and z coordinates at which to evaluate the points
- sigma_x: width of the gaussian
- sigma_z: height of the gaussian
- amp: amplitude of the gaussian
- theta: angle of the gaussian (in radian)
- x0: shift in x of the gaussian
- z0: shift in z of the gaussian
- y0: shift in y of the gaussian
x,y = np.meshgrid(np.linspace(-10,10,100),np.linspace(-10,10,100))
fig = plt.figure(figsize=(10,4))

ax = fig.add_subplot(121, projection='3d')
z = gaussian_2D((x,y), sigma_x=1, sigma_z=1, amp=1, theta=0, x0=0, z0=0, y0=0).reshape(100,100)
ax.plot_wireframe(x, y, z, rstride=4, cstride=4, label="sigma:1 amp:1")
_ = ax.legend()

ax = fig.add_subplot(122, projection='3d')
z = gaussian_2D((x,y), sigma_x=2, sigma_z=2, amp=1, theta=0, x0=0, z0=0, y0=0).reshape(100,100)
ax.plot_wireframe(x, y, z, rstride=4, cstride=4, label="sigma:2 amp:1", color='r')
_ = ax.legend()

sum_of_2D_gaussian[source]

sum_of_2D_gaussian(xz, sigma_x_1, sigma_z_1, amp_1, theta_1, x0_1, z0_1, sigma_x_2, sigma_z_2, amp_2, theta_2, x0_2, z0_2, y0)

Sum of two 2D Gaussian function. For the params, see gaussian_2D. However, both share the y0 parameter.

x,y = np.meshgrid(np.linspace(-10,10,100),np.linspace(-10,10,100))
fig = plt.figure(figsize=(5,4))

ax = fig.add_subplot(111, projection='3d')
z = sum_of_2D_gaussian((x,y), sigma_x_1=1, sigma_z_1=1, amp_1=1, theta_1=0, x0_1=0, z0_1=0,
                           sigma_x_2=2, sigma_z_2=2, amp_2=-0.5, theta_2=0, x0_2=0, z0_2=0, y0=0).reshape(100,100)
ax.plot_wireframe(x, y, z, rstride=3, cstride=3, label="sigma:1 amp:1")
_ = ax.legend()

exponential_decay[source]

exponential_decay(x, tau, baseline, amplitude)

Exponential decay function params:

- x: 1D numpy array at which to evaluate the points
- tau: decay speed
- baseline: baseline after the decay. Must be lower than the amplitude
- amplitude: Initial amplitude before the decay
x = np.linspace(0,10,100)
plt.figure()
plt.plot(x, exponential_decay(x, tau=1, baseline=0, amplitude=1), label="tau:1 baseline:0 amp:1")
plt.plot(x, exponential_decay(x, tau=2, baseline=0, amplitude=1), label="tau:2 baseline:0 amp:1")
plt.plot(x, exponential_decay(x, tau=1, baseline=.5, amplitude=1), label="tau:1 baseline:0.5 amp:1")
plt.plot(x, exponential_decay(x, tau=1, baseline=0, amplitude=2), label="tau:1 baseline:0 amp:2")

_ = plt.legend()

sin_exponent[source]

sin_exponent(x, amp, phi, freq, exp)

Sine raised to an exponent power. To the power 2, the sine is going between 0 and 1 while its period is halfed. Subsequent round powers have the effect of narowing the peaks and making the rest of the sine approach zero. (graph it to see it) params:

- x: 1D numpy array at which to evaluate the points
- amp: amplitude of the sine
- phi: phase of the sine in radian
- freq: frequency of the sine in Hz
- exp: Power the sine is raised to
x = np.linspace(0,5,100)
plt.figure(figsize=(10,4))
plt.subplot(1,2,1)
plt.plot(x, sin_exponent(x, amp=1, phi=0, freq=0.5, exp=1), label="freq:0.5 exp:1")
plt.plot(x, sin_exponent(x, amp=1, phi=0, freq=0.5, exp=2), label="freq:0.5 exp:2")
_ = plt.legend()

plt.subplot(1,2,2)
plt.plot(x, sin_exponent(x, amp=1, phi=0, freq=0.5, exp=2), label=    "phi:0      freq:1 exp:2")
plt.plot(x, sin_exponent(x, amp=1, phi=0, freq=0.5, exp=8), label=    "phi:0      freq:1 exp:8")
plt.plot(x, sin_exponent(x, amp=1, phi=np.pi, freq=0.5, exp=8), label="phi:3.14 freq:1 exp:8")
_ = plt.legend()

sinexp_sigm[source]

sinexp_sigm(x, sigma, x0, y0, amp, phi, freq, exp)

A sine exponent weighted by a sigmoid. For parameters, see sigmoid and sin_exponent. However the amp and y0 (baseline) parameters are shared between the two sub-functions

x = np.linspace(0,10,200)
plt.figure()
plt.plot(x, sinexp_sigm(x, sigma=1, x0=4, y0=0, amp=1, phi=0, freq=1, exp=2), label="phi=0")
plt.plot(x, sinexp_sigm(x, sigma=1, x0=4, y0=0, amp=1, phi=np.pi, freq=1, exp=2), label="phi=3.14")
_ = plt.legend()

Fitting functions

fit_sigmoid[source]

fit_sigmoid(nonlin, t=None)

Fit a sigmoid to a 1D input. Used to fit the nonlinearities of the cells responses. Also, the data must be increasing (from left to right). params:

- nonlin: 1D data to fit
- t     : timepoints of the data to fit. Same lenght as nonlin

return:

- A dictionary containing the parameters of the fit for a [`sigmoid`](/theonerig/modelling.html#sigmoid)
- A quality index of the fit (explained variance of the model)

fit_spatial_sta[source]

fit_spatial_sta(sta)

Fit a sum of 2D gaussian to a 2D matrix. Used to fit the STA obtained from checkerboard like stimulus. params:

- sta: 2D data to fit

return:

- A dictionary containing the parameters of the fit for a [`sum_of_2D_gaussian`](/theonerig/modelling.html#sum_of_2D_gaussian)
- A quality index of the fit (explained variance of the model)

fit_temporal_sta[source]

fit_temporal_sta(sta, frame_rate=60)

Fit a sum of gaussian to a 1D input. Used to fit the STA obtained from fullfield flicker like stimulus. params:

- sta: 1D data to fit

return:

- A dictionary containing the parameters of the fit for a [`diff_of_gaussian`](/theonerig/modelling.html#diff_of_gaussian). The parameter with _1 suffix
correspond to the closest gaussian to the zero timepoint.
- A quality index of the fit (explained variance of the model)

fit_chirp_am[source]

fit_chirp_am(cell_mean, start=420, stop=960, freq=1.5, frame_rate=60)

Fit a sinexp_sigm to the mean response of a cell to chirp_am stimulus. params:

- cell_mean: Average response of the cell to the chirp_am stimulus
- start:     Index in the cell mean where to start the fitting
- stop:      Index in the cell mean where to stop the fitting
- freq:      Frequency of the amplitude modulation of the stimulus in Hz

return :

- A dictionary containing the parameters of the fit for a [`sinexp_sigm`](/theonerig/modelling.html#sinexp_sigm)
- A quality index of the fit (explained variance of the model)

fit_chirp_freq_epoch[source]

fit_chirp_freq_epoch(cell_mean, start=360, freqs=[1.875, 3.75, 7.5, 15, 30], durations=[2, 2, 2, 1, 1], frame_rate=60)

Fit multiple sinexp_sigm to the mean response of a cell to chirp_freq_epoch stimulus. Each epoch is fitted independantly by a sin_exponent params:

- cell_mean: Average response of the cell to the chirp_am stimulus
- start:     Index in the cell mean where to start the fitting
- freqs:     Frequencies of each epoch in Hz
- durations: Duration of each epoch in seconds
- freq:      Frequency of the amplitude modulation of the stimulus

return :

- List of dictionary containing the parameters of each epoch fit for a [`sin_exponent`](/theonerig/modelling.html#sin_exponent)
- List of quality index of each epoch fit (explained variance of the model)

fit_transiency[source]

fit_transiency(pref_response, frame_rate=60)

Fit an exponential decay to a 1D input. Used to fit the response of cells to their preferred sustained stimulus. params:

- pref_response: 1D data to fit

return:

- A dictionary containing the parameters of the fit for a [`exponential_decay`](/theonerig/modelling.html#exponential_decay)
- A quality index of the fit (explained variance of the model)

repetition_quality_index[source]

repetition_quality_index(cell_response)

Return a quality index of cell response to a repeated stimulus. params:

- cell_response: response of a cell of shape (n_rep, time)

onoff_transient_index[source]

onoff_transient_index(cell_response, start_on=120, stop_on=240, start_off=240, stop_off=360)

Return both on-off and transient indexes of cell response. Transiency index is not always ideal to describe the cells response. Consider fitting an exponential decay with fit_transiency. params:

- cell_response: response of a cell of shape (time,...)
- start_on:  starting index of ON stimulation
- stop_on:   stop index of ON stimulation
- start_off: starting index of OFF stimulation
- stop_off:  stop index of OFF stimulation

return:

- ON-OFF index
- Transiency index