Activation thresholds with arbitrary settings
This script shows how to use the axon_AP_threshold()
function to evaluate axon thresholds with various stimulation waveforms.
[1]:
import nrv
import matplotlib.pyplot as plt
# axon def
y = 0 # axon y position, in [um]
z = 0 # axon z position, in [um]
d = 10 # axon diameter, in [um]
L = 50_000 # axon length, along x axis, in [um]
axon1 = nrv.myelinated(y,z,d,L,T=37,rec='nodes',dt=0.001)
#electrode
n_node = len(axon1.x_nodes)
x_elec = axon1.x_nodes[n_node//2]
y_elec = 1000
z_elec = 0
PSA = nrv.point_source_electrode(x_elec,y_elec,z_elec)
stim1 = nrv.stimulus() #dummy stim
### define extra cellular stimulation
extra_stim = nrv.stimulation('endoneurium_bhadra')
extra_stim.add_electrode(PSA, stim1)
axon1.attach_extracellular_stimulation(extra_stim)
#Update function for cathodic (monopolar) pulse
def cath_pulse_update(axon,amp, pw, start_p = 1, elec_id = 0):
stim_1 = nrv.stimulus()
stim_1.pulse(start=start_p, duration=pw, value = -amp)
axon.change_stimulus_from_electrode(elec_id, stim_1)
#Update function for biphasic pulse
def biphasic_pulse_update(axon,amp, pw, start_p = 1, elec_id = 0,t_inter=50e-3):
stim_1 = nrv.stimulus()
stim_1.biphasic_pulse(start = start_p, s_anod=amp,t_stim=pw,s_cathod=amp,t_inter=t_inter)
axon.change_stimulus_from_electrode(elec_id, stim_1)
#Update function for cathodic sine pulse
def cath_sine_pulse_update(axon,amp, pw, start_p = 1, elec_id = 0):
stim_1 = nrv.stimulus()
freq = 1/(2*pw)
stim_1.sinus(start_p, pw, amp, freq, offset=0, phase=180, dt=0)
axon.change_stimulus_from_electrode(elec_id, stim_1)
#parameters for the waveforms
arg_stim = {'pw':50e-3, 'elec_id':0, 'start_p':1}
max_amp = 300 #maximum search boundary
threshold_pulse = nrv.axon_AP_threshold(axon = axon1,amp_max = max_amp,
update_func = cath_pulse_update, args_update=arg_stim)
threshold_biphasic = nrv.axon_AP_threshold(axon = axon1,amp_max = max_amp,
update_func = biphasic_pulse_update, args_update=arg_stim)
threshold_sine = nrv.axon_AP_threshold(axon = axon1,amp_max = max_amp,
update_func = cath_sine_pulse_update, args_update=arg_stim)
NRV INFO: Iteration 1, Amp is 300µA ...
NRV INFO: Iteration 1, Amp is 300µA (100.0%)... AP Detected! (in 0.932s)
NRV INFO: Iteration 2, Amp is 0µA (100.0%)... AP Not Detected! (in 0.934s)
NRV INFO: Iteration 3, Amp is 150.0µA (100.0%)... AP Not Detected! (in 0.928s)
NRV INFO: Iteration 4, Amp is 225.0µA (33.33%)... AP Detected! (in 0.924s)
NRV INFO: Iteration 5, Amp is 187.5µA (16.67%)... AP Detected! (in 0.921s)
NRV INFO: Iteration 6, Amp is 168.75µA (10.0%)... AP Not Detected! (in 0.916s)
NRV INFO: Iteration 7, Amp is 178.12µA (5.26%)... AP Not Detected! (in 0.923s)
NRV INFO: Iteration 8, Amp is 182.81µA (2.56%)... AP Not Detected! (in 0.92s)
NRV INFO: Iteration 9, Amp is 185.16µA (1.27%)... AP Not Detected! (in 0.927s)
NRV INFO: Iteration 10, Amp is 186.33µA (0.63%)... AP Not Detected! (in 0.921s)
NRV INFO: Activation threshold is 186.91µA (0.63%), found in 10 iterations (9.38s).
NRV INFO: Iteration 1, Amp is 300µA ...
NRV INFO: Iteration 1, Amp is 300µA (100.0%)... AP Detected! (in 0.921s)
NRV INFO: Iteration 2, Amp is 0µA (100.0%)... AP Not Detected! (in 0.92s)
NRV INFO: Iteration 3, Amp is 150.0µA (100.0%)... AP Not Detected! (in 0.916s)
NRV INFO: Iteration 4, Amp is 225.0µA (33.33%)... AP Detected! (in 0.915s)
NRV INFO: Iteration 5, Amp is 187.5µA (16.67%)... AP Not Detected! (in 0.925s)
NRV INFO: Iteration 6, Amp is 206.25µA (9.09%)... AP Detected! (in 0.934s)
NRV INFO: Iteration 7, Amp is 196.88µA (4.55%)... AP Detected! (in 0.926s)
NRV INFO: Iteration 8, Amp is 192.19µA (2.38%)... AP Not Detected! (in 0.931s)
NRV INFO: Iteration 9, Amp is 194.53µA (1.2%)... AP Not Detected! (in 0.926s)
NRV INFO: Iteration 10, Amp is 195.7µA (0.6%)... AP Not Detected! (in 0.922s)
NRV INFO: Activation threshold is 196.29µA (0.6%), found in 10 iterations (9.37s).
NRV INFO: Iteration 1, Amp is 300µA ...
NRV INFO: Iteration 1, Amp is 300µA (100.0%)... AP Detected! (in 0.94s)
NRV INFO: Iteration 2, Amp is 0µA (100.0%)... AP Not Detected! (in 0.939s)
NRV INFO: Iteration 3, Amp is 150.0µA (100.0%)... AP Not Detected! (in 0.941s)
NRV INFO: Iteration 4, Amp is 225.0µA (33.33%)... AP Not Detected! (in 0.939s)
NRV INFO: Iteration 5, Amp is 262.5µA (14.29%)... AP Not Detected! (in 0.933s)
NRV INFO: Iteration 6, Amp is 281.25µA (6.67%)... AP Detected! (in 0.935s)
NRV INFO: Iteration 7, Amp is 271.88µA (3.33%)... AP Detected! (in 0.985s)
NRV INFO: Iteration 8, Amp is 267.19µA (1.72%)... AP Detected! (in 0.946s)
NRV INFO: Iteration 9, Amp is 264.84µA (0.88%)... AP Not Detected! (in 0.937s)
NRV INFO: Activation threshold is 266.02µA (0.88%), found in 9 iterations (8.66s).
[42]:
plt.plot([0.25], [threshold_pulse], "s", markersize=10, label="pulse", color="c")
plt.plot([1.25], [threshold_biphasic], "s", markersize=10, label="biphasic", color="orange")
plt.plot([2.25], [threshold_sine], "s", markersize=10, label="sine", color="limegreen")
plt.fill_between([0, 0.5], [threshold_pulse, threshold_pulse], 0, color="c")
plt.fill_between([1,1.5], [threshold_biphasic, threshold_biphasic], 0, color="orange")
plt.fill_between([2, 2.5], [threshold_sine, threshold_sine], 0, color="limegreen")
plt.xticks([])
plt.xlabel("stimulus type")
plt.ylabel("amplitude $\mu m$")
plt.legend()
[42]:
<matplotlib.legend.Legend at 0x30c8aa600>