pyfibers package

Initializer for core PyFibers code.

enum pyfibers.BisectionMean(value)

Bases: str, Enum

Mean type used during bisection search phase of finding threshold.

Member Type:

str

Valid values are as follows:

GEOMETRIC = <BisectionMean.GEOMETRIC: 'geometric'>
ARITHMETIC = <BisectionMean.ARITHMETIC: 'arithmetic'>
enum pyfibers.BoundsSearchMode(value)

Bases: str, Enum

Modes for adjusting bounds in the bounds search phase of finding threshold.

Member Type:

str

Valid values are as follows:

PERCENT_INCREMENT = <BoundsSearchMode.PERCENT_INCREMENT: 'percent'>
ABSOLUTE_INCREMENT = <BoundsSearchMode.ABSOLUTE_INCREMENT: 'absolute'>
class pyfibers.Fiber(fiber_model, diameter, temperature=37, passive_end_nodes=True, is_3d=False)

Bases: object

Base class for model fibers.

The Fiber class provides functionality for constructing, configuring, and simulating fiber models. It encapsulates key methods for:

  • Generating sections specified by a fiber model subclass

  • Recording membrane voltage, current, and gating variables

  • Calculating extracellular potentials and single fiber action potentials

  • Measuring conduction velocity along the fiber

  • Handling 3D or 1D fiber geometry

_create_sections(function_list)

Create and connect NEURON sections for each node or internode in the fiber.

The provided function_list starts with a function for a node, followed by each internodal section in order. Each node is optionally converted to a passive node if it is within the range of Fiber.passive_end_nodes.

Parameters:

function_list (list[Callable]) – A list of functions that each return a new NEURON h.Section.

Return type:

Fiber

Returns:

The updated Fiber instance.

_make_passive(node)

Convert a node section to passive by removing all active mechanisms.

For more info, see Implementations of Fiber Models.

Parameters:

node (Section) – The node h.Section to be made passive.

Return type:

Section

Returns:

The modified section with a passive mechanism inserted.

Raises:

AssertionError – If the node’s name does not contain ‘passive’.

add_intrinsic_activity(loc=0.1, loc_index=None, avg_interval=1, num_stims=1, start_time=1, noise=0, synapse_tau=0.1, synapse_reversal_potential=0, netcon_weight=0.1)

Add a spike-generating synapse to the fiber for intrinsic activity.

The synapse is generated via a NetStim (spike generator), which is connected to an ExpSyn (exponential synapse) on the chosen node. A NetCon object links them together, injecting an exponentially decaying current upon each spike event.

Parameters:
  • loc (float) – Normalized location along the fiber where the synapse is placed ([0,1]).

  • loc_index (int) – Alternatively, specify an integer index of the node.

  • avg_interval (float) – Average interval between NetStim spikes (ms).

  • num_stims (int) – Number of spikes to deliver.

  • start_time (float) – Time to start delivering spikes (ms).

  • noise (float) – Noise parameter for spike intervals (0 = regular, 1 = Poisson).

  • synapse_tau (float) – Time constant (ms) for synaptic current decay.

  • synapse_reversal_potential (float) – Reversal potential (mV) of the synapse.

  • netcon_weight (float) – Weight of the NetCon between the spike generator and synapse.

Raises:

AssertionError – If neither loc nor loc_index is specified, or if both are specified.

Return type:

None

apcounts(thresh=-30)

Create NEURON APCount objects at each node to detect action potentials.

Parameters:

thresh (float) – Threshold voltage (mV) for an action potential.

Return type:

None

static calculate_periaxonal_current(from_sec, to_sec, vext_from, vext_to)

Compute the periaxonal current between two compartments for myelinated models.

Parameters:
  • from_sec (Section) – The NEURON h.Section from which current flows.

  • to_sec (Section) – The NEURON h.Section receiving current.

  • vext_from (float) – Extracellular (periaxonal) potential at from_sec (mV).

  • vext_to (float) – Extracellular (periaxonal) potential at to_sec (mV).

Return type:

float

Returns:

The periaxonal current in mA.

generate(function_list, n_nodes=None, n_sections=None, length=None, enforce_odd_nodecount=True)

Build the fiber model sections in NEURON according to a specified generation strategy.

This method either uses the desired number of nodes (n_nodes), the number of sections (n_sections), or determines the total number of sections from the fiber length and the spacing Fiber.delta_z.

Parameters:
  • function_list (list[Callable]) – List of callable functions that each create a section. Typically structured as [node_func, myelin_func].

  • n_nodes (int) – Number of nodes of Ranvier.

  • n_sections (int) – Total number of sections in the fiber.

  • length (float) – Total length of the fiber (µm). Overrides n_sections if given.

  • enforce_odd_nodecount (bool) – If True, ensures that the fiber has an odd number of nodes.

Return type:

Fiber

Returns:

The updated Fiber instance after generation.

Raises:

AssertionError – If the computed number of sections does not align with the function_list-based pattern.

is_3d()

Check if the fiber is using 3D coordinates.

Return type:

bool

Returns:

True if 3D, otherwise False.

loc(loc, target='nodes')

Retrieve a node or section at a given normalized position along the fiber.

Returns either the node or section nearest to loc * (len - 1) from the fiber. Note that the indexing is performed on the nodes or sections list. This means that the number represents the proportion along the list of nodes or sections, not necessarily along the physical fiber (though these are generally the same).

Parameters:
  • loc (float) – Location in the range [0, 1].

  • target (str) – Specifies whether to retrieve from 'nodes' or 'sections'.

Raises:

AssertionError – If loc is not in [0, 1] or if target is not 'nodes' or 'sections'.

Return type:

Section

Returns:

The chosen node or section as a h.Section.

loc_index(loc, target='nodes')

Convert a normalized location [0, 1] into an integer index for nodes or sections.

Parameters:
  • loc (float) – Location in the fiber (from 0 to 1).

  • target (str) – Indicates whether to index into 'nodes' or 'sections'.

Raises:

AssertionError – If loc is not in [0, 1] or if target is invalid.

Return type:

int

Returns:

The integer index corresponding to the node or section.

measure_cv(start=0.25, end=0.75, tolerance=0.005)

Estimate conduction velocity (m/s) by measuring AP times at two points (start and end).

This method calculates the conduction velocity by comparing the action potential times at two specified normalized locations (using NEURON indexing). It also checks for linear conduction between the two points, within a specified tolerance.

Parameters:
  • start (float) – Starting position for conduction velocity measurement (from 0 to 1).

  • end (float) – Ending position for conduction velocity measurement (from 0 to 1).

  • tolerance (float) – Tolerance (ms) for checking linearity of AP times.

Raises:
  • ValueError – If conduction is not approximately linear between start and end.

  • AssertionError – If no APs are detected at one or both of the measurement nodes.

Return type:

float

Returns:

Conduction velocity in meters per second (m/s).

membrane_currents(downsample=1)

Compute membrane currents at each section for each time point in the simulation.

Uses the methods described in Pena et. al 2024: http://dx.doi.org/10.1371/journal.pcbi.1011833

For myelinated fibers, this calculation includes periaxonal currents between adjacent sections (based on h.Section.xraxial <neuron.Section.xraxial>). The result is a matrix of shape: [num_timepoints, num_sections].

This method returns a tuple consisting of: 1. A 2D array (time steps x number of sections) of membrane currents in mA. 2. The array of time points corresponding to those currents (downsampled by the specified factor).

Parameters:

downsample (int) – Factor to reduce the temporal resolution (e.g., downsample=2 takes every 2nd time step).

Return type:

ndarray

Returns:

A tuple (i_membrane_matrix, downsampled_time). The matrix contains total currents (mA). The time array contains the corresponding simulation times (ms).

Raises:

RuntimeError – If membrane currents or extracellular potentials were not recorded properly.

nodebuilder(ind, node_type)

Create a generic node of Ranvier.

This method sets the node length, diameter, and inserts the extracellular mechanism, but does not add active ion channel mechanisms. Subclasses or external code should add the relevant channels to define the node’s electrophysiological behavior.

Parameters:
  • ind (int) – Index of this node in the overall fiber construction.

  • node_type (str) – A string identifying the node as ‘active’ or ‘passive’.

Return type:

Section

Returns:

The newly created node as a NEURON h.Section.

point_source_potentials(x, y, z, i0, sigma, inplace=False)

Compute extracellular potentials at the fiber’s coordinates due to a point source stimulus.

See also

Documentation on extracellular potentials <extracellular potentials.md> in PyFibers.

Parameters:
  • x (float) – x-coordinate of the source in µm.

  • y (float) – y-coordinate of the source in µm.

  • z (float) – z-coordinate of the source in µm.

  • i0 (float) – Magnitude of the point-source current (mA).

  • sigma (float | tuple) – Conductivity (S/m). A float for isotropic or a tuple (sigma_x, sigma_y, sigma_z) for anisotropic.

  • inplace (bool) – If True, update Fiber.potentials in-place.

Return type:

ndarray

Returns:

Extracellular potentials at each fiber coordinate, in mV.

record_gating(**kwargs)

Record gating parameters (ion channel states) from axon nodes.

The gating variables must be declared in Fiber.gating_variables within the fiber model class.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

dict[str, list[h.Vector | None]]

Returns:

A dictionary mapping gating variable names to lists of recorded NEURON Vector objects.

Raises:

AssertionError – If Fiber.gating_variables is empty.

record_im(**kwargs)

Record membrane current (nA) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane current.

record_sfap(rec_potentials, downsample=1)

Compute the SFAP time course at a given electrode location.

Parameters:
  • rec_potentials (list | ndarray) – 1D array of precomputed potentials (mV) at each fiber section due to the electrode placement.

  • downsample (int) – Downsampling factor for the time vector (applies to the current matrix).

Return type:

ndarray

Returns:

A tuple (sfap_trace, downsampled_time) where sfap_trace is the computed single-fiber action potential in microvolts (µV) and downsampled_time is the corresponding time array (ms).

record_values(ref_attr, allsec=False, indices=None, allow_missing=True, recording_dt=None, recording_tvec=None)

Record NEURON variable references (e.g. membrane voltage) along the fiber.

Note that recording_dt and recording_tvec are mutually exclusive. If both are None, the variable is recorded at every simulation timestep. For more info, see the NEURON docs: https://nrn.readthedocs.io/en/latest/python/programming/math/vector.html#Vector.record

Parameters:
  • ref_attr (str) – The NEURON attribute to record (e.g. '_ref_v').

  • allsec (bool) – If True, record from sections (including nodes). Otherwise, only record from nodes.

  • indices (list[int]) – Specific indices to record from (if None, record from all).

  • allow_missing (bool) – If True, allows missing attributes without raising an error (returns None).

  • recording_dt (float) – The time step [ms] for recording the values (separate from simulation dt). Should be larger than the simulation dt.

  • recording_tvec (h.Vector) –

    A NEURON Vector of time points at which to record the values (ms). Note that the user MUST keep this Vector in memory for the duration of the simulation. This means you must assign it to a variable that is not overwritten or deleted. For example, to record at time points 0, 1, 2, and 3 ms:

    recording_tvec = h.Vector([0, 1, 2, 3])  # store times in a Vector
    fiber.record_values("_ref_v", recording_tvec=recording_tvec)  # pass Vector to record
    stimulation.find_threshold(fiber)  # run the simulation
    plt.plot(recording_tvec, fiber.vm[0])  # plot the recorded values
    

Raises:

ValueError – If indices is an empty list.

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects or None (if allow_missing=True and the requested attribute is missing).

record_vext()

Record extracellular potential (mV) from each section along the fiber.

Return type:

list[Vector]

Returns:

A list of NEURON Vector objects recording extracellular potential.

record_vm(**kwargs)

Record membrane voltage (mV) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane voltage.

resample_potentials(potentials, potential_coords, center=False, inplace=False)

Use linear interpolation to resample external potentials onto the fiber’s coordinate system (1D).

This is used when extracellular potentials are calculated from an external source, such as a finite element model. The potentials provided by the user should be sampled at high resolution along the fiber’s path and provided alongside the corresponding arc-length coordinates.

If center=True, both the input coordinates and the fiber’s coordinates will be shifted such that their midpoints align.

Parameters:
  • potentials (ndarray) – 1D array of external potential values.

  • potential_coords (ndarray) – 1D array of coordinates corresponding to potentials.

  • center (bool) – If True, center the potentials around the midpoint of each domain.

  • inplace (bool) – If True, update Fiber.potentials with the resampled values.

Return type:

ndarray

Returns:

Interpolated potential values aligned with Fiber.longitudinal_coordinates.

Raises:

AssertionError – If input array sizes or monotonicity checks fail.

resample_potentials_3d(potentials, potential_coords, center=False, inplace=False)

Interpolate external potentials onto a 3D fiber coordinate system.

A wrapper around Fiber.resample_potentials() that handles 3D coordinates by first computing the arc length of the provided coordinate array. As with the 1D version, this method is used to resample external potentials (e.g., from a finite element model) onto the fiber.

At present, this does not check that the input coordinates lie exactly along the 3D fiber path. Therefore, it is recommended to use the same coordinates to construct the fiber as to use here. Alternatively, you can create a 1D fiber and calculate the coordinate arc lengths. For more information, see Supplying Extracellular Potentials.

Parameters:
  • potentials (ndarray) – 1D array of external potential values.

  • potential_coords (ndarray) – 2D array of shape (N, 3) representing the (x, y, z) coordinates where the potentials are measured or computed.

  • center (bool) – If True, center the potentials around the midpoint of each domain.

  • inplace (bool) – If True, update Fiber.potentials with the resampled values.

Return type:

ndarray

Returns:

Interpolated potential values aligned with the fiber’s 3D arc-length coordinates.

Raises:

AssertionError – If called on a non-3D fiber or if input coordinate shapes are invalid.

set_save_gating(**kwargs)

Record gating parameters (ion channel states) from axon nodes.

The gating variables must be declared in Fiber.gating_variables within the fiber model class.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

dict[str, list[h.Vector | None]]

Returns:

A dictionary mapping gating variable names to lists of recorded NEURON Vector objects.

Raises:

AssertionError – If Fiber.gating_variables is empty.

set_save_im(**kwargs)

Record membrane current (nA) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane current.

set_save_vext()

Record extracellular potential (mV) from each section along the fiber.

Return type:

list[Vector]

Returns:

A list of NEURON Vector objects recording extracellular potential.

set_save_vm(**kwargs)

Record membrane voltage (mV) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane voltage.

set_xyz(x=0, y=0, z=0)

Assign new (x, y, z) shifts to a straight (1D) fiber.

The fiber is assumed to be along the z-axis initially, with x=0 and y=0. This method sets the x and y coordinates to the specified values and shifts the z coordinate by the given amount.

Parameters:
  • x (float) – Shift in the x-direction (µm).

  • y (float) – Shift in the y-direction (µm).

  • z (float) – Shift in the z-direction (µm).

Raises:

AssertionError – If this fiber is a 3D fiber (since this method is for 1D only).

Return type:

None

static sfap(current_matrix, potentials)

Compute the Single-Fiber Action Potential (SFAP) by multiplying currents with recording potentials.

Parameters:
  • current_matrix (ndarray) – 2D array of shape [timepoints, sections], containing currents in mA.

  • potentials (ndarray) – 1D array of potentials (mV) at each section, length = number of sections.

Raises:

AssertionError – If the number of columns in the current matrix does not match the length of potentials.

Return type:

ndarray

Returns:

The computed SFAP in microvolts (µV).

enum pyfibers.FiberModel(value)

Bases: NoAliasEnum

An enumeration.

Valid values are as follows:

MRG_DISCRETE = <FiberModel.MRG_DISCRETE: <class 'pyfibers.models.mrg.MRGFiber'>>
MRG_INTERPOLATION = <FiberModel.MRG_INTERPOLATION: <class 'pyfibers.models.mrg.MRGFiber'>>
RATTAY = <FiberModel.RATTAY: <class 'pyfibers.models.rattay.RattayFiber'>>
SCHILD94 = <FiberModel.SCHILD94: <class 'pyfibers.models.schild.SchildFiber'>>
SCHILD97 = <FiberModel.SCHILD97: <class 'pyfibers.models.schild.SchildFiber'>>
SMALL_MRG_INTERPOLATION = <FiberModel.SMALL_MRG_INTERPOLATION: <class 'pyfibers.models.mrg.MRGFiber'>>
SUNDT = <FiberModel.SUNDT: <class 'pyfibers.models.sundt.SundtFiber'>>
THIO_AUTONOMIC = <FiberModel.THIO_AUTONOMIC: <class 'pyfibers.models.thio.ThioFiber'>>
THIO_CUTANEOUS = <FiberModel.THIO_CUTANEOUS: <class 'pyfibers.models.thio.ThioFiber'>>
TIGERHOLM = <FiberModel.TIGERHOLM: <class 'pyfibers.models.tigerholm.TigerholmFiber'>>
class pyfibers.IntraStim(dt=0.001, tstop=50, t_init_ss=-200, dt_init_ss=5, istim_ind=None, istim_loc=None, clamp_kws=None)

Bases: Stimulation

Manage intracellular stimulation of model fibers.

This class extends the Stimulation class to provide intracellular stimulation capabilities. The intracellular stimulation is managed via a custom h.trainIClamp mechanism. This mechanism allows for repeated square pulses of current to be injected into a fiber. Its arguments are provided as clamp_kws when creating an instance of this class.

Example Usage

# Assuming you have already created a fiber "my_fiber"
clamp_kws = {
    "delay": 1,  # ms
    "pw": 0.1,  # ms
    "dur": 50,  # ms
    "freq": 100,  # Hz
    "amp": 1,  # nA, recommended to set 1 for scaling purposes
}
istim_ind = 0  # Stimulate the first section
dt, tstop = 0.001, 50
istim = IntraStim(dt=dt, tstop=tstop, istim_ind=istim_ind, clamp_kws=clamp_kws)

# Run a simulation with a given stimulation amplitude
stimamp = 2
istim.run_sim(2, my_fiber)

# Run a threshold search
istim.find_threshold(my_fiber, condition="activation")
_add_istim(fiber)

Create instance of h.trainIClamp for intracellular stimulation.

This method is not called by the user, see :class:IntraStim. Note that trainIClamp is a mod file included in this package. It is an extension of NEURON’s built-in IClamp that allows repeated square pulses.

Parameters:

fiber (Fiber) – The Fiber object to attach intracellular stimulation to.

Return type:

IntraStim

Returns:

The Stimulation instance (self).

run_sim(stimamp, fiber, ap_detect_location=0.9, exit_func=<function IntraStim.<lambda>>, exit_func_interval=100, exit_func_kws=None, use_exit_t=False, fail_on_end_excitation=True, ap_detect_threshold=-30)

Run a simulation for a single stimulation amplitude.

Parameters:
  • stimamp (float) – Amplitude to be applied to extracellular stimulation - Should be a single float for one source - If stimamp is a single float and there are multiple sources, the same stimamp is applied to all sources - If stimamp is a list of floats, each float is applied to the corresponding source

  • fiber (Fiber) – The Fiber to be stimulated.

  • ap_detect_location (float) – Location to detect action potentials (percent along fiber)

  • exit_func (Callable) – Function to call to check if simulation should be exited

  • exit_func_interval (int) – Interval to call exit_func

  • exit_func_kws (dict) – Keyword arguments to pass to exit_func

  • use_exit_t (bool) – If True, use the time returned by exit_func as the simulation end time

  • fail_on_end_excitation (bool) – Behavior for end excitation detection - if True, raise error if end excitation is detected - if False, continue simulation if end excitation is detected - if None, do not check for end excitation

  • ap_detect_threshold (float) – Threshold for detecting action potentials (default: -30 mV)

Raises:

RuntimeError – If NaNs are detected in fiber potentials

Return type:

tuple[int, float]

Returns:

Number of detected APs and time of last detected AP.

class pyfibers.ScaledStim(waveform, dt=0.001, tstop=50, t_init_ss=-200, dt_init_ss=5, pad_waveform=True, truncate_waveform=True)

Bases: Stimulation

Manage extracellular stimulation of model fibers.

This class takes one or more waveforms, each of which is expected to match one corresponding set of fiber potentials (i.e., one source) in the fiber being stimulated. Therefore, if you have N potential sets on the fiber, you must provide N waveforms, each describing the time course of stimulation for that source.

The waveform can be either:
  • A callable that accepts a single float argument for the stimulation time (in ms) and returns

the waveform value at that time.
  • A list of such callables with length N if there are N sources

See also

For more information on defining potentials, see Supplying Extracellular Potentials.

** Example Usage **

# Assuming you have already created a fiber "my_fiber"
# Add potentials to the fiber
my_fiber.potentials = electrical_potentials_array


# Create function which takes in a time and returns a waveform value
def my_square_pulse(t: float) -> float:
    if t > 0.1 and t <= 0.2:  # if time is during the pulse
        return 1  # on
    else:  # if time is outside the pulse
        return 0  # off


# Same waveform as above but using scipy.interpolate.interp1d
my_square_pulse = interp1d(
    [0, 0.1, 0.2, 50],  # start, on, off, stop, in ms
    [0, 1, 0, 0],  # waveform values at those times
    kind="previous",
)

# Create a :class:`ScaledStim` object
dt, tstop = 0.001, 50
stim = ScaledStim(waveform=my_square_pulse, dt=dt, tstop=tstop)

# Run the simulation
stim.run_sim(-1, my_fiber)

# Calculate threshold
stim.find_threshold(my_fiber, condition="activation")
run_sim(stimamp, fiber, ap_detect_location=0.9, exit_func=<function ScaledStim.<lambda>>, exit_func_interval=100, exit_func_kws=None, use_exit_t=False, fail_on_end_excitation=True, ap_detect_threshold=-30)

Run a simulation with a given stimulus amplitude(s).

See also

For more information on the underlying math, see the Algorithms docs page.

Parameters:
  • stimamp (float | list[float]) – Amplitude to scale the product of extracellular potentials and waveform. - Should be a single float for one source - If stimamp is a single float and there are multiple sources, the same stimamp is applied to all sources - If stimamp is a list of floats, each float is applied to the corresponding source

  • fiber (Fiber) – The Fiber object to stimulate.

  • ap_detect_location (float) – Normalized location in [0,1] to check for APs.

  • exit_func (Callable[..., bool]) – Callback to check if the simulation can be ended early (e.g., upon detection of an AP).

  • exit_func_interval (int) – How often (in time steps) to call exit_func.

  • exit_func_kws (dict | None) – Additional arguments for exit_func.

  • use_exit_t (bool) – If True, simulation will stop after self._exit_t (if set).

  • fail_on_end_excitation (bool | None) – Behavior for end excitation detection - If True, raise error if end excitation is detected - If False, continue simulation if end excitation is detected - If None, do not check for end excitation

  • ap_detect_threshold (float) – Threshold for detecting action potentials (default: -30 mV)

Raises:

RuntimeError – If NaNs are detected in membrane potentials or if required setup (e.g., istim) is missing.

Return type:

tuple[int, float]

Returns:

Tuple (number_of_APs, time_of_last_AP).

class pyfibers.Stimulation(dt=0.001, tstop=50, t_init_ss=-200, dt_init_ss=5, custom_run_sim=None)

Bases: object

Manage stimulation of NEURON simulations.

Provides methods to configure time stepping, run simulations, and perform threshold search. This class is not meant to be used as is; subclasses should override the Stimulation.run_sim() method. For example, the ScaledStim subclass provides extracellular stimulation capabilities.

For more details on using this class to write custom stimulation routines, see Custom Simulation Code.

static ap_checker(fiber, ap_detect_location=0.9, precision=3, check_all_apc=True)

Check how many action potentials occurred at a specified node location.

Parameters:
  • fiber (Fiber) – The Fiber object to evaluate for APs.

  • ap_detect_location (float) – Normalized location along the fiber in [0,1] to check for APs.

  • precision (int) – Decimal places to round the detected AP time.

  • check_all_apc (bool) – If True, raise a warning if APs occur elsewhere but not at the detect location.

Return type:

tuple[int, float | None]

Returns:

A tuple (num_aps, last_ap_time). If no APs are detected, last_ap_time is None.

Raises:

AssertionError – If the detected AP time is non-positive.

static end_excitation_checker(fiber, multi_site_check=True, fail_on_end_excitation=True)

Check for end-excitation.

Determines activation sites by finding local minima in each node’s AP time. If an AP is detected near in or adjacent to the passive end nodes, raise an error or issue a warning based on fail_on_end_excitation.

Parameters:
  • fiber (Fiber) – The Fiber object to check.

  • multi_site_check (bool) – If True, warn if multiple activation sites are detected.

  • fail_on_end_excitation (bool | None) – Controls handling of end-excitation: - True: Raise RuntimeError if end-excitation is detected. - False: Only warn if end-excitation is detected. - None: Skip the check entirely.

Return type:

bool

Returns:

True if end excitation is detected, False otherwise.

Raises:

RuntimeError – If end excitation is detected and fail_on_end_excitation is True.

find_threshold(fiber, condition=ThresholdCondition.ACTIVATION, bounds_search_mode=BoundsSearchMode.PERCENT_INCREMENT, bounds_search_step=10, termination_mode=TerminationMode.PERCENT_DIFFERENCE, termination_tolerance=1, stimamp_top=-1, stimamp_bottom=-0.01, max_iterations=50, exit_t_shift=5, bisection_mean=BisectionMean.ARITHMETIC, block_delay=0, thresh_num_aps=1, silent=False, **kwargs)

Perform a bisection search to find the threshold stimulus amplitude.

This method first locates bounds where one amplitude is subthreshold and

another is suprathreshold (bounds search phase). Then, it repeatedly narrows the bounds until they converge based on the specified termination mode and tolerance (bisection search phase). Ideally, the initial bounds should be chosen such that stimamp_top is supra-threshold and stimamp_bottom is sub-threshold.

Note that enums (ThresholdCondition, BoundsSearchMode,

TerminationMode, BisectionMean) can be provided as an enum member (e.g. ThresholdCondition.ACTIVATION) or as the member’s string value (e.g. “activation”).

See also

For more details on the threshold search process, see Algorithms in PyFibers.

Parameters:
  • fiber (Fiber) – Instance of pyfibers.fiber.Fiber to apply stimulation to.

  • condition (ThresholdCondition) – The threshold condition (ThresholdCondition.ACTIVATION or ThresholdCondition.BLOCK).

  • bounds_search_mode (BoundsSearchMode) – The bounds search mode (BoundsSearchMode.PERCENT_INCREMENT or BoundsSearchMode.ABSOLUTE_INCREMENT).

  • bounds_search_step (float) – The iterative increase/decrease of the upper/lower bound during bounds search - if bounds_search_mode is “percent” this is the percentage increase/decrease - if bounds_search_mode is “absolute” this is the absolute increase/decrease

  • termination_mode (TerminationMode) – The termination mode (TerminationMode.PERCENT_DIFFERENCE or TerminationMode.ABSOLUTE_DIFFERENCE).

  • termination_tolerance (float) – Difference between upper and lower bounds that indicates convergence - absolute difference if termination_mode is “absolute” - percentage difference if termination_mode is “percent”

  • stimamp_top (float) – Initial upper-bound stimulus amplitude to test.

  • stimamp_bottom (float) – Initial lower-bound stimulus amplitude to test.

  • max_iterations (int) – Maximum attempts to find bounding amplitudes before bisection.

  • exit_t_shift (float) – Extra time (ms) after an AP is detected, beyond which the simulation can be cut short.

  • bisection_mean (BisectionMean) – The bisection mean type (BisectionMean.ARITHMETIC or BisectionMean.GEOMETRIC).

  • block_delay (float) – Time (ms) after start to check for a blocked AP, used in block searches.

  • thresh_num_aps (int) – Number of action potentials for threshold search - if activation, suprathreshold requires detected aps >= thresh_num_aps - if block, suprathreshold requires detected aps < thresh_num_aps

  • silent (bool) – If True, suppress print statements for the search process.

  • kwargs – Additional arguments passed to the run_sim method.

Return type:

tuple[float, tuple[int, float | None]]

Returns:

A tuple (threshold_amplitude, (num_detected_aps, last_detected_ap_time)).

Raises:

RuntimeError – If contradictory bounding conditions occur or if the search fails to converge.

pre_run_setup(fiber, ap_detect_threshold=-30)

Prepare the simulation environment before running.

This method sets the temperature, initializes the membrane potential, configures AP detection, and optionally balances the membrane currents for certain fiber models (e.g., Tigerholm). It also applies any intracellular stimulation parameters if provided.

Parameters:
  • fiber (Fiber) – The Fiber object for which the simulation will be configured.

  • ap_detect_threshold (float) – The voltage threshold for detecting action potentials (mV).

Return type:

None

run_sim(*args, **kwargs)

Run a simulation using either a custom_run_sim method or a subclass override.

Return type:

tuple[int, float | None]

supra_exit(fiber, ap_detect_location, thresh_num_aps=1)

Determine if simulation can be exited early, for activation threshold searches only.

Parameters:
  • fiber (Fiber) – The Fiber object to check for an action potential.

  • ap_detect_location (float) – Normalized location to check in [0,1].

  • thresh_num_aps (int) – Number of APs required to consider it suprathreshold.

Return type:

bool

Returns:

True if the specified number of APs has occurred at or before this time.

static threshold_checker(fiber, block=False, ap_detect_location=0.9, block_delay=0, thresh_num_aps=1, check_all_apc=True)

Determine whether a stimulation was above or below threshold, for activation or block.

Parameters:
  • fiber (Fiber) – The Fiber object to evaluate.

  • ap_detect_location (float) – Normalized location in [0,1] where APs are detected.

  • block (bool) – If True, check for block threshold; otherwise, check for activation threshold.

  • block_delay (float) – Time after simulation start to check for block (ms).

  • thresh_num_aps (int) – Number of APs that constitutes a suprathreshold response.

  • check_all_apc (bool) – Passed to Stimulation.ap_checker() for additional warnings.

Return type:

bool

Returns:

True if stimulation is suprathreshold; False if subthreshold.

Raises:
threshsim(stimamp, fiber, condition=ThresholdCondition.ACTIVATION, block_delay=0, thresh_num_aps=1, **kwargs)

Run a single stimulation trial at a given amplitude and check for threshold.

Parameters:
Return type:

tuple[bool, tuple[int, float | None]]

Returns:

A tuple (is_suprathreshold, (num_aps, last_ap_time)).

enum pyfibers.TerminationMode(value)

Bases: str, Enum

Modes for determining when to terminate bisection search phase of finding threshold.

Member Type:

str

Valid values are as follows:

PERCENT_DIFFERENCE = <TerminationMode.PERCENT_DIFFERENCE: 'percent'>
ABSOLUTE_DIFFERENCE = <TerminationMode.ABSOLUTE_DIFFERENCE: 'absolute'>
enum pyfibers.ThresholdCondition(value)

Bases: str, Enum

Different threshold search conditions.

Member Type:

str

Valid values are as follows:

ACTIVATION = <ThresholdCondition.ACTIVATION: 'activation'>
BLOCK = <ThresholdCondition.BLOCK: 'block'>
pyfibers.build_fiber(fiber_model, diameter, length=None, n_sections=None, n_nodes=None, enforce_odd_nodecount=True, **kwargs)

Generate a 1D (straight) fiber model in NEURON.

This function creates a model fiber as an instance of the Fiber class using the specific subclass specified from the FiberModel enumerator. with user-specified diameter and length (from one of: number of sections, number of nodes, or length in microns). Additional keyword arguments are forwarded to the fiber model class constructor.

By default, the first section of the fiber is located at the origin (0, 0, 0), and the fiber extends along the z-axis in the positive direction. To change the fiber’s location, the method Fiber.set_xyz() can be used to translate the fiber along the x, y, or z axes. To create fibers along a custom path in 3D space, use build_fiber_3d() instead.

Parameters:
  • fiber_model (FiberModel) – A FiberModel enumerator specifying the type of fiber to instantiate.

  • diameter (float) – The fiber diameter in micrometers (µm).

  • length (float) – The total length of the fiber in micrometers (µm), if defining by length.

  • n_sections (int) – The total number of sections for discretizing the fiber, if defining by sections.

  • n_nodes (int) – The total number of nodes along the fiber, if defining by nodes.

  • enforce_odd_nodecount (bool) – If True, ensure that the number of nodes is odd.

  • kwargs – Additional arguments forwarded to the underlying fiber model class.

Raises:

ValueError – If more than one among length, n_sections, or n_nodes is specified.

Return type:

Fiber

Returns:

A Fiber class instance.

Example Usage

from PyFibers import build_fiber, FiberModel

fiber = build_fiber(fiber_model=FiberModel.MRG_INTERPOLATION, diameter=10, n_nodes=25)
pyfibers.build_fiber_3d(fiber_model, diameter, path_coordinates, shift=0, shift_ratio=None, center_shift=False, **kwargs)

Generate a 3D fiber model in NEURON based on a specified path.

This function calculates the fiber’s length from the user-supplied path_coordinates and uses it internally to instantiate a 3D fiber model. The coordinates are a 2D numpy array of shape (number_of_points, 3), where each row represents a point in 3D space (x, y, z).

The fiber model will be created by repeating sections along the path until no more nodes can be added without exceeding the path length. By default, the center of the first section is placed at the origin (0, 0, 0), and the fiber extends along the path.

Parameters:
  • fiber_model (FiberModel) – A FiberModel enumerator specifying the type of fiber to instantiate.

  • diameter (float) – The fiber diameter in micrometers (µm).

  • path_coordinates (ndarray) – A numpy array of shape (N, 3) specifying the 3D coordinates (x, y, z) of the fiber path.

  • shift (float) – A shift in microns to apply to the fiber coordinates.

  • shift_ratio (float) – Ratio of the internodal length to shift the fiber coordinates.

  • center_shift (bool) – If True, center the fiber before applying the shift.

  • kwargs – Additional arguments forwarded to the underlying fiber model class.

Raises:

ValueError – If path_coordinates is not provided, or if n_sections, n_nodes, or length is specified (these are invalid in 3D mode).

Return type:

Fiber

Returns:

A fully instantiated 3D fiber model Fiber instance.

Example:

import numpy as np
from PyFibers import build_fiber_3d, FiberModel

coords = np.array(
    [
        [0.0, 0.0, 0.0],
        [0.0, 1.0, 3.0],
        [0.0, 2.0, 7.0],
        # ...
    ]
)
fiber = build_fiber_3d(
    fiber_model=FiberModel.MRG_INTERPOLATION, diameter=10, path_coordinates=coords
)
print(fiber)

Subpackages

Submodules

pyfibers.compile module

Install script for PyFibers.

pyfibers.compile.main()

Compile NEURON MOD files.

Raises:

RuntimeError – If nrnivmodl is not found or fails.

Return type:

None

pyfibers.fiber module

Defines the Fiber class and helper functions for building fiber models.

This module provides functionality for building and simulating both 1D and 3D fiber models in the NEURON environment.

class pyfibers.fiber.Fiber(fiber_model, diameter, temperature=37, passive_end_nodes=True, is_3d=False)

Bases: object

Base class for model fibers.

The Fiber class provides functionality for constructing, configuring, and simulating fiber models. It encapsulates key methods for:

  • Generating sections specified by a fiber model subclass

  • Recording membrane voltage, current, and gating variables

  • Calculating extracellular potentials and single fiber action potentials

  • Measuring conduction velocity along the fiber

  • Handling 3D or 1D fiber geometry

_create_sections(function_list)

Create and connect NEURON sections for each node or internode in the fiber.

The provided function_list starts with a function for a node, followed by each internodal section in order. Each node is optionally converted to a passive node if it is within the range of Fiber.passive_end_nodes.

Parameters:

function_list (list[Callable]) – A list of functions that each return a new NEURON h.Section.

Return type:

Fiber

Returns:

The updated Fiber instance.

_make_passive(node)

Convert a node section to passive by removing all active mechanisms.

For more info, see Implementations of Fiber Models.

Parameters:

node (Section) – The node h.Section to be made passive.

Return type:

Section

Returns:

The modified section with a passive mechanism inserted.

Raises:

AssertionError – If the node’s name does not contain ‘passive’.

add_intrinsic_activity(loc=0.1, loc_index=None, avg_interval=1, num_stims=1, start_time=1, noise=0, synapse_tau=0.1, synapse_reversal_potential=0, netcon_weight=0.1)

Add a spike-generating synapse to the fiber for intrinsic activity.

The synapse is generated via a NetStim (spike generator), which is connected to an ExpSyn (exponential synapse) on the chosen node. A NetCon object links them together, injecting an exponentially decaying current upon each spike event.

Parameters:
  • loc (float) – Normalized location along the fiber where the synapse is placed ([0,1]).

  • loc_index (int) – Alternatively, specify an integer index of the node.

  • avg_interval (float) – Average interval between NetStim spikes (ms).

  • num_stims (int) – Number of spikes to deliver.

  • start_time (float) – Time to start delivering spikes (ms).

  • noise (float) – Noise parameter for spike intervals (0 = regular, 1 = Poisson).

  • synapse_tau (float) – Time constant (ms) for synaptic current decay.

  • synapse_reversal_potential (float) – Reversal potential (mV) of the synapse.

  • netcon_weight (float) – Weight of the NetCon between the spike generator and synapse.

Raises:

AssertionError – If neither loc nor loc_index is specified, or if both are specified.

Return type:

None

apcounts(thresh=-30)

Create NEURON APCount objects at each node to detect action potentials.

Parameters:

thresh (float) – Threshold voltage (mV) for an action potential.

Return type:

None

static calculate_periaxonal_current(from_sec, to_sec, vext_from, vext_to)

Compute the periaxonal current between two compartments for myelinated models.

Parameters:
  • from_sec (Section) – The NEURON h.Section from which current flows.

  • to_sec (Section) – The NEURON h.Section receiving current.

  • vext_from (float) – Extracellular (periaxonal) potential at from_sec (mV).

  • vext_to (float) – Extracellular (periaxonal) potential at to_sec (mV).

Return type:

float

Returns:

The periaxonal current in mA.

generate(function_list, n_nodes=None, n_sections=None, length=None, enforce_odd_nodecount=True)

Build the fiber model sections in NEURON according to a specified generation strategy.

This method either uses the desired number of nodes (n_nodes), the number of sections (n_sections), or determines the total number of sections from the fiber length and the spacing Fiber.delta_z.

Parameters:
  • function_list (list[Callable]) – List of callable functions that each create a section. Typically structured as [node_func, myelin_func].

  • n_nodes (int) – Number of nodes of Ranvier.

  • n_sections (int) – Total number of sections in the fiber.

  • length (float) – Total length of the fiber (µm). Overrides n_sections if given.

  • enforce_odd_nodecount (bool) – If True, ensures that the fiber has an odd number of nodes.

Return type:

Fiber

Returns:

The updated Fiber instance after generation.

Raises:

AssertionError – If the computed number of sections does not align with the function_list-based pattern.

is_3d()

Check if the fiber is using 3D coordinates.

Return type:

bool

Returns:

True if 3D, otherwise False.

loc(loc, target='nodes')

Retrieve a node or section at a given normalized position along the fiber.

Returns either the node or section nearest to loc * (len - 1) from the fiber. Note that the indexing is performed on the nodes or sections list. This means that the number represents the proportion along the list of nodes or sections, not necessarily along the physical fiber (though these are generally the same).

Parameters:
  • loc (float) – Location in the range [0, 1].

  • target (str) – Specifies whether to retrieve from 'nodes' or 'sections'.

Raises:

AssertionError – If loc is not in [0, 1] or if target is not 'nodes' or 'sections'.

Return type:

Section

Returns:

The chosen node or section as a h.Section.

loc_index(loc, target='nodes')

Convert a normalized location [0, 1] into an integer index for nodes or sections.

Parameters:
  • loc (float) – Location in the fiber (from 0 to 1).

  • target (str) – Indicates whether to index into 'nodes' or 'sections'.

Raises:

AssertionError – If loc is not in [0, 1] or if target is invalid.

Return type:

int

Returns:

The integer index corresponding to the node or section.

measure_cv(start=0.25, end=0.75, tolerance=0.005)

Estimate conduction velocity (m/s) by measuring AP times at two points (start and end).

This method calculates the conduction velocity by comparing the action potential times at two specified normalized locations (using NEURON indexing). It also checks for linear conduction between the two points, within a specified tolerance.

Parameters:
  • start (float) – Starting position for conduction velocity measurement (from 0 to 1).

  • end (float) – Ending position for conduction velocity measurement (from 0 to 1).

  • tolerance (float) – Tolerance (ms) for checking linearity of AP times.

Raises:
  • ValueError – If conduction is not approximately linear between start and end.

  • AssertionError – If no APs are detected at one or both of the measurement nodes.

Return type:

float

Returns:

Conduction velocity in meters per second (m/s).

membrane_currents(downsample=1)

Compute membrane currents at each section for each time point in the simulation.

Uses the methods described in Pena et. al 2024: http://dx.doi.org/10.1371/journal.pcbi.1011833

For myelinated fibers, this calculation includes periaxonal currents between adjacent sections (based on h.Section.xraxial <neuron.Section.xraxial>). The result is a matrix of shape: [num_timepoints, num_sections].

This method returns a tuple consisting of: 1. A 2D array (time steps x number of sections) of membrane currents in mA. 2. The array of time points corresponding to those currents (downsampled by the specified factor).

Parameters:

downsample (int) – Factor to reduce the temporal resolution (e.g., downsample=2 takes every 2nd time step).

Return type:

ndarray

Returns:

A tuple (i_membrane_matrix, downsampled_time). The matrix contains total currents (mA). The time array contains the corresponding simulation times (ms).

Raises:

RuntimeError – If membrane currents or extracellular potentials were not recorded properly.

nodebuilder(ind, node_type)

Create a generic node of Ranvier.

This method sets the node length, diameter, and inserts the extracellular mechanism, but does not add active ion channel mechanisms. Subclasses or external code should add the relevant channels to define the node’s electrophysiological behavior.

Parameters:
  • ind (int) – Index of this node in the overall fiber construction.

  • node_type (str) – A string identifying the node as ‘active’ or ‘passive’.

Return type:

Section

Returns:

The newly created node as a NEURON h.Section.

point_source_potentials(x, y, z, i0, sigma, inplace=False)

Compute extracellular potentials at the fiber’s coordinates due to a point source stimulus.

See also

Documentation on extracellular potentials <extracellular potentials.md> in PyFibers.

Parameters:
  • x (float) – x-coordinate of the source in µm.

  • y (float) – y-coordinate of the source in µm.

  • z (float) – z-coordinate of the source in µm.

  • i0 (float) – Magnitude of the point-source current (mA).

  • sigma (float | tuple) – Conductivity (S/m). A float for isotropic or a tuple (sigma_x, sigma_y, sigma_z) for anisotropic.

  • inplace (bool) – If True, update Fiber.potentials in-place.

Return type:

ndarray

Returns:

Extracellular potentials at each fiber coordinate, in mV.

record_gating(**kwargs)

Record gating parameters (ion channel states) from axon nodes.

The gating variables must be declared in Fiber.gating_variables within the fiber model class.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

dict[str, list[h.Vector | None]]

Returns:

A dictionary mapping gating variable names to lists of recorded NEURON Vector objects.

Raises:

AssertionError – If Fiber.gating_variables is empty.

record_im(**kwargs)

Record membrane current (nA) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane current.

record_sfap(rec_potentials, downsample=1)

Compute the SFAP time course at a given electrode location.

Parameters:
  • rec_potentials (list | ndarray) – 1D array of precomputed potentials (mV) at each fiber section due to the electrode placement.

  • downsample (int) – Downsampling factor for the time vector (applies to the current matrix).

Return type:

ndarray

Returns:

A tuple (sfap_trace, downsampled_time) where sfap_trace is the computed single-fiber action potential in microvolts (µV) and downsampled_time is the corresponding time array (ms).

record_values(ref_attr, allsec=False, indices=None, allow_missing=True, recording_dt=None, recording_tvec=None)

Record NEURON variable references (e.g. membrane voltage) along the fiber.

Note that recording_dt and recording_tvec are mutually exclusive. If both are None, the variable is recorded at every simulation timestep. For more info, see the NEURON docs: https://nrn.readthedocs.io/en/latest/python/programming/math/vector.html#Vector.record

Parameters:
  • ref_attr (str) – The NEURON attribute to record (e.g. '_ref_v').

  • allsec (bool) – If True, record from sections (including nodes). Otherwise, only record from nodes.

  • indices (list[int]) – Specific indices to record from (if None, record from all).

  • allow_missing (bool) – If True, allows missing attributes without raising an error (returns None).

  • recording_dt (float) – The time step [ms] for recording the values (separate from simulation dt). Should be larger than the simulation dt.

  • recording_tvec (h.Vector) –

    A NEURON Vector of time points at which to record the values (ms). Note that the user MUST keep this Vector in memory for the duration of the simulation. This means you must assign it to a variable that is not overwritten or deleted. For example, to record at time points 0, 1, 2, and 3 ms:

    recording_tvec = h.Vector([0, 1, 2, 3])  # store times in a Vector
    fiber.record_values("_ref_v", recording_tvec=recording_tvec)  # pass Vector to record
    stimulation.find_threshold(fiber)  # run the simulation
    plt.plot(recording_tvec, fiber.vm[0])  # plot the recorded values
    

Raises:

ValueError – If indices is an empty list.

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects or None (if allow_missing=True and the requested attribute is missing).

record_vext()

Record extracellular potential (mV) from each section along the fiber.

Return type:

list[Vector]

Returns:

A list of NEURON Vector objects recording extracellular potential.

record_vm(**kwargs)

Record membrane voltage (mV) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane voltage.

resample_potentials(potentials, potential_coords, center=False, inplace=False)

Use linear interpolation to resample external potentials onto the fiber’s coordinate system (1D).

This is used when extracellular potentials are calculated from an external source, such as a finite element model. The potentials provided by the user should be sampled at high resolution along the fiber’s path and provided alongside the corresponding arc-length coordinates.

If center=True, both the input coordinates and the fiber’s coordinates will be shifted such that their midpoints align.

Parameters:
  • potentials (ndarray) – 1D array of external potential values.

  • potential_coords (ndarray) – 1D array of coordinates corresponding to potentials.

  • center (bool) – If True, center the potentials around the midpoint of each domain.

  • inplace (bool) – If True, update Fiber.potentials with the resampled values.

Return type:

ndarray

Returns:

Interpolated potential values aligned with Fiber.longitudinal_coordinates.

Raises:

AssertionError – If input array sizes or monotonicity checks fail.

resample_potentials_3d(potentials, potential_coords, center=False, inplace=False)

Interpolate external potentials onto a 3D fiber coordinate system.

A wrapper around Fiber.resample_potentials() that handles 3D coordinates by first computing the arc length of the provided coordinate array. As with the 1D version, this method is used to resample external potentials (e.g., from a finite element model) onto the fiber.

At present, this does not check that the input coordinates lie exactly along the 3D fiber path. Therefore, it is recommended to use the same coordinates to construct the fiber as to use here. Alternatively, you can create a 1D fiber and calculate the coordinate arc lengths. For more information, see Supplying Extracellular Potentials.

Parameters:
  • potentials (ndarray) – 1D array of external potential values.

  • potential_coords (ndarray) – 2D array of shape (N, 3) representing the (x, y, z) coordinates where the potentials are measured or computed.

  • center (bool) – If True, center the potentials around the midpoint of each domain.

  • inplace (bool) – If True, update Fiber.potentials with the resampled values.

Return type:

ndarray

Returns:

Interpolated potential values aligned with the fiber’s 3D arc-length coordinates.

Raises:

AssertionError – If called on a non-3D fiber or if input coordinate shapes are invalid.

set_save_gating(**kwargs)

Record gating parameters (ion channel states) from axon nodes.

The gating variables must be declared in Fiber.gating_variables within the fiber model class.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

dict[str, list[h.Vector | None]]

Returns:

A dictionary mapping gating variable names to lists of recorded NEURON Vector objects.

Raises:

AssertionError – If Fiber.gating_variables is empty.

set_save_im(**kwargs)

Record membrane current (nA) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane current.

set_save_vext()

Record extracellular potential (mV) from each section along the fiber.

Return type:

list[Vector]

Returns:

A list of NEURON Vector objects recording extracellular potential.

set_save_vm(**kwargs)

Record membrane voltage (mV) along the fiber.

Parameters:

kwargs – Additional arguments passed to Fiber.record_values().

Return type:

list[h.Vector | None]

Returns:

A list of NEURON Vector objects recording membrane voltage.

set_xyz(x=0, y=0, z=0)

Assign new (x, y, z) shifts to a straight (1D) fiber.

The fiber is assumed to be along the z-axis initially, with x=0 and y=0. This method sets the x and y coordinates to the specified values and shifts the z coordinate by the given amount.

Parameters:
  • x (float) – Shift in the x-direction (µm).

  • y (float) – Shift in the y-direction (µm).

  • z (float) – Shift in the z-direction (µm).

Raises:

AssertionError – If this fiber is a 3D fiber (since this method is for 1D only).

Return type:

None

static sfap(current_matrix, potentials)

Compute the Single-Fiber Action Potential (SFAP) by multiplying currents with recording potentials.

Parameters:
  • current_matrix (ndarray) – 2D array of shape [timepoints, sections], containing currents in mA.

  • potentials (ndarray) – 1D array of potentials (mV) at each section, length = number of sections.

Raises:

AssertionError – If the number of columns in the current matrix does not match the length of potentials.

Return type:

ndarray

Returns:

The computed SFAP in microvolts (µV).

pyfibers.fiber.build_fiber(fiber_model, diameter, length=None, n_sections=None, n_nodes=None, enforce_odd_nodecount=True, **kwargs)

Generate a 1D (straight) fiber model in NEURON.

This function creates a model fiber as an instance of the Fiber class using the specific subclass specified from the FiberModel enumerator. with user-specified diameter and length (from one of: number of sections, number of nodes, or length in microns). Additional keyword arguments are forwarded to the fiber model class constructor.

By default, the first section of the fiber is located at the origin (0, 0, 0), and the fiber extends along the z-axis in the positive direction. To change the fiber’s location, the method Fiber.set_xyz() can be used to translate the fiber along the x, y, or z axes. To create fibers along a custom path in 3D space, use build_fiber_3d() instead.

Parameters:
  • fiber_model (FiberModel) – A FiberModel enumerator specifying the type of fiber to instantiate.

  • diameter (float) – The fiber diameter in micrometers (µm).

  • length (float) – The total length of the fiber in micrometers (µm), if defining by length.

  • n_sections (int) – The total number of sections for discretizing the fiber, if defining by sections.

  • n_nodes (int) – The total number of nodes along the fiber, if defining by nodes.

  • enforce_odd_nodecount (bool) – If True, ensure that the number of nodes is odd.

  • kwargs – Additional arguments forwarded to the underlying fiber model class.

Raises:

ValueError – If more than one among length, n_sections, or n_nodes is specified.

Return type:

Fiber

Returns:

A Fiber class instance.

Example Usage

from PyFibers import build_fiber, FiberModel

fiber = build_fiber(fiber_model=FiberModel.MRG_INTERPOLATION, diameter=10, n_nodes=25)
pyfibers.fiber.build_fiber_3d(fiber_model, diameter, path_coordinates, shift=0, shift_ratio=None, center_shift=False, **kwargs)

Generate a 3D fiber model in NEURON based on a specified path.

This function calculates the fiber’s length from the user-supplied path_coordinates and uses it internally to instantiate a 3D fiber model. The coordinates are a 2D numpy array of shape (number_of_points, 3), where each row represents a point in 3D space (x, y, z).

The fiber model will be created by repeating sections along the path until no more nodes can be added without exceeding the path length. By default, the center of the first section is placed at the origin (0, 0, 0), and the fiber extends along the path.

Parameters:
  • fiber_model (FiberModel) – A FiberModel enumerator specifying the type of fiber to instantiate.

  • diameter (float) – The fiber diameter in micrometers (µm).

  • path_coordinates (ndarray) – A numpy array of shape (N, 3) specifying the 3D coordinates (x, y, z) of the fiber path.

  • shift (float) – A shift in microns to apply to the fiber coordinates.

  • shift_ratio (float) – Ratio of the internodal length to shift the fiber coordinates.

  • center_shift (bool) – If True, center the fiber before applying the shift.

  • kwargs – Additional arguments forwarded to the underlying fiber model class.

Raises:

ValueError – If path_coordinates is not provided, or if n_sections, n_nodes, or length is specified (these are invalid in 3D mode).

Return type:

Fiber

Returns:

A fully instantiated 3D fiber model Fiber instance.

Example:

import numpy as np
from PyFibers import build_fiber_3d, FiberModel

coords = np.array(
    [
        [0.0, 0.0, 0.0],
        [0.0, 1.0, 3.0],
        [0.0, 2.0, 7.0],
        # ...
    ]
)
fiber = build_fiber_3d(
    fiber_model=FiberModel.MRG_INTERPOLATION, diameter=10, path_coordinates=coords
)
print(fiber)

pyfibers.model_enum module

Dynamically create an enum of all fiber models, including built-in and plugin models.

This module imports fiber model classes from the local models package, discovers additional plugin models via Python entry points, and aggregates all submodels into a single FiberModel enum. This allows users to refer to fiber models in a uniform way and makes it easier to extend the codebase with new fiber models or external plugins.

Classes:

FiberModel: A dynamically generated enum of fiber models.

enum pyfibers.model_enum.FiberModel(value)

Bases: NoAliasEnum

An enumeration.

Valid values are as follows:

MRG_DISCRETE = <FiberModel.MRG_DISCRETE: <class 'pyfibers.models.mrg.MRGFiber'>>
MRG_INTERPOLATION = <FiberModel.MRG_INTERPOLATION: <class 'pyfibers.models.mrg.MRGFiber'>>
RATTAY = <FiberModel.RATTAY: <class 'pyfibers.models.rattay.RattayFiber'>>
SCHILD94 = <FiberModel.SCHILD94: <class 'pyfibers.models.schild.SchildFiber'>>
SCHILD97 = <FiberModel.SCHILD97: <class 'pyfibers.models.schild.SchildFiber'>>
SMALL_MRG_INTERPOLATION = <FiberModel.SMALL_MRG_INTERPOLATION: <class 'pyfibers.models.mrg.MRGFiber'>>
SUNDT = <FiberModel.SUNDT: <class 'pyfibers.models.sundt.SundtFiber'>>
THIO_AUTONOMIC = <FiberModel.THIO_AUTONOMIC: <class 'pyfibers.models.thio.ThioFiber'>>
THIO_CUTANEOUS = <FiberModel.THIO_CUTANEOUS: <class 'pyfibers.models.thio.ThioFiber'>>
TIGERHOLM = <FiberModel.TIGERHOLM: <class 'pyfibers.models.tigerholm.TigerholmFiber'>>
pyfibers.model_enum.discover_plugins()

Discover plugin classes using entry points under the ‘pyfibers.fiber_plugins’ group.

This function looks for any package registering an entry point with the group ‘pyfibers.fiber_plugins’. Each plugin class must define a submodels attribute (list) that enumerates the specific fiber model(s) identifiers the plugin provides.

Raises:

ValueError – If a discovered plugin class does not contain a submodels attribute.

Return type:

dict[str, type]

Returns:

A dictionary mapping each submodel name (converted to uppercase) to the plugin class.

pyfibers.stimulation module

Defines classes for running simulations using model fibers.

This module provides classes and functionalities to manage stimulation of model fibers. It includes enumerations for different threshold, termination, and bounds search modes, as well as a base Stimulation class and the more specialized ScaledStim class for extracellular stimulation and IntraStim class for intracellular stimulation.

enum pyfibers.stimulation.BisectionMean(value)

Bases: str, Enum

Mean type used during bisection search phase of finding threshold.

Member Type:

str

Valid values are as follows:

GEOMETRIC = <BisectionMean.GEOMETRIC: 'geometric'>
ARITHMETIC = <BisectionMean.ARITHMETIC: 'arithmetic'>
enum pyfibers.stimulation.BoundsSearchMode(value)

Bases: str, Enum

Modes for adjusting bounds in the bounds search phase of finding threshold.

Member Type:

str

Valid values are as follows:

PERCENT_INCREMENT = <BoundsSearchMode.PERCENT_INCREMENT: 'percent'>
ABSOLUTE_INCREMENT = <BoundsSearchMode.ABSOLUTE_INCREMENT: 'absolute'>
class pyfibers.stimulation.IntraStim(dt=0.001, tstop=50, t_init_ss=-200, dt_init_ss=5, istim_ind=None, istim_loc=None, clamp_kws=None)

Bases: Stimulation

Manage intracellular stimulation of model fibers.

This class extends the Stimulation class to provide intracellular stimulation capabilities. The intracellular stimulation is managed via a custom h.trainIClamp mechanism. This mechanism allows for repeated square pulses of current to be injected into a fiber. Its arguments are provided as clamp_kws when creating an instance of this class.

Example Usage

# Assuming you have already created a fiber "my_fiber"
clamp_kws = {
    "delay": 1,  # ms
    "pw": 0.1,  # ms
    "dur": 50,  # ms
    "freq": 100,  # Hz
    "amp": 1,  # nA, recommended to set 1 for scaling purposes
}
istim_ind = 0  # Stimulate the first section
dt, tstop = 0.001, 50
istim = IntraStim(dt=dt, tstop=tstop, istim_ind=istim_ind, clamp_kws=clamp_kws)

# Run a simulation with a given stimulation amplitude
stimamp = 2
istim.run_sim(2, my_fiber)

# Run a threshold search
istim.find_threshold(my_fiber, condition="activation")
_add_istim(fiber)

Create instance of h.trainIClamp for intracellular stimulation.

This method is not called by the user, see :class:IntraStim. Note that trainIClamp is a mod file included in this package. It is an extension of NEURON’s built-in IClamp that allows repeated square pulses.

Parameters:

fiber (Fiber) – The Fiber object to attach intracellular stimulation to.

Return type:

IntraStim

Returns:

The Stimulation instance (self).

run_sim(stimamp, fiber, ap_detect_location=0.9, exit_func=<function IntraStim.<lambda>>, exit_func_interval=100, exit_func_kws=None, use_exit_t=False, fail_on_end_excitation=True, ap_detect_threshold=-30)

Run a simulation for a single stimulation amplitude.

Parameters:
  • stimamp (float) – Amplitude to be applied to extracellular stimulation - Should be a single float for one source - If stimamp is a single float and there are multiple sources, the same stimamp is applied to all sources - If stimamp is a list of floats, each float is applied to the corresponding source

  • fiber (Fiber) – The Fiber to be stimulated.

  • ap_detect_location (float) – Location to detect action potentials (percent along fiber)

  • exit_func (Callable) – Function to call to check if simulation should be exited

  • exit_func_interval (int) – Interval to call exit_func

  • exit_func_kws (dict) – Keyword arguments to pass to exit_func

  • use_exit_t (bool) – If True, use the time returned by exit_func as the simulation end time

  • fail_on_end_excitation (bool) – Behavior for end excitation detection - if True, raise error if end excitation is detected - if False, continue simulation if end excitation is detected - if None, do not check for end excitation

  • ap_detect_threshold (float) – Threshold for detecting action potentials (default: -30 mV)

Raises:

RuntimeError – If NaNs are detected in fiber potentials

Return type:

tuple[int, float]

Returns:

Number of detected APs and time of last detected AP.

class pyfibers.stimulation.ScaledStim(waveform, dt=0.001, tstop=50, t_init_ss=-200, dt_init_ss=5, pad_waveform=True, truncate_waveform=True)

Bases: Stimulation

Manage extracellular stimulation of model fibers.

This class takes one or more waveforms, each of which is expected to match one corresponding set of fiber potentials (i.e., one source) in the fiber being stimulated. Therefore, if you have N potential sets on the fiber, you must provide N waveforms, each describing the time course of stimulation for that source.

The waveform can be either:
  • A callable that accepts a single float argument for the stimulation time (in ms) and returns

the waveform value at that time.
  • A list of such callables with length N if there are N sources

See also

For more information on defining potentials, see Supplying Extracellular Potentials.

** Example Usage **

# Assuming you have already created a fiber "my_fiber"
# Add potentials to the fiber
my_fiber.potentials = electrical_potentials_array


# Create function which takes in a time and returns a waveform value
def my_square_pulse(t: float) -> float:
    if t > 0.1 and t <= 0.2:  # if time is during the pulse
        return 1  # on
    else:  # if time is outside the pulse
        return 0  # off


# Same waveform as above but using scipy.interpolate.interp1d
my_square_pulse = interp1d(
    [0, 0.1, 0.2, 50],  # start, on, off, stop, in ms
    [0, 1, 0, 0],  # waveform values at those times
    kind="previous",
)

# Create a :class:`ScaledStim` object
dt, tstop = 0.001, 50
stim = ScaledStim(waveform=my_square_pulse, dt=dt, tstop=tstop)

# Run the simulation
stim.run_sim(-1, my_fiber)

# Calculate threshold
stim.find_threshold(my_fiber, condition="activation")
run_sim(stimamp, fiber, ap_detect_location=0.9, exit_func=<function ScaledStim.<lambda>>, exit_func_interval=100, exit_func_kws=None, use_exit_t=False, fail_on_end_excitation=True, ap_detect_threshold=-30)

Run a simulation with a given stimulus amplitude(s).

See also

For more information on the underlying math, see the Algorithms docs page.

Parameters:
  • stimamp (float | list[float]) – Amplitude to scale the product of extracellular potentials and waveform. - Should be a single float for one source - If stimamp is a single float and there are multiple sources, the same stimamp is applied to all sources - If stimamp is a list of floats, each float is applied to the corresponding source

  • fiber (Fiber) – The Fiber object to stimulate.

  • ap_detect_location (float) – Normalized location in [0,1] to check for APs.

  • exit_func (Callable[..., bool]) – Callback to check if the simulation can be ended early (e.g., upon detection of an AP).

  • exit_func_interval (int) – How often (in time steps) to call exit_func.

  • exit_func_kws (dict | None) – Additional arguments for exit_func.

  • use_exit_t (bool) – If True, simulation will stop after self._exit_t (if set).

  • fail_on_end_excitation (bool | None) – Behavior for end excitation detection - If True, raise error if end excitation is detected - If False, continue simulation if end excitation is detected - If None, do not check for end excitation

  • ap_detect_threshold (float) – Threshold for detecting action potentials (default: -30 mV)

Raises:

RuntimeError – If NaNs are detected in membrane potentials or if required setup (e.g., istim) is missing.

Return type:

tuple[int, float]

Returns:

Tuple (number_of_APs, time_of_last_AP).

class pyfibers.stimulation.Stimulation(dt=0.001, tstop=50, t_init_ss=-200, dt_init_ss=5, custom_run_sim=None)

Bases: object

Manage stimulation of NEURON simulations.

Provides methods to configure time stepping, run simulations, and perform threshold search. This class is not meant to be used as is; subclasses should override the Stimulation.run_sim() method. For example, the ScaledStim subclass provides extracellular stimulation capabilities.

For more details on using this class to write custom stimulation routines, see Custom Simulation Code.

static ap_checker(fiber, ap_detect_location=0.9, precision=3, check_all_apc=True)

Check how many action potentials occurred at a specified node location.

Parameters:
  • fiber (Fiber) – The Fiber object to evaluate for APs.

  • ap_detect_location (float) – Normalized location along the fiber in [0,1] to check for APs.

  • precision (int) – Decimal places to round the detected AP time.

  • check_all_apc (bool) – If True, raise a warning if APs occur elsewhere but not at the detect location.

Return type:

tuple[int, float | None]

Returns:

A tuple (num_aps, last_ap_time). If no APs are detected, last_ap_time is None.

Raises:

AssertionError – If the detected AP time is non-positive.

static end_excitation_checker(fiber, multi_site_check=True, fail_on_end_excitation=True)

Check for end-excitation.

Determines activation sites by finding local minima in each node’s AP time. If an AP is detected near in or adjacent to the passive end nodes, raise an error or issue a warning based on fail_on_end_excitation.

Parameters:
  • fiber (Fiber) – The Fiber object to check.

  • multi_site_check (bool) – If True, warn if multiple activation sites are detected.

  • fail_on_end_excitation (bool | None) – Controls handling of end-excitation: - True: Raise RuntimeError if end-excitation is detected. - False: Only warn if end-excitation is detected. - None: Skip the check entirely.

Return type:

bool

Returns:

True if end excitation is detected, False otherwise.

Raises:

RuntimeError – If end excitation is detected and fail_on_end_excitation is True.

find_threshold(fiber, condition=ThresholdCondition.ACTIVATION, bounds_search_mode=BoundsSearchMode.PERCENT_INCREMENT, bounds_search_step=10, termination_mode=TerminationMode.PERCENT_DIFFERENCE, termination_tolerance=1, stimamp_top=-1, stimamp_bottom=-0.01, max_iterations=50, exit_t_shift=5, bisection_mean=BisectionMean.ARITHMETIC, block_delay=0, thresh_num_aps=1, silent=False, **kwargs)

Perform a bisection search to find the threshold stimulus amplitude.

This method first locates bounds where one amplitude is subthreshold and

another is suprathreshold (bounds search phase). Then, it repeatedly narrows the bounds until they converge based on the specified termination mode and tolerance (bisection search phase). Ideally, the initial bounds should be chosen such that stimamp_top is supra-threshold and stimamp_bottom is sub-threshold.

Note that enums (ThresholdCondition, BoundsSearchMode,

TerminationMode, BisectionMean) can be provided as an enum member (e.g. ThresholdCondition.ACTIVATION) or as the member’s string value (e.g. “activation”).

See also

For more details on the threshold search process, see Algorithms in PyFibers.

Parameters:
  • fiber (Fiber) – Instance of pyfibers.fiber.Fiber to apply stimulation to.

  • condition (ThresholdCondition) – The threshold condition (ThresholdCondition.ACTIVATION or ThresholdCondition.BLOCK).

  • bounds_search_mode (BoundsSearchMode) – The bounds search mode (BoundsSearchMode.PERCENT_INCREMENT or BoundsSearchMode.ABSOLUTE_INCREMENT).

  • bounds_search_step (float) – The iterative increase/decrease of the upper/lower bound during bounds search - if bounds_search_mode is “percent” this is the percentage increase/decrease - if bounds_search_mode is “absolute” this is the absolute increase/decrease

  • termination_mode (TerminationMode) – The termination mode (TerminationMode.PERCENT_DIFFERENCE or TerminationMode.ABSOLUTE_DIFFERENCE).

  • termination_tolerance (float) – Difference between upper and lower bounds that indicates convergence - absolute difference if termination_mode is “absolute” - percentage difference if termination_mode is “percent”

  • stimamp_top (float) – Initial upper-bound stimulus amplitude to test.

  • stimamp_bottom (float) – Initial lower-bound stimulus amplitude to test.

  • max_iterations (int) – Maximum attempts to find bounding amplitudes before bisection.

  • exit_t_shift (float) – Extra time (ms) after an AP is detected, beyond which the simulation can be cut short.

  • bisection_mean (BisectionMean) – The bisection mean type (BisectionMean.ARITHMETIC or BisectionMean.GEOMETRIC).

  • block_delay (float) – Time (ms) after start to check for a blocked AP, used in block searches.

  • thresh_num_aps (int) – Number of action potentials for threshold search - if activation, suprathreshold requires detected aps >= thresh_num_aps - if block, suprathreshold requires detected aps < thresh_num_aps

  • silent (bool) – If True, suppress print statements for the search process.

  • kwargs – Additional arguments passed to the run_sim method.

Return type:

tuple[float, tuple[int, float | None]]

Returns:

A tuple (threshold_amplitude, (num_detected_aps, last_detected_ap_time)).

Raises:

RuntimeError – If contradictory bounding conditions occur or if the search fails to converge.

pre_run_setup(fiber, ap_detect_threshold=-30)

Prepare the simulation environment before running.

This method sets the temperature, initializes the membrane potential, configures AP detection, and optionally balances the membrane currents for certain fiber models (e.g., Tigerholm). It also applies any intracellular stimulation parameters if provided.

Parameters:
  • fiber (Fiber) – The Fiber object for which the simulation will be configured.

  • ap_detect_threshold (float) – The voltage threshold for detecting action potentials (mV).

Return type:

None

run_sim(*args, **kwargs)

Run a simulation using either a custom_run_sim method or a subclass override.

Return type:

tuple[int, float | None]

supra_exit(fiber, ap_detect_location, thresh_num_aps=1)

Determine if simulation can be exited early, for activation threshold searches only.

Parameters:
  • fiber (Fiber) – The Fiber object to check for an action potential.

  • ap_detect_location (float) – Normalized location to check in [0,1].

  • thresh_num_aps (int) – Number of APs required to consider it suprathreshold.

Return type:

bool

Returns:

True if the specified number of APs has occurred at or before this time.

static threshold_checker(fiber, block=False, ap_detect_location=0.9, block_delay=0, thresh_num_aps=1, check_all_apc=True)

Determine whether a stimulation was above or below threshold, for activation or block.

Parameters:
  • fiber (Fiber) – The Fiber object to evaluate.

  • ap_detect_location (float) – Normalized location in [0,1] where APs are detected.

  • block (bool) – If True, check for block threshold; otherwise, check for activation threshold.

  • block_delay (float) – Time after simulation start to check for block (ms).

  • thresh_num_aps (int) – Number of APs that constitutes a suprathreshold response.

  • check_all_apc (bool) – Passed to Stimulation.ap_checker() for additional warnings.

Return type:

bool

Returns:

True if stimulation is suprathreshold; False if subthreshold.

Raises:
threshsim(stimamp, fiber, condition=ThresholdCondition.ACTIVATION, block_delay=0, thresh_num_aps=1, **kwargs)

Run a single stimulation trial at a given amplitude and check for threshold.

Parameters:
Return type:

tuple[bool, tuple[int, float | None]]

Returns:

A tuple (is_suprathreshold, (num_aps, last_ap_time)).

enum pyfibers.stimulation.TerminationMode(value)

Bases: str, Enum

Modes for determining when to terminate bisection search phase of finding threshold.

Member Type:

str

Valid values are as follows:

PERCENT_DIFFERENCE = <TerminationMode.PERCENT_DIFFERENCE: 'percent'>
ABSOLUTE_DIFFERENCE = <TerminationMode.ABSOLUTE_DIFFERENCE: 'absolute'>
enum pyfibers.stimulation.ThresholdCondition(value)

Bases: str, Enum

Different threshold search conditions.

Member Type:

str

Valid values are as follows:

ACTIVATION = <ThresholdCondition.ACTIVATION: 'activation'>
BLOCK = <ThresholdCondition.BLOCK: 'block'>