PyFibers Documentation

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:

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:

Returns:

The updated Fiber instance after generation.

Raises:

ValueError – 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.

property last_shift_amount: float¶

Amount of coordinate shift from last resampling operation.

Returns:: Shift amount in micrometers (µm). Positive values indicate shift in the positive coordinate direction.

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:

ValueError – 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:

ValueError – 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.
RuntimeError – If no APs are detected at one or both of the measurement nodes.

Return type:

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:

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:: RuntimeError – 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. These potentials should be calculated assuming a 1 mA unit current source from the recording electrode (e.g., using Fiber.point_source_potentials() with i0=1).
downsample (int) – Downsampling factor for the time vector (applies to the current matrix).

Return type:

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, shift=0, shift_ratio=None)¶

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.

Shifting Behavior: For 1D fibers, shifting is handled during potential resampling. This allows you to test various fiber positions relative to a given potential distribution without recreating the fiber. You can apply shifts using either: - shift: a shift in microns, OR - shift_ratio: a fraction of delta_z (the internodal length).

If center=True, both the input coordinates and the fiber’s coordinates will be shifted such that their midpoints align. If a shift is also specified, the shift is applied after centering.

Note: Unlike 3D fibers, 1D fiber shifting is temporary and only affects the potential resampling. The underlying fiber remains unchanged.

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. If a shift is also specified, the shift is applied after centering.
inplace (bool) – If True, update Fiber.potentials with the resampled values.
shift (float) – Shift distance in microns (µm).
shift_ratio (float) – Shift as a ratio of delta_z.

Return type:

Returns:

Interpolated potential values aligned with Fiber.longitudinal_coordinates.

Raises:

ValueError – If input array sizes or monotonicity checks fail, or if potential coordinates don’t span the fiber coordinates.

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.

Note on shifting: Unlike the 1D version, this method does not support shifting parameters (shift, shift_ratio). For 3D fibers, shifting is handled during fiber creation via build_fiber_3d() parameters. If you need to test different fiber positions, create a new 3D fiber with the desired shift parameters.

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:

Returns:

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

Raises:

ValueError – 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:: RuntimeError – 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:

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

Return type:

static sfap(current_matrix, potentials)¶

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

This method uses the principle of reciprocity to calculate the SFAP.

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. These potentials should be for a 1 mA unit current source from the recording electrode.

Raises:

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

Return type:

Returns:

The computed SFAP in microvolts (µV).

property shifted_coordinates: ndarray¶

Coordinates after applying last resampling shift.

These represent the effective position of the fiber sections after the last call to resample_potentials().

Returns:: 1D array of shifted longitudinal coordinates in micrometers (µm).

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'>>¶

PENA = <FiberModel.PENA: <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'>>¶

SWEENEY = <FiberModel.SWEENEY: <class 'pyfibers.models.sweeney.SweeneyFiber'>>¶

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=lambda *x, **y: ..., 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:

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=lambda *x, **y: ..., 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:

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:

RuntimeError – 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:

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, **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
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:

ValueError – If invalid enum values are provided for condition, bounds_search_mode, termination_mode, or bisection_mean.
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:

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:

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:

Returns:

True if stimulation is suprathreshold; False if subthreshold.

Raises:

ValueError – If thresh_num_aps is not positive.
NotImplementedError – If block is True and thresh_num_aps != 1.
RuntimeError – If no APs are detected at all in a block threshold search.

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:

stimamp (float) – Stimulus amplitude to apply.
fiber (Fiber) – The Fiber object to stimulate.
condition (ThresholdCondition) – Threshold condition (ThresholdCondition.ACTIVATION or ThresholdCondition.BLOCK).
block_delay (float) – If condition=BLOCK, time after which AP detection is considered blocked (ms).
thresh_num_aps (int) – Number of APs required to be considered suprathreshold.
kwargs – Additional arguments for the run_sim method.

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.

TerminationMode.PERCENT_DIFFERENCE: Convergence is based on percentage difference between bounds.
TerminationMode.ABSOLUTE_DIFFERENCE: Convergence is based on the absolute difference between bounds.

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.

ThresholdCondition.ACTIVATION:
Search for the minimal stimulus amplitude required to generate an action potential.
ThresholdCondition.BLOCK:
Search for the stimulus amplitude required to block propagation (i.e., prevent action potentials).

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
ValueError – If is_3d is specified in kwargs.
RuntimeError – If node count does not match number of nodes.

Return type:

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=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.

Shifting Behavior: For 3D fibers, shifting is handled during fiber creation (not during potential resampling). This is because 3D fiber geometry represents the physical position of the fiber in space, so shifting affects the actual fiber coordinates. You can apply shifts using either: - shift: a shift in microns, OR - shift_ratio: a fraction of delta_z (the internodal length).

If center=True, the shift is applied after the fiber is first centered about the midpoint of the 3D path. To test different 3D fiber positions, create a new fiber with the desired shift parameters.

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 (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:

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.enable_logging(level=logging.INFO, format_string=None, stream=None)¶

Enable logging output for pyfibers.

Parameters:

level (int) – Logging level. Defaults to logging.INFO. Common values: logging.DEBUG, logging.INFO, logging.WARNING, logging.ERROR
format_string (str) – Custom format string for log messages. Defaults to a standard format showing level, logger name, and message.
stream (TextIO) – Output stream for logging. Defaults to stderr (standard for logging). Use sys.stdout for notebooks or when you want logging mixed with regular output. If None and PYFIBERS_DOCS_BUILD=1, automatically uses stdout for documentation builds.

Return type:

import pyfibers

pyfibers.enable_logging()  # Enable INFO level logging (stderr in normal use, stdout in docs)
pyfibers.enable_logging(stream=sys.stdout)  # Force stdout
pyfibers.enable_logging(level=pyfibers.logging.DEBUG)  # Enable DEBUG level
pyfibers.enable_logging(level=pyfibers.logging.WARNING)  # Only show warnings and errors

pyfibers.register_custom_fiber(fiber_class)¶

Register a custom fiber model class with the FiberModel enum at runtime.

This function allows users to dynamically add custom fiber models to the FiberModel enum without needing to modify the main package or create a plugin.

Parameters:: fiber_class (type) – The fiber model class to register. Must have a ‘submodels’ attribute.
Return type:: None

Subpackages¶

pyfibers.models package

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:: ValueError – 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:

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

Return type:

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:

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:

Returns:

The updated Fiber instance after generation.

Raises:

ValueError – 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.

property last_shift_amount: float¶

Amount of coordinate shift from last resampling operation.

Returns:: Shift amount in micrometers (µm). Positive values indicate shift in the positive coordinate direction.

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:

ValueError – 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:

ValueError – 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.
RuntimeError – If no APs are detected at one or both of the measurement nodes.

Return type:

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:

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:: RuntimeError – 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. These potentials should be calculated assuming a 1 mA unit current source from the recording electrode (e.g., using Fiber.point_source_potentials() with i0=1).
downsample (int) – Downsampling factor for the time vector (applies to the current matrix).

Return type:

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, shift=0, shift_ratio=None)¶

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.

Shifting Behavior: For 1D fibers, shifting is handled during potential resampling. This allows you to test various fiber positions relative to a given potential distribution without recreating the fiber. You can apply shifts using either: - shift: a shift in microns, OR - shift_ratio: a fraction of delta_z (the internodal length).

If center=True, both the input coordinates and the fiber’s coordinates will be shifted such that their midpoints align. If a shift is also specified, the shift is applied after centering.

Note: Unlike 3D fibers, 1D fiber shifting is temporary and only affects the potential resampling. The underlying fiber remains unchanged.

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. If a shift is also specified, the shift is applied after centering.
inplace (bool) – If True, update Fiber.potentials with the resampled values.
shift (float) – Shift distance in microns (µm).
shift_ratio (float) – Shift as a ratio of delta_z.

Return type:

Returns:

Interpolated potential values aligned with Fiber.longitudinal_coordinates.

Raises:

ValueError – If input array sizes or monotonicity checks fail, or if potential coordinates don’t span the fiber coordinates.

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.

Note on shifting: Unlike the 1D version, this method does not support shifting parameters (shift, shift_ratio). For 3D fibers, shifting is handled during fiber creation via build_fiber_3d() parameters. If you need to test different fiber positions, create a new 3D fiber with the desired shift parameters.

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:

Returns:

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

Raises:

ValueError – 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:: RuntimeError – 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:

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

Return type:

static sfap(current_matrix, potentials)¶

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

This method uses the principle of reciprocity to calculate the SFAP.

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. These potentials should be for a 1 mA unit current source from the recording electrode.

Raises:

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

Return type:

Returns:

The computed SFAP in microvolts (µV).

property shifted_coordinates: ndarray¶

Coordinates after applying last resampling shift.

These represent the effective position of the fiber sections after the last call to resample_potentials().

Returns:: 1D array of shifted longitudinal coordinates in micrometers (µm).

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
ValueError – If is_3d is specified in kwargs.
RuntimeError – If node count does not match number of nodes.

Return type:

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=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.

Shifting Behavior: For 3D fibers, shifting is handled during fiber creation (not during potential resampling). This is because 3D fiber geometry represents the physical position of the fiber in space, so shifting affects the actual fiber coordinates. You can apply shifts using either: - shift: a shift in microns, OR - shift_ratio: a fraction of delta_z (the internodal length).

If center=True, the shift is applied after the fiber is first centered about the midpoint of the 3D path. To test different 3D fiber positions, create a new fiber with the desired shift parameters.

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 (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:

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.

pyfibers.model_enum.register_custom_fiber(fiber_class)¶

Register a custom fiber model class with the FiberModel enum at runtime.

This function allows users to dynamically add custom fiber models to the FiberModel enum without needing to modify the main package or create a plugin.

Parameters:: fiber_class (type) – The fiber model class to register. Must have a ‘submodels’ attribute.
Return type:: None

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.

BisectionMean.GEOMETRIC: Use geometric mean (sqrt(bottom * top)).
BisectionMean.ARITHMETIC: Use arithmetic mean ((bottom + top) / 2).

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.

BoundsSearchMode.PERCENT_INCREMENT: Adjust bounds by multiplying/dividing by a percentage factor.
BoundsSearchMode.ABSOLUTE_INCREMENT: Adjust bounds by adding/subtracting a fixed increment.

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=lambda *x, **y: ..., 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:

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=lambda *x, **y: ..., 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:

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:

RuntimeError – 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:

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, **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
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:

ValueError – If invalid enum values are provided for condition, bounds_search_mode, termination_mode, or bisection_mean.
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:

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:

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: