diff --git a/PointNICE/plt/pltutils.py b/PointNICE/plt/pltutils.py index f68c8d4..6406703 100644 --- a/PointNICE/plt/pltutils.py +++ b/PointNICE/plt/pltutils.py @@ -1,472 +1,510 @@ # -*- coding: utf-8 -*- # @Author: Theo Lemaire # @Date: 2017-08-23 14:55:37 # @Last Modified by: Theo Lemaire -# @Last Modified time: 2017-08-24 18:14:46 +# @Last Modified time: 2017-08-24 18:47:04 ''' Plotting utilities ''' import pickle import ntpath import re +import logging import inspect import tkinter as tk from tkinter import filedialog import numpy as np import matplotlib.pyplot as plt from matplotlib.patches import Rectangle from .. import channels from ..bls import BilayerSonophore from .pltvars import pltvars +# Get package logger +logger = logging.getLogger('PointNICE') # Define global variables neuron = None bls = None # Regular expression for input files rgxp = re.compile('([E,A])STIM_([A-Za-z]*)_(.*).pkl') +# Figure naming conventions +ESTIM_CW_title = '{} neuron: E-STIM {:.2f} mA/m2, {:.0f} ms)' +ASTIM_CW_title = '{} neuron: CW A-STIM {:.0f}kHz, {:.0f}kPa, {:.0f}ms' +ASTIM_PW_title = '{} neuron: PW A-STIM {:.0f}kHz, {:.0f}kPa, {:.0f}ms, {:.2f}kHz PRF, {:.0f}% DC' + class InteractiveLegend(object): """ Class defining an interactive matplotlib legend, where lines visibility can be toggled by simply clicking on the corresponding legend label. Other graphic objects can also be associated to the toggle of a specific line Adapted from: http://stackoverflow.com/questions/31410043/hiding-lines-after-showing-a-pyplot-figure """ def __init__(self, legend, aliases): self.legend = legend self.fig = legend.axes.figure self.lookup_artist, self.lookup_handle = self._build_lookups(legend) self._setup_connections() self.handles_aliases = aliases self.update() def _setup_connections(self): for artist in self.legend.texts + self.legend.legendHandles: artist.set_picker(10) # 10 points tolerance self.fig.canvas.mpl_connect('pick_event', self.on_pick) def _build_lookups(self, legend): ''' Method of the InteractiveLegend class building the legend lookups. ''' labels = [t.get_text() for t in legend.texts] handles = legend.legendHandles label2handle = dict(zip(labels, handles)) handle2text = dict(zip(handles, legend.texts)) lookup_artist = {} lookup_handle = {} for artist in legend.axes.get_children(): if artist.get_label() in labels: handle = label2handle[artist.get_label()] lookup_handle[artist] = handle lookup_artist[handle] = artist lookup_artist[handle2text[handle]] = artist lookup_handle.update(zip(handles, handles)) lookup_handle.update(zip(legend.texts, handles)) return lookup_artist, lookup_handle def on_pick(self, event): handle = event.artist if handle in self.lookup_artist: artist = self.lookup_artist[handle] artist.set_visible(not artist.get_visible()) self.update() def update(self): for artist in self.lookup_artist.values(): handle = self.lookup_handle[artist] if artist.get_visible(): handle.set_visible(True) if artist in self.handles_aliases: for al in self.handles_aliases[artist]: al.set_visible(True) else: handle.set_visible(False) if artist in self.handles_aliases: for al in self.handles_aliases[artist]: al.set_visible(False) self.fig.canvas.draw() def show(self): ''' showing the interactive legend ''' plt.show() def getPatchesLoc(t, states): ''' Determine the location of stimulus patches. :param t: simulation time vector (s). :param states: a vector of stimulation state (ON/OFF) at each instant in time. :return: 3-tuple with number of patches, timing of STIM-ON an STIM-OFF instants. ''' # Compute states derivatives and identify bounds indexes of pulses dstates = np.diff(states) ipatch_on = np.insert(np.where(dstates > 0.0)[0] + 1, 0, 0) ipatch_off = np.where(dstates < 0.0)[0] # Get time instants for pulses ON and OFF npatches = ipatch_on.size tpatch_on = t[ipatch_on] tpatch_off = t[ipatch_off] # return 3-tuple with #patches, pulse ON and pulse OFF instants return (npatches, tpatch_on, tpatch_off) def SaveFigDialog(dirname, filename): """ Open a FileSaveDialogBox to set the directory and name of the figure to be saved. The default directory and filename are given, and the default extension is ".pdf" :param dirname: default directory :param filename: default filename :return: full path to the chosen filename """ root = tk.Tk() root.withdraw() filename_out = filedialog.asksaveasfilename(defaultextension=".pdf", initialdir=dirname, initialfile=filename) return filename_out def plotComp(yvars, filepaths, fs=15, show_patches=True): ''' Compare profiles of several specific output variables of NICE simulations. :param yvars: list of variables names to extract and compare :param filepaths: list of full paths to output data files to be compared :param fs: labels fontsize :param show_patches: boolean indicating whether to indicate periods of stimulation with colored rectangular patches ''' nvars = len(yvars) # x and y variables plotting information t_plt = pltvars['t'] y_pltvars = [pltvars[key] for key in yvars] # Dictionary of neurons neurons = {} for _, obj in inspect.getmembers(channels): if inspect.isclass(obj) and isinstance(obj.name, str): neurons[obj.name] = obj # Initialize figure and axes if nvars == 1: _, ax = plt.subplots(figsize=(11, 4)) axes = [ax] else: _, axes = plt.subplots(nvars, 1, figsize=(11, min(3 * nvars, 9))) - labels = [ntpath.basename(fp)[4:-4].replace('_', ' ') for fp in filepaths] + + for i in range(nvars): ax = axes[i] pltvar = y_pltvars[i] if 'min' in pltvar and 'max' in pltvar: ax.set_ylim(pltvar['min'], pltvar['max']) if pltvar['unit']: ax.set_ylabel('${}\ ({})$'.format(pltvar['label'], pltvar['unit']), fontsize=fs) else: ax.set_ylabel('${}$'.format(pltvar['label']), fontsize=fs) if i < nvars - 1: ax.get_xaxis().set_ticklabels([]) else: ax.set_xlabel('${}\ ({})$'.format(t_plt['label'], t_plt['unit']), fontsize=fs) for tick in ax.xaxis.get_major_ticks(): tick.label.set_fontsize(fs) ax.locator_params(axis='y', nbins=2) for tick in ax.yaxis.get_major_ticks(): tick.label.set_fontsize(fs) # Loop through data files tstim_ref = 0.0 nstim = 0 j = 0 aliases = {} for filepath in filepaths: pkl_filename = ntpath.basename(filepath) # Retrieve neuron name mo = rgxp.fullmatch(pkl_filename) if not mo: print('Error: PKL file does not match regular expression pattern') quit() sim_type = mo.group(1) neuron_name = mo.group(2) # Load data print('Loading data from "' + pkl_filename + '"') with open(filepath, 'rb') as pkl_file: data = pickle.load(pkl_file) # Extract useful variables t = data['t'] states = data['states'] tstim = data['tstim'] # Initialize BLS and channels mechanism global neuron neuron = neurons[neuron_name]() neuron_states = [data[sn] for sn in neuron.states_names] # Initialize BLS if sim_type == 'A': global bls params = data['params'] Fdrive = data['Fdrive'] a = data['a'] d = data['d'] geom = {"a": a, "d": d} Qm0 = neuron.Cm0 * neuron.Vm0 * 1e-3 bls = BilayerSonophore(geom, params, Fdrive, neuron.Cm0, Qm0) # Get data of variables to plot vrs = [] for i in range(nvars): pltvar = y_pltvars[i] if 'alias' in pltvar: var = eval(pltvar['alias']) elif 'key' in pltvar: var = data[pltvar['key']] else: var = data[yvars[i]] vrs.append(var) # Determine patches location npatches, tpatch_on, tpatch_off = getPatchesLoc(t, states) # Adding onset to all signals if t_plt['onset'] > 0.0: t = np.insert(t + t_plt['onset'], 0, 0.0) for i in range(nvars): vrs[i] = np.insert(vrs[i], 0, vrs[i][0]) tpatch_on += t_plt['onset'] tpatch_off += t_plt['onset'] + # Legend label + if sim_type == 'E': + label = ESTIM_CW_title.format(neuron.name, data['Astim'], data['tstim'] * 1e3) + elif sim_type == 'A': + if data['DF'] == 1.0: + label = ASTIM_CW_title.format(neuron.name, Fdrive * 1e-3, + data['Adrive'] * 1e-3, data['tstim'] * 1e3) + else: + label = ASTIM_PW_title.format(neuron.name, Fdrive * 1e-3, + data['Adrive'] * 1e-3, data['tstim'] * 1e3, + data['PRF'] * 1e-3, data['DF'] * 1e2) + # Plotting handles = [axes[i].plot(t * t_plt['factor'], vrs[i] * y_pltvars[i]['factor'], - linewidth=2, label=labels[j]) for i in range(nvars)] - plt.tight_layout() + linewidth=2, label=label) for i in range(nvars)] if show_patches: k = 0 # stimulation patches for ax in axes: handle = handles[k] (ybottom, ytop) = ax.get_ylim() la = [] for i in range(npatches): la.append(ax.add_patch(Rectangle((tpatch_on[i] * t_plt['factor'], ybottom), (tpatch_off[i] - tpatch_on[i]) * t_plt['factor'], ytop - ybottom, color=handle[0].get_color(), alpha=0.2))) aliases[handle[0]] = la k += 1 if tstim != tstim_ref: if nstim == 0: nstim += 1 tstim_ref = tstim else: print('Warning: comparing different stimulation durations') j += 1 + plt.tight_layout() iLegends = [] for k in range(nvars): axes[k].legend(loc='upper left', fontsize=fs) iLegends.append(InteractiveLegend(axes[k].legend_, aliases)) plt.show() def plotBatch(vars_dict, directory, filepaths, plt_show=True, plt_save=False, - ask_before_save=1, fig_ext='png', tag='fig', fs=15, lw=4, show_patches=True): + ask_before_save=1, fig_ext='png', tag='fig', fs=15, lw=4, title=True, + show_patches=True): ''' Plot a figure with profiles of several specific NICE output variables, for several NICE simulations. :param vars_dict: dict of lists of variables names to extract and plot together :param positions: subplot indexes of each variable :param filepaths: list of full paths to output data files to be compared :param fs: labels fontsize :param show_patches: boolean indicating whether to indicate periods of stimulation with colored rectangular patches ''' labels = list(vars_dict.keys()) naxes = len(vars_dict) # x and y variables plotting information t_plt = pltvars['t'] # Dictionary of neurons neurons = {} for _, obj in inspect.getmembers(channels): if inspect.isclass(obj) and isinstance(obj.name, str): neurons[obj.name] = obj # Loop through data files for filepath in filepaths: # Get code from file name pkl_filename = ntpath.basename(filepath) filecode = pkl_filename[0:-4] # Retrieve neuron name mo = rgxp.fullmatch(pkl_filename) if not mo: print('Error: PKL file does not match regular expression pattern') quit() sim_type = mo.group(1) neuron_name = mo.group(2) + assert sim_type in ['A', 'E'], 'invalid stimulation type (should be "ESTIM" or "ASTIM")' # Load data print('Loading data from "' + pkl_filename + '"') with open(filepath, 'rb') as pkl_file: data = pickle.load(pkl_file) # Extract variables print('Extracting variables') t = data['t'] states = data['states'] nsamples = t.size # Initialize channel mechanism global neuron neuron = neurons[neuron_name]() neuron_states = [data[sn] for sn in neuron.states_names] # Initialize BLS if sim_type == 'A': global bls params = data['params'] Fdrive = data['Fdrive'] a = data['a'] d = data['d'] geom = {"a": a, "d": d} Qm0 = neuron.Cm0 * neuron.Vm0 * 1e-3 bls = BilayerSonophore(geom, params, Fdrive, neuron.Cm0, Qm0) # Determine patches location npatches, tpatch_on, tpatch_off = getPatchesLoc(t, states) # Adding onset to time vector and patches if t_plt['onset'] > 0.0: t = np.insert(t + t_plt['onset'], 0, 0.0) tpatch_on += t_plt['onset'] tpatch_off += t_plt['onset'] # Plotting if naxes == 1: _, ax = plt.subplots(figsize=(11, 4)) axes = [ax] else: _, axes = plt.subplots(naxes, 1, figsize=(11, min(3 * naxes, 9))) for i in range(naxes): ax = axes[i] ax_pltvars = [pltvars[j] for j in vars_dict[labels[i]]] nvars = len(ax_pltvars) # X-axis if i < naxes - 1: ax.get_xaxis().set_ticklabels([]) else: ax.set_xlabel('${}\ ({})$'.format(t_plt['label'], t_plt['unit']), fontsize=fs) for tick in ax.xaxis.get_major_ticks(): tick.label.set_fontsize(fs) # Y-axis if ax_pltvars[0]['unit']: ax.set_ylabel('${}\ ({})$'.format(labels[i], ax_pltvars[0]['unit']), fontsize=fs) else: ax.set_ylabel('${}$'.format(labels[i]), fontsize=fs) if 'min' in ax_pltvars[0] and 'max' in ax_pltvars[0]: ax_min = min([ap['min'] for ap in ax_pltvars]) ax_max = max([ap['max'] for ap in ax_pltvars]) ax.set_ylim(ax_min, ax_max) ax.locator_params(axis='y', nbins=2) for tick in ax.yaxis.get_major_ticks(): tick.label.set_fontsize(fs) # Time series icolor = 0 for j in range(nvars): # Extract variable pltvar = ax_pltvars[j] if 'alias' in pltvar: var = eval(pltvar['alias']) elif 'key' in pltvar: var = data[pltvar['key']] elif 'constant' in pltvar: var = eval(pltvar['constant']) * np.ones(nsamples) else: var = data[vars_dict[labels[i]][j]] if t_plt['onset'] > 0.0: var = np.insert(var, 0, var[0]) # Plot variable if 'constant' in pltvar: ax.plot(t * t_plt['factor'], var * pltvar['factor'], '--', c='black', lw=lw, label='${}$'.format(pltvar['label'])) else: ax.plot(t * t_plt['factor'], var * pltvar['factor'], c='C{}'.format(icolor), lw=lw, label='${}$'.format(pltvar['label'])) icolor += 1 # Patches if show_patches == 1: (ybottom, ytop) = ax.get_ylim() for j in range(npatches): ax.add_patch(Rectangle((tpatch_on[j] * t_plt['factor'], ybottom), (tpatch_off[j] - tpatch_on[j]) * t_plt['factor'], ytop - ybottom, color='#8A8A8A', alpha=0.1)) # Legend if nvars > 1: ax.legend(fontsize=fs, loc=7, ncol=nvars // 4 + 1) + + # Title + if title: + if sim_type == 'E': + fig_title = ESTIM_CW_title.format(neuron.name, data['Astim'], data['tstim'] * 1e3) + elif sim_type == 'A': + if data['DF'] == 1.0: + fig_title = ASTIM_CW_title.format(neuron.name, Fdrive * 1e-3, + data['Adrive'] * 1e-3, data['tstim'] * 1e3) + else: + fig_title = ASTIM_PW_title.format(neuron.name, Fdrive * 1e-3, + data['Adrive'] * 1e-3, data['tstim'] * 1e3, + data['PRF'] * 1e-3, data['DF'] * 1e2) + axes[0].set_title(fig_title, fontsize=fs) + plt.tight_layout() # Save figure if needed (automatic or checked) if plt_save == 1: if ask_before_save == 1: plt_filename = SaveFigDialog(directory, '{}_{}.{}'.format(filecode, tag, fig_ext)) else: plt_filename = '{}/{}_{}.{}'.format(directory, filecode, tag, fig_ext) if plt_filename: plt.savefig(plt_filename) print('Saving figure as "{}"'.format(plt_filename)) plt.close() # Show all plots if needed if plt_show == 1: plt.show() diff --git a/sim/ASTIM_batch.py b/sim/ASTIM_batch.py index 014d2c8..8a36c70 100644 --- a/sim/ASTIM_batch.py +++ b/sim/ASTIM_batch.py @@ -1,97 +1,97 @@ #!/usr/bin/env python # -*- coding: utf-8 -*- # @Author: Theo Lemaire # @Date: 2017-02-13 18:16:09 # @Email: theo.lemaire@epfl.ch # @Last Modified by: Theo Lemaire -# @Last Modified time: 2017-08-24 18:16:23 +# @Last Modified time: 2017-08-24 18:37:50 """ Run batch acoustic simulations of the NICE model. """ import os import logging import numpy as np from PointNICE.solvers import checkBatchLog, runAStimBatch from PointNICE.channels import * from PointNICE.utils import LoadParams from PointNICE.plt import plotBatch # Set logging options logging.basicConfig(format='%(asctime)s %(message)s', datefmt='%d/%m/%Y %H:%M:%S:') logger = logging.getLogger('PointNICE') logger.setLevel(logging.DEBUG) # BLS parameters bls_params = LoadParams() # Geometry of BLS structure a = 32e-9 # in-plane radius (m) d = 0.0e-6 # embedding tissue thickness (m) geom = {"a": a, "d": d} # Channels mechanisms -neurons = [ThalamoCortical()] +neurons = [CorticalLTS()] # Stimulation parameters stim_params = { 'freqs': [3.5e5], # Hz 'amps': [100e3], # Pa 'durations': [50e-3], # s 'PRFs': [1e2], # Hz - 'DFs': [1.0] + 'DFs': [0.5] } stim_params['offsets'] = [30e-3] * len(stim_params['durations']) # s # Select output directory try: (batch_dir, log_filepath) = checkBatchLog('A-STIM') except AssertionError as err: logger.error(err) quit() # Run A-STIM batch pkl_filepaths = runAStimBatch(batch_dir, log_filepath, neurons, bls_params, geom, stim_params) pkl_dir, _ = os.path.split(pkl_filepaths[0]) vars_RS_FS = { 'Z': ['Z'], 'Q_m': ['Qm'], 'i_{Na}\ kin.': ['m', 'h'], 'i_K\ kin.': ['n'], 'i_M\ kin.': ['p'], 'I': ['iNa', 'iK', 'iM', 'iL', 'iNet'] } vars_LTS = { 'Z': ['Z'], 'Q_m': ['Qm'], 'i_{Na}\ kin.': ['m', 'h'], 'i_K\ kin.': ['n'], 'i_M\ kin.': ['p'], 'i_T\ kin.': ['s', 'u'], 'I': ['iNa', 'iK', 'iM', 'iT', 'iL', 'iNet'] } vars_RE = { 'Z': ['Z'], 'Q_m': ['Qm'], 'i_{Na}\ kin.': ['m', 'h'], 'i_K\ kin.': ['n'], 'i_{TS}\ kin.': ['s', 'u'], 'I': ['iNa', 'iK', 'iTs', 'iL', 'iNet'] } vars_TC = { 'Z': ['Z'], 'Q_m': ['Qm'], 'i_{Na}\ kin.': ['m', 'h'], 'i_K\ kin.': ['n'], 'i_{T}\ kin.': ['s', 'u'], 'i_{H}\ kin.': ['O', 'OL', 'O + 2OL'], 'I': ['iNa', 'iK', 'iT', 'iH', 'iKL', 'iL', 'iNet'] } -plotBatch(vars_TC, pkl_dir, pkl_filepaths) +plotBatch(vars_LTS, pkl_dir, pkl_filepaths) diff --git a/sim/ESTIM_sim_batch.py b/sim/ESTIM_sim_batch.py deleted file mode 100644 index bdc74ff..0000000 --- a/sim/ESTIM_sim_batch.py +++ /dev/null @@ -1,315 +0,0 @@ -#!/usr/bin/env python -# -*- coding: utf-8 -*- -# @Author: Theo Lemaire -# @Date: 2016-10-11 20:35:38 -# @Email: theo.lemaire@epfl.ch -# @Last Modified by: Theo Lemaire -# @Last Modified time: 2017-08-22 18:08:12 - -""" Run simulations of the HH system with injected electric current, -and plot resulting dynamics. """ - -import matplotlib.pyplot as plt -import matplotlib.patches as patches - -from PointNICE.solvers import SolverElec -from PointNICE.channels import * - - -# -------------- SIMULATION ----------------- - - -# Create channels mechanism -neuron = LeechTouch() - -print('initial states:') -print('Vm0 = {:.2f}'.format(neuron.Vm0)) -for i in range(len(neuron.states_names)): - if neuron.states_names[i] == 'C_Ca': - print('{}0 = {:.3f} uM'.format(neuron.states_names[i], neuron.states0[i] * 1e6)) - else: - print('{}0 = {:.2f}'.format(neuron.states_names[i], neuron.states0[i])) - - -# Set pulse parameters -tstim = 1.5 # s -toffset = 0.5 # s -Astim = 4.1 # mA/m2 - -show_currents = False - -# Run simulation -print('stimulating {} neuron ({:.2f} mA/m2, {:.0f} ms)'.format(neuron.name, Astim, tstim * 1e3)) -solver = SolverElec() -(t, y) = solver.runSim(neuron, Astim, tstim, toffset) - - -# -------------- VARIABLES SEPARATION ----------------- - - -# Membrane potential and states -Vm = y[:, 0] -states = y[:, 1:].T - -# Final states -statesf = y[-1, 1:] -print('final states:') -print('Vmf = {:.2f}'.format(Vm[-1])) -for i in range(len(neuron.states_names)): - if len(neuron.states_names[i]) == 1: # channel state - print('{}f = {:.2f}'.format(neuron.states_names[i], statesf[i])) - else: # other state - print('{}f = {:f}'.format(neuron.states_names[i], statesf[i])) - - -# Leakage current and net current -iL = neuron.currL(Vm) -iNet = neuron.currNet(Vm, states) - -# Sodium and Potassium gating dynamics and currents -m = y[:, 1] -h = y[:, 2] -n = y[:, 3] -iNa = neuron.currNa(m, h, Vm) -iK = neuron.currK(n, Vm) - - -corticals = ['RS', 'FS', 'LTS'] -thalamics = ['RE', 'TC'] -leeches = ['LeechT'] - -# Cortical neurons -if neuron.name in corticals: - p = y[:, 4] - iM = neuron.currM(p, Vm) - - # Special case: LTS neuron - if neuron.name == 'LTS': - s = y[:, 5] - u = y[:, 6] - iCa = neuron.currCa(s, u, Vm) - - -# Thalamic neurons -if neuron.name in thalamics: - s = y[:, 4] - u = y[:, 5] - iCa = neuron.currCa(s, u, Vm) - - # Special case: TC neuron - if neuron.name == 'TC': - O = y[:, 6] - C = y[:, 7] - P0 = y[:, 8] - C_Ca = y[:, 9] - OL = 1 - O - C - P1 = 1 - P0 - Ih = neuron.currH(O, C, Vm) - IKL = neuron.currKL(Vm) - -# Leech neurons -if neuron.name in leeches: - s = y[:, 4] - C_Na = y[:, 5] - A_Na = y[:, 6] - C_Ca = y[:, 7] - A_Ca = y[:, 8] - iCa = neuron.currCa(s, Vm) - iPumpNa = neuron.currPumpNa(C_Na, Vm) - iKCa = neuron.currKCa(C_Ca, Vm) - - -# -------------- PLOTTING ----------------- - -fs = 12 -if neuron.name == 'TC': - naxes = 7 -if neuron.name in ['LTS', 'RE']: - naxes = 5 -if neuron.name in ['RS', 'FS']: - naxes = 4 -if neuron.name in leeches: - naxes = 7 - -if not show_currents: - naxes -= 1 - -height = 5.5 -if neuron.name == 'TC': - height = 7 -fig, axes = plt.subplots(naxes, 1, figsize=(10, height)) - -# Membrane potential -i = 0 -ax = axes[i] -ax.plot(t * 1e3, Vm, linewidth=2) -ax.set_ylabel('$V_m\ (mV)$', fontsize=fs) -if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) -ax.locator_params(axis='y', nbins=2) -for item in ax.get_yticklabels(): - item.set_fontsize(fs) -(ybottom, ytop) = ax.get_ylim() -ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) - -# iNa dynamics -i += 1 -ax = axes[i] -ax.set_ylim([-0.1, 1.1]) -ax.set_ylabel('$Na^+ \ kin.$', fontsize=fs) -ax.plot(t * 1e3, m, color='blue', linewidth=2, label='$m$') -ax.plot(t * 1e3, h, color='red', linewidth=2, label='$h$') -ax.plot(t * 1e3, m**2 * h, '--', color='black', linewidth=2, label='$m^2h$') -(ybottom, ytop) = ax.get_ylim() -ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) -ax.legend(fontsize=fs, loc=7) -if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) -ax.locator_params(axis='y', nbins=2) -for item in ax.get_yticklabels(): - item.set_fontsize(fs) - - -# iK & iM dynamics -i += 1 -ax = axes[i] -ax.set_ylim([-0.1, 1.1]) -ax.set_ylabel('$K^+ \ kin.$', fontsize=fs) -ax.plot(t * 1e3, n, color='#734d26', linewidth=2, label='$n$') -if neuron.name in ['RS', 'FS', 'LTS']: - ax.plot(t * 1e3, p, color='#660099', linewidth=2, label='$p$') -(ybottom, ytop) = ax.get_ylim() -ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) -ax.legend(fontsize=fs, loc=7) -if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) -ax.locator_params(axis='y', nbins=2) -for item in ax.get_yticklabels(): - item.set_fontsize(fs) - -# iCa dynamics -if neuron.name in ['LTS', 'RE', 'TC', 'LeechT']: - i += 1 - ax = axes[i] - ax.set_ylim([-0.1, 1.1]) - ax.set_ylabel('$Ca^{2+} \ kin.$', fontsize=fs) - ax.plot(t * 1e3, s, color='#2d862d', linewidth=2, label='$s$') - if neuron.name in ['LTS', 'RE', 'TC']: - ax.plot(t * 1e3, u, color='#e68a00', linewidth=2, label='$u$') - ax.plot(t * 1e3, s**2 * u, '--', color='black', linewidth=2, label='$s^2u$') - (ybottom, ytop) = ax.get_ylim() - ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) - ax.legend(fontsize=fs, loc=7) - if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) - ax.locator_params(axis='y', nbins=2) - for item in ax.get_yticklabels(): - item.set_fontsize(fs) - - -# iH dynamics -if neuron.name == 'TC': - i += 1 - ax = axes[i] - ax.set_ylim([-0.1, 2.1]) - ax.set_ylabel('$i_H\ kin.$', fontsize=fs) - # ax.plot(t * 1e3, C, linewidth=2, label='$C$') - ax.plot(t * 1e3, O, linewidth=2, label='$O$') - ax.plot(t * 1e3, OL, linewidth=2, label='$O_L$') - ax.plot(t * 1e3, O + 2 * OL, '--', color='black', linewidth=2, label='$O + 2O_L$') - (ybottom, ytop) = ax.get_ylim() - ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) - ax.legend(fontsize=fs, ncol=2, loc=7) - if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) - ax.locator_params(axis='y', nbins=2) - for item in ax.get_yticklabels(): - item.set_fontsize(fs) - -# submembrane [Ca2+] dynamics -if neuron.name in ['TC', 'LeechT']: - i += 1 - ax = axes[i] - if neuron.name == 'TC': - ax.set_ylabel('$[Ca^{2+}_i]\ (uM)$', fontsize=fs) - ax.plot(t * 1e3, C_Ca * 1e6, linewidth=2, label='$[Ca^{2+}_i]$') - if neuron.name == 'LeechT': - ax.set_ylabel('$[Ca^{2+}_i]\ (arb.)$', fontsize=fs) - ax.plot(t * 1e3, C_Ca, linewidth=2, label='$[Ca^{2+}_i]$') - ax.plot(t * 1e3, A_Ca, linewidth=2, label='$A_{Ca}$') - ax.legend(fontsize=fs, loc=7) - (ybottom, ytop) = ax.get_ylim() - ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) - if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) - ax.locator_params(axis='y', nbins=2) - for item in ax.get_yticklabels(): - item.set_fontsize(fs) - -# submembrane [Na+] dynamics -if neuron.name == 'LeechT': - i += 1 - ax = axes[i] - ax.set_ylabel('$[Na^{+}_i]\ (arb.)$', fontsize=fs) - ax.plot(t * 1e3, C_Na, linewidth=2, label='$[Na^{+}_i]$') - ax.plot(t * 1e3, A_Na, linewidth=2, label='$A_{Na}$') - (ybottom, ytop) = ax.get_ylim() - ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) - ax.legend(fontsize=fs, loc=7) - if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) - ax.locator_params(axis='y', nbins=2) - for item in ax.get_yticklabels(): - item.set_fontsize(fs) - -# currents -if show_currents: - i += 1 - ax = axes[i] - ax.set_ylabel('$I\ (A/m^2)$', fontsize=fs) - ax.set_xlabel('$time\ (ms)$', fontsize=fs) - ax.plot(t * 1e3, iNa * 1e-3, linewidth=2, label='$i_{Na}$') - ax.plot(t * 1e3, iK * 1e-3, linewidth=2, label='$i_K$') - if neuron.name in ['RS', 'FS', 'LTS']: - ax.plot(t * 1e3, iM * 1e-3, linewidth=2, label='$i_M$') - if neuron.name in ['LTS', 'TC', 'LeechT']: - ax.plot(t * 1e3, iCa * 1e-3, linewidth=2, label='$i_{T}$') - if neuron.name == 'RE': - ax.plot(t * 1e3, iCa * 1e-3, linewidth=2, label='$i_{TS}$') - if neuron.name == 'TC': - ax.plot(t * 1e3, Ih * 1e-3, linewidth=2, label='$i_{H}$') - ax.plot(t * 1e3, IKL * 1e-3, linewidth=2, label='$i_{KL}$') - if neuron.name == 'LeechT': - ax.plot(t * 1e3, iKCa * 1e-3, linewidth=2, label='$i_{K,Ca}$') - ax.plot(t * 1e3, iPumpNa * 1e-3, linewidth=2, label='$i_{Na\ pump}$') - ax.plot(t * 1e3, iL * 1e-3, linewidth=2, label='$i_L$') - ax.plot(t * 1e3, iNet * 1e-3, '--', linewidth=2, color='black', label='$i_{Net}$') - ax.legend(fontsize=fs, ncol=2, loc=7) - ax.locator_params(axis='y', nbins=2) - for item in ax.get_yticklabels(): - item.set_fontsize(fs) - if i < naxes - 1: - ax.get_xaxis().set_ticklabels([]) - (ybottom, ytop) = ax.get_ylim() - ax.add_patch(patches.Rectangle((0.0, ybottom), tstim * 1e3, ytop - ybottom, - color='#8A8A8A', alpha=0.1)) - - -axes[-1].set_xlabel('$time\ (ms)$', fontsize=fs) -for item in axes[-1].get_xticklabels(): - item.set_fontsize(fs) - -if tstim > 0.0: - title = '{} neuron ({:.2f} mA/m2, {:.0f} ms)'.format(neuron.name, Astim, tstim * 1e3) -else: - title = '{} neuron (free, {:.0f} ms)'.format(neuron.name, toffset * 1e3) -fig.suptitle(title, fontsize=fs) - -plt.show()