Patch-clamp data analysis in Python: passive membrane properties

The analysis of neurophysiological data is usually time-consuming and programming languages such as R and Python are powerful tools to automatize the analysis of large datasets. Moreover, learning how to code is becoming a valuable and needed skill for neuroscientists, as nicely pointed out by Ashley Juavinett. I will focus on Python because it is an open-source programming language and is widely used for the analysis of neurophysiological data.

In this mini-tutorial, I explain how to analyze some passive membrane and sag properties using the eFEL (Electrophys Feature Extraction Library) Python package. Read the Clampfit tutorial for a brief description of each property. In this blog post series, you can also find tutorials for the analysis of action potentials in Clampfit (post) and Python (post).

I am also a beginner in Python so I will just describe the basics to use the packages, and I will refer to resources to know more about specific topics. Of course, all the credit goes to the creators of the packages and I just hope this brief tutorial can be useful for new users of Python. Please contact me if you have suggestions or corrections to improve this tutorial.


Python installation

This is a brief explanation of how to install Python (there are other ways) and use the packages in Jupyter notebooks.
Important notes:
▪ Replace path/to/directory, filename, or myenv with your actual directory, name of your file, and name of your Python environment, respectively.
▪ Read about the difference between absolute and relative paths.
▪ The # symbol in Python means a comment.

  • Install a Python distribution: Miniconda or Anaconda. The main difference is that Anaconda has more packages pre-installed and has a graphical user interface called Navigator.
  • Run the Anaconda Prompt and type the commands in bold:
    • To change the path to your working directory: cd Path/to/directory.
    • To create a Conda environment (libraries, packages, etc. are isolated from other environments) with some useful packages (Anaconda might have some pre-installed): conda create -n myenv pip jupyterlab matplotlib ipympl numpy pandas scipy
    • Then, activate your new environment: conda activate myenv
    • Other packages can be installed in your Python environment using pip or conda install: pip install name-of-the-package or conda install name-of-the-package
    • Run the Jupyter notebook (more about notebooks here): jupyter lab


Data example

For the scripts below, I used a whole-cell patch-clamp recording (in current-clamp mode) of a pyramidal cell from an acute brain slice (see video) containing the Prefrontal Cortex (post). You can download the ABF file here. Then, paste the scripts below (without the line numbers) into your Jupyter Notebook and play around with them.

pyABF

You need this easy and up-to-date package to read ABF files in Python. The package was created by Scott W Harden and you can read tutorials, documentation, and examples on his website and GitHub. The scripts also include some examples to nicely visualize traces. To learn more about plotting in Python, visit the pyABF and matplotlib websites.
Install pyABF (info): pip install pyabf


eFEL by the Blue Brain Project

EFEL package was created by the Blue Brain Project to extract electrophysiological features recorded from neurons. This package is more versatile than the IPFX from the Allen Institute which, I think, you can only use to analyze action potentials (post) using your own data. The eFEL package also allows you to analyze some passive membrane properties and voltage sag features. You can find documentation, examples, and notebooks on the eFEL GitHub and website.

Install the package (info). It should work in any environment. Activate the environment and install the package in the Anaconda Prompt: pip install efel

Load packages and define paths
This is the common part for the rest of the eFEL features so you only need to write it at the beginning of your notebook. Note: Change path/to/directory and filename for your actual directory and the name of your file, respectively.

# Import the packages 
import efel
import pyabf
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# Define the path to your file
# data = '/Path/to/file/filename.abf'
abf = pyabf.ABF(data)
print (abf)  # Optional: extract information from the ABF file. 

# List of available features in the eFEL package
# efel.getFeatureNames() 


Passive membrane properties and voltage sag characteristics

I used some of the output features of the eFEL package to calculate the input resistance and the capacitance myself because the feature ohmic_input_resistance did not work (let me know if you made it work). I recommend reading the description of the features to understand how they are calculated. You can also read the Clampfit tutorial for a brief description of each property.

# Optional: columns of the output table
table = pd.DataFrame(columns=['voltage','steady_state_voltage_stimend',
                              'voltage_delta',
                              'current_pA',
                              'input_resistance_Gohm',
                              'time_constant',
                              'capacitance_pF',
                              'sag_amplitude', 'sag_ratio1'])  

# Loop function
for sweep in abf.sweepList[0:12]:  # To select a range of traces
    abf.setSweep(sweep)
    # Defines a trace and region of analysis
    trace = {'T': abf.sweepX*1000, # convert to ms
             'V': abf.sweepY,
             'I': abf.sweepC,
             'stim_start' : [230],
             'stim_end' : [700]}
    traces = [trace]
      
    # Output features
    feature_values = efel.getFeatureValues(traces,
                                           ['voltage', 'steady_state_voltage_stimend', 
                                            'current', 'time_constant', 
                                            'sag_amplitude', 'sag_ratio1'],
                                           raise_warnings=None)[0] # If true, returns warnings

    # Optional: create table from the results
    # Use [0] to extract the values from the Spikecount array
    length = len(table)
    table.loc[length, 'voltage'] = feature_values['voltage'][0]
    table.loc[length, 'steady_state_voltage_stimend'] = feature_values['steady_state_voltage_stimend'][0]
    table.loc[length, 'voltage_delta'] = (feature_values['steady_state_voltage_stimend']-feature_values['voltage'])[0]
    table.loc[length, 'current_pA'] = feature_values['current'][4000:6000].mean() # [] to select the rows of current step
    table.loc[length, 'sag_amplitude'] = feature_values['sag_amplitude']
    table.loc[length, 'sag_ratio1'] = feature_values['sag_ratio1']
    
    if (table.loc[length, 'current_pA']  >= -80) & (table.loc[length, 'current_pA']  < 0):
        table.loc[length, 'input_resistance_Gohm'] = (table.loc[length, 'voltage_delta']/table.loc[length, 'current_pA'])
        if (feature_values['time_constant'] is not None):
            table.loc[length, 'time_constant'] = feature_values['time_constant'][0]
            table.loc[length, 'capacitance_pF'] = (feature_values['time_constant']/table.loc[length, 'input_resistance_Gohm'])[0]

# Optional plotting: example with three subplots
fig = plt.figure(figsize=(12, 4))
fig.tight_layout()

# I-V curve
ax1 = fig.add_subplot(1, 2, 1)
currents = [] # Current value between t1 and t2 (ms) for each step
for sweep in abf.sweepList[0:12]:
    abf.setSweep(sweep)
    t1 = int(250*abf.dataPointsPerMs) 
    t2 = int(600*abf.dataPointsPerMs)
    currents.append(np.average(abf.sweepC[t1:t2]))
ax1.scatter(x = currents, y = table.loc[:,'steady_state_voltage_stimend'])
ax1.set_xlabel('Current (pA)')
ax1.set_ylabel('Voltage (mV)')

# Voltage traces
ax2 = fig.add_subplot(2, 2, 2)
for sweep in abf.sweepList[0:12]:
        abf.setSweep(sweep)
        ax2.plot(abf.sweepX*1000, abf.sweepY, alpha=0.3)
ax2.set_xlabel('Time (ms)')
ax2.set_ylabel('Membrane voltage (mV)')
ax2.axes.set_xlim(0, 150)

# Current steps
ax3 = fig.add_subplot(2, 2, 4, sharex=ax2)
for sweep in abf.sweepList[0:12]:
    abf.setSweep(sweep)
    ax3.plot(abf.sweepX*1000, abf.sweepC, alpha=0.3)
ax3.set_xlabel('Time (ms)')
ax3.set_ylabel('Current (pA)')
ax3.axes.set_xlim(0, 1500)

# Optional: Export the graph and the table
# fig.savefig('Path/to/file/filename.png', dpi=300)
# table.to_csv('Path/to/file/filename.csv')
          
# Display the graph and the table
plt.show()
table

[Output]:

You can calculate the input resistance from the slope of the I-V graph using the first hyperpolarizing steps (although it is better to use a recording protocol with smaller hyperpolarizing steps, see Clampfit post).

# I-V graph
fig = plt.figure(figsize=(4, 3))
fig.tight_layout()
ax = fig.add_subplot()
current = [-60, -40, -20, 0] 
voltage = [-77.1, -74.7, -72.1, -68.7]
ax.scatter(current, voltage)
ax.set_xlabel('Current (pA)')
ax.set_ylabel('Voltage (mV)')

# Linear regression using the function polyfit (degree 1)
m, b = np.polyfit(current, voltage, 1)
# Plot regression line
line = np.polyval([m, b], current)
plt.plot(current, line, color='r')

# Print the graph and the input resistance (GOhm)
print (m) # Input resistance (GOhm)
plt.show()

[Output]:

Further reading and resources

Patch-clamp data analysis in Clampfit: passive membrane properties

Patch-clamp is still the gold-standard technique to study the intrinsic electrophysiological properties of individual neurons as well as the synaptic communication between neurons. This series of mini-tutorials explains the analysis of passive membrane properties (the current post), action potentials (post), and synaptic events (upcoming post) using Clampfit (pClamp, Molecular Devices), a user-friendly electrophysiological software used in a lot of laboratories. The basics of the analysis should apply to the analysis of action potentials (post) and passive membrane properties (post) in Python or other software packages, so I hope these brief tutorials can be useful for any newcomer to patch-clamp data analysis. Please contact me if you have suggestions or corrections to improve this tutorial.


Recording protocols in Clampex: passive membrane properties

To keep it simple, I consider that the neurons from brain acute slices (see video) were patched in whole-cell mode using artificial cerebrospinal fluid (ACSF) as the external solution and potassium gluconate as the internal (pipette) solution. Two types of recording protocols can be used depending on the properties to measure: current clamp, the change in membrane potential is measured, and/or voltage clamp, the change in current is measured. I will mainly focus on the current clamp protocols.

Current clamp
You can measure the main three passive membrane properties (input resistance, capacitance, and time constant) and the resting membrane potential using the same current-clamp protocol. You just need a few small current steps around the resting membrane potentials. For example, 6 square current pulses of 500 ms duration from -30 pA to 30 pA. You can see the previous post and video from Molecular Devices on how to set a Clampex protocol.

Voltage clamp
You can use a short and small voltage step (for example, a -5 mV of 100 ms duration) to measure the passive membrane properties. Average the protocol 5-10 times to reduce the noise.


Data analysis in Clampfit: general options

  • Do NOT overwrite the raw data. Save the processed files with another name: File > Save as.
  • Cursors measure voltage and time differences. You can bring the cursors with Shift + Left mouse click.
  • Lock the cursors to each other. Tools > Cursors > Cursor Properties (OR double click on cursor): Lock to partner
  • Select sweeps (traces). View > Select sweeps (or mouse right-click).
  • Visualize traces. View > Data display: Sweeps, Continue, Concatenated.
  • Transfer traces. Edit > Transfer traces: Results window (to export data) or Graph window (display traces).
  • Average traces. Analyze > Average files. A file will be generated with the average traces from the selected files.
  • View > Windows default: Save selected Analysis Window settings as default. This is a handy option when you have to analyze the same type of recordings.
  • Analyze > Power Spectrum. You can leave the default options and check the graph for sources of noise in your recording. Yes, that 50 (Europe) or 60 (US) Hz noise peak (see the fantastic Neurogig posts about electrical noise).
  • Analyze > Filter. See the pClamp user guide on how to filter the noise.
  • Read more about Clampfit options in the pClamp user guide (Help > PDF Manual).


Resting membrane potentials (RMP)

The membrane potential is defined as the potential of the inner side of the membrane relative to the potential of the outer side. The RMP defines the steady-state of the cell when no current is applied. The RMP depends on the intracellular and extracellular concentrations of ions. I include it here but it is not a passive membrane property. Keep in mind that the voltage measured in current and voltage-clamp does not account for the liquid junction potential (LJP) (see Nat Blair’s post to know more). You can correct the LJP offline using the Clampex LJP calculator.

  • Select a region of the pre-pulse membrane potential using cursors 1 and 2.
  • Analyze > Statistics.
    • Trace Selection: Active signal (voltage).
    • Peak polarity: Absolute.
    • Baseline Region: 0 units.
    • Search region. Range: Cursors 1, 2.
    • Measurements: Mean (units).


The three main passive membrane properties are the input/membrane resistance, the time constant, and the input resistance. These three interconnected properties influence synaptic integration, neuronal excitability and determine the speed of signal propagation. For example, inhibitory neurons such as parvalbumin+ interneurons (post), which include Chandelier and basket cells, have fast time constants and low input resistances that contribute to fast synaptic responses.

Input resistance

The input resistance or membrane resistance is the total resistance measured by the amplifier, which the main component is the membrane resistance, i.e., the resistance of the cell membrane to the ions flowing across ion channels. Neurons with high membrane resistance, such as VIP+ interneurons (Tremblay et al., 2016), have high excitability according to Ohm’s law (V = IR) since small current injections will produce high voltage responses. The input resistance is affected by the size of the cell (larger cells, lower input resistance) and the number of open ion channels (more channels, also lower input resistance). In Clampfit you can calculate the membrane resistance using any of the below methods.

Current-clamp

  • Put cursors 1 and 2 on a steady-state region of the voltage response and cursors 3 and 4 on a baseline voltage region (pre-pulse).
  • Analyze > Resistance
    • Experiment type: Current clamp
    • Waveform: current signal (pA)
    • Conversion factor: 0.001 (for pA and mV). Note: change it if needed.
  • Resistance (MOhms) is calculated from Ohm’s law (R = V (mV)/ I (pA)) for each trace.
  • Put cursors 1 and 2 in a steady-state region of the voltage response.
  • Analyze > Quick graph  > I-V. Note: In Current clamp, it is the opposite V-I relationship.
    • X-Axis (current signal). Waveform: Current signal. Waveform at cursor 1.
    • Y-Axis (voltage signal). Region: Cursors 1, 2 (Mean).
  • Analyze > Fit.
    • Predefined function: Select Straight line, standard. Click OK.
    • A regression line is fitted to the I-V data.
  • Results window. The slope of the line (m) is the Input resistance (GΩ) from y = m x + b, where R (m) is calculated in GΩ = mV (b) / pA (y).

Voltage-clamp
Use the first method above to calculate the membrane resistance from Ohm’s law.


Membrane time constant

The change in the membrane potential after the injection of a small current is mainly characterized by a single exponential with a time constant (tau) equals to Rm * Cm, where R is the membrane resistance and C is the membrane capacitance. The time constant is the time to reach 63% of the final voltage. In Clampfit you can calculate the membrane resistance using any of the below methods.

Exponential fitting to the decay voltage response to current steps from -20 pA to +20 pA in 10 pA steps. Recording from a somatostatin interneuron (post) characterized by a slow time constant and high input resistance.

Current-clamp

  • Put cursor 1 at the end of the voltage response and cursor 2 on a steady-state voltage region, between 50 to 300 ms apart from cursor 1 (this depends on the type of neuron).
  • Analyze > Fit
    • Predefined function: Select Exponential, standard. Select Fitting methods automatically. Note: Change to Levenberg-Marquardt if the Chebyshev method gives strange results.
    • Data/Options. Check if the Active Signal and the Traces you want to measure are selected. Optional: Check Compare Models: Automatic. This is useful to check if 2 terms are needed in the model, in case there is a slow inactivation component.
    • Exponential lines are fitted to the voltage traces.
  • Results window. Tau results (ms) and other parameters of the exponential fitting.
  • Put cursor 1 at the end of the voltage response and cursor 2 approx. 100 ms apart from cursor 1, on a steady-state voltage region.
  • Analyze > Statistics
    • Peak polarity: Absolute
    • Baseline region: Cursors 1, 2.
    • Search region: Cursors 3, 4.
    • Measurements: Decay time (ms): From 100% to 37%. Note: depending on the sampling frequency, this method might give unreliable results.
  • Results window.

Voltage-clamp
You can use the same methods as in current-clamp mode.


Membrane capacitance

The lipid membrane of the cell is an insulator and therefore works as a capacitator because it can store charges at each side of the membrane. The amount of charges stored is proportional to the area of the membrane defined by Capacitance =ε0 A/d, where ε0 is the electric constant, d is the distance between the layers (broadly similar across neurons, usually 8-10 nm), and A is the area of the cell. Hence, capacitance is an estimation of neuronal size.

There are many nuances in measuring the capacitance of non-ideal isopotential neurons, both in voltage and current clamp modes (I recommend reading Golowasch et al., 2009).

Current-clamp

  • Calculate the membrane time constant and the input resistance (see above).
  • Calculate the Capacitance (pF) = tau (ms) / Rm (GOhm)


Voltage-clamp
You can calculate Capacitance (pF) = tau (ms) / Rm (GOhm), as in Current-clamp mode.
Alternatively, you can calculate from the total charge (Q), i.e. area of the current transient (Q = pA x ms), accumulated on the membrane during the initial current transient. Then, C (pF) = Q (pA x ms) / V (mV)

  • Put the cursors 1 and 2 on a pre-stimulus region and cursors 3 and 4 between the initial current transient.
  • Analyze > Statistics:
    • Peak polarity: Absolute
    • Baseline region: Cursors 1, 2.
    • Search region: Cursors 3, 4.
    • Measurements: Area.
  • In Results windows, decay time (ms), i.e., tau, can be used to calculate the capacitance: C (nF) = Tau (ms) / R (MOhm) x 1000 = pF


Voltage sag

The initial negative voltage peak in response to a hyperpolarizing step is associated with the hyperpolarization-activated current (Ih), an inward current mediated by HCN channels that start to open at approx. -60 mV (Biel et al., 2008). The Ih currents play an essential role in action potential firing frequency and HCN channels influence differences in membrane properties between humans and mice (Kalmbach et al., 2008). Note that the direct measurement of the Ih current needs to be done in voltage-clamp. Although voltage sag is not a membrane passive property, it also influences the resting membrane potential and sag voltage characteristics are commonly reported as subthreshold properties of the neurons.

Recording protocols in Clampex: voltage sag
Hyperpolarizing current pulses of 500-1000 ms duration from -200 pA to 0 pA in 20 pA steps from resting membrane potentials or a holding potential (for, example -65 mV).
Depolarizing current steps can be added to this protocol to also measure firing properties (see post). For example, from -200 pA to 200 pA (500 ms duration) in 20 pA steps.

Voltage traces in response to hyperpolarizing current steps (from -200 pA to 0 pA) show the characteristic voltage sag of pyramidal neurons. The neuron was recorded from layer 3 of a mouse prefrontal cortex (post).
  • Put cursor 3 at the negative peak of the hyperpolarizing voltage response and cursor 4 on a steady-state point of the same trace. The voltage difference (mV) is the voltage sag. Leave cursors 1 and 2 in a pre-stimulus region.
  • Analyze > Statistics
    • Peak polarity: Absolute
    • Baseline region: Cursors 1, 2.
    • Search region: Cursors 3, 4.
    • Measurements. Click on Peak amplitude (mV) and Antipeak amplitude (mV).
  • Results window. You can subtract the Antipeak amplitude from the Peak amplitude to obtain the voltage sag of each trace. To normalize the values, you can also calculate the sag ratio by dividing the voltage by the negative peak amplitude.


Bibliography

Text edited on May 14, 2022

Patch-clamp data analysis in Clampfit: action potentials

Patch-clamp is still the gold-standard technique to study the passive and active electrophysiological properties of individual neurons as well as the synaptic communication between neurons. This series of mini-tutorials explains the analysis of action potentials (current post), passive membrane properties (post), and synaptic events (upcoming post) using Clampfit (pClamp, Molecular Devices), a user-friendly electrophysiological software used in a lot of laboratories. The basics of the analysis should apply to the analysis of action potentials (post) and passive membrane properties (post) in Python or other software packages, so I hope these brief tutorials can be useful for any newcomer to patch-clamp data analysis. Please contact me if you have suggestions or corrections to improve this tutorial.


Recording protocol in Clampex for action potentials

To keep things simple, I consider that the neurons from brain acute slices (see video) were patched in whole-cell mode using artificial cerebrospinal fluid (ACSF) as the external solution and potassium gluconate as the internal (pipette) solution. Two types of recording protocols can be used depending on the properties to measure: current clamp, the change in membrane potential is measured, and/or voltage clamp, the change in current is measured.

In whole-cell patch clamp (video), you can record action potentials in current-clamp mode. If the neurons do not fire spontaneously, you can either use long square current pulses (500 to 1000 ms steps) to evoke firing patterns or short square current pulses (3 to 5 ms) that are more suitable to study properties (reducing the effect of calcium loading) of single action potentials. Alternatively, a ramp protocol, where you gradually and slowly increase the injected current (25-100 pA/s), can be used to study excitability and action potentials near the action potential threshold. Note that the amplitude and intervals of the current injection depend on each type of neuron. For example, 0 to 200 pA in 20 pA current steps for most of the pyramidal cortical neurons.

An example of how to configure a simple square protocol in Clampex is given below. To know more, see the video from Molecular Devices.

  • Acquire > New Protocol.
    • Mode/Rate.
      • Acquisition Mode: Episodic stimulation.
      • Trial Hierarchy. Trial delay: 0, Runs/trial: 1, Sweeps/run: 16, Sweeps duration (s): 2.5
      • Sampling rate: 20000 Hz
    • Inputs. Channel #0: Vm_prime (voltage signal), Channel #1: Im_sec (current signal). Note: the names of the channels depend on how you have configured MultiClamp. See the guide from Molecular Devices.
    • Outputsmult. Channel #0: I_clamp (current clamp command signal).
    • Waveform. Analog Waveform. Epochs.
      • Epoch A. Type: Step. First level: 0 mV. Delta level: 0 mV. First duration: 200 ms.
      • Epoch B. Type: Step. First level: 0 mV. Delta level: 20 mV. First duration: 500 ms.
      • Epoch C. Type: Step. First level: 0 mV. Delta level: 0 mV. First duration: 1700 ms.
    • Set other parameters as default.
  • Acquire > Save protocol as. Save the protocol.


Data analysis in Clampfit: general options

  • Do NOT overwrite the raw data. Save the processed files with another name: File > Save as.
  • Cursors measure voltage and time differences. You can bring the cursors with Shift + Left mouse click.
  • Lock the cursors to each other. Tools > Cursors > Cursor Properties (OR double click on cursor): Lock to partner
  • Select sweeps (traces). View > Select sweeps (or mouse right-click).
  • Visualize traces. View > Data display: Sweeps, Continue, Concatenated.
  • Transfer traces. Edit > Transfer traces: Results window (to export data) or Graph window (display traces).
  • Average traces. Analyze > Average files. A file will be generated with the average traces from the selected files.
  • View > Windows default: Save selected Analysis Window settings as default. This is a handy option when you have to analyze the same type of recordings.
  • Analyze > Power Spectrum. You can leave the default options and check the graph for sources of noise in your recording. Yes, that 50 (Europe) or 60 (US) Hz noise peak (see the fantastic Neurogig posts about electrical noise).
  • Analyze > Filter. See the pClamp user guide on how to filter the noise.
  • Read more about Clampfit options in the pClamp user guide (Help > PDF Manual).


Action potentials and firing patterns

Neurons encode and transfer information by action potentials (see post). Firing frequency and patterns are also two of the most distinctive electrophysiological features of neurons and these features can be used to classify neuronal types (see post). In addition, the shape of the action potential contains information about the underlying ion channels (see the great review by Bruce Bean). Note that further characterization of channels and currents may require voltage-clamp recordings combined with pharmacological and/or optogenetic (see post) manipulations.


Analysis of action potentials in Clampfit

Analysis window of Clampfit 11 showing the action potentials evoked in a Parvalbumin+ interneuron (post and see also post on Cre-lox mice) from the prefrontal cortex (post). The top panel shows the action potentials detected (blue arrows) between cursors 1 and 2. The window on the right displays the detected events aligned in time. The bottom panel shows the injected current: from 0 to 300 pA in 20 pA steps of 500 ms duration.
  • Select the region with action potentials between cursors 1 and 2.
  • Go to Event detection > Threshold search. Note: Clampfit 11 has optimized the event analysis (see video) but the new features are not freely available yet.
    • Event Polarity: Positive-going.
    • Move horizontal cursor B to select the Baseline and horizontal cursor 0 to select the Trigger (action potential threshold can be between -10 and 0 mV).
    • Optionally, you have the option to add Re-arm to indicate the lower threshold at the end of the event (to remove trunked action potentials, for example). You can also add Rejection to discard a big voltage artifact.
    • Pre-trigger and Post-trigger lengths. Extract a segment of the trace before and after the event.
    • Event rejection. A Min. and max. duration of the action potential can be set to reject artifacts and truncated action potentials.
    • Search Region: Cursors 1, 2. To limit the analysis to the selected region above.
    • SelectActive signal, All visible traces. To limit the search to all traces of the voltage signal.
    • Click OK.
  • Then, click on Event detection > Non-stop (also the flash-forward button in the Toolbar) OR Event detection > Accept the entire category (also the flash-forward with trace button in the Toolbar) to count all the action potentials. All action potentials are identified with a blue target symbol. See the figure above.
  • Results from all voltage traces are displayed in the Event viewer (Event detection > Event viewer). Here you can reject individual events (artifacts) if needed.
  • To see the numerical results, Event detection > Event statistics. For each or all traces, there is a summary of statistics such as peak amplitude, half-width, rise and decay slopes, etc.
  • Go to the Results window and Events tab to see the result for each action potential in all traces. With this and the above results, you get the values or can calculate some of the below features, among many others.
Anatomy of an action potential from the review by Bruce Bean.

Afterhyperpolarization (AHP). You can use the cursors to measure the voltage difference between the threshold and the most negative voltage of the AHP. AHP can also last longer (medium and slow AHP). AHP can slow the rate of firing (increasing the spike frequency adaptation) and, generally, Ca2+-activated K+ currents (SK and BK channels) underlie AHP.

Amplitude or spike height. Event detection > Event statistics or measured with cursors from peak to threshold. The peak amplitude (mV) relative to the resting membrane potential or the afterhyperpolarization (AHP).

Firing rate is the number of action potentials per second. This is one of the key features of neurons and interneurons such as Parvalbumin+ can fire more than 150 spikes per second (see video). You can create Input-output curves by plotting the number of action potentials per current step. One fast way is by copying into Excel the column Trace from the Events tab in the Results window. In Excel, create a COUNTIF formula for each trace number (e.g. COUNTIF($A1:$A$10000,1), etc.) and you get the number of action potentials. Then, just plot those values (y-axis) with the corresponding current value (x-axis) of each trace.

Frequency adaptation can be calculated as the ratio between the time (ms) between the first two action potentials and the last two action potentials in the same trace (choosing a criterion such as a trace +100 pA above the rheobase) using the cursors in Clampfit. Adaptation in the firing pattern is a characteristic feature of neurons such as pyramidal neurons or some interneurons (see the paper by Ascoli and colleagues) and plays a role in filtering input signals. Spike frequency adaptation is caused by several mechanisms, including the M-type potassium current and the AHP currents (see article in Scholarpedia).

Interspike interval (ISI) can be used to calculate the coefficient of variation (CVISI = SDISI /meanISI from the Event statistics) to differentiate between neurons with regular (low CVISI, e.g.: parvalbumin+ interneurons) and irregular or bursting (high CVISI, e.g.: some somatostatin+ interneurons) firing patterns.

Latency is the time between the start of the current stimulus and the time of the first evoked action potential. the cursors. This can be a differential feature between neurons (see Gouwens et al., 2019).

Rheobase is the minimum current that elicits action potential. It gives an idea of how excitable is the cell. Of course, rheobase will depend on the resting potential of the neuron (mV) so, when direct comparisons are needed, it is also calculated from a holding current to keep the neurons at a similar membrane voltage (e.g. -70 mV).

Threshold is the voltage (mV) at which the action potential is triggered. It can be calculated as the voltage point where steep increases suddenly, for example, at dV/dt > 20 mV/s. Go to Edit > Transfer Traces: Results window, Region to transfer: Cursors 1, 2, Active signal. Then, calculate the derivative of the voltage response. In Clampfit, Analyze > Column Arithmetic: cC = diff(cB). Finally, you can create a phase-plane plot by plotting dV/dt (mV/ms) on the y-axis versus the V (mV) on the x-axis (see video, min. 27 onwards).

Upstroke/Downstroke or rise (depolarizing) and decay (repolarizing) phase of the action potential. Event detection > Event statistics. The maximum value of the rise slope (dV/dt) gives information about the maximum sodium current, whereas the maximum decay slope (dV/dt) is slower than the rising phase and it reflects the activation of potassium channels.

Width is the full width (in ms) of the action potential at half height. Event detection > Event statistics. Spike width is mainly dependent on the density and type of potassium channels as well as the inactivation of the sodium channels. For example, Kv3 voltage-gated channels are highly expressed in fast-spiking interneurons (Parvalbumin+) which have short action potentials and short refractory periods, allowing high-frequency firing patterns (Hu et al., 2014).


BIBLIOGRAPHY

— Text updated on May 15, 2022

Patch-clamp data analysis in Python: action potentials

The analysis of neurophysiological data is usually time-consuming and programming languages such as R and Python are powerful tools to automatize the analysis of large datasets. Moreover, learning how to code is becoming a valuable and needed skill for neuroscientists, as nicely pointed out by Ashley Juavinett.

I will focus on Python because it is an open-source programming language and is widely used for the analysis of neurophysiological data. An accessible introduction to Python applied to electrophysiology can be the great online notebooks (web-based interactive files) designed by Ashley Juavinett to explore The Allen Cell Types Database. Check out the course website.

There are also several packages available for the analysis of whole-cell patch-clamp data in Python, which can be a good complement or alternative to data analysis in Clampfit (see posts). However, most of the packages require an intermediate knowledge of Python, whereas others were a standalone effort with no further resources to keep the package up to date. For these reasons, I will focus on two libraries/packages for patch-clamp analysis with current support: IPFX (Intrinsic Physiology Feature eXtractor) from the Allen Institute and eFEL (Electrophys Feature Extraction Library) from the Blue Brain Project. In this blog post series, you can also find tutorials for the analysis of passive membrane properties in Clampfit (post) and Python (post).

I am also a beginner in Python so I will just describe the basics to use the packages, and I will refer to resources to know more about specific topics. Of course, all the credit goes to the creators of the packages and I just hope this brief tutorial can be useful for new users of Python. Please contact me if you have suggestions or corrections to improve this tutorial.


Python installation

This is a brief explanation of how to install Python (there are other ways) and use the packages in Jupyter notebooks.
Important notes:
▪ Replace path/to/directory, filename, or myenv with your actual directory, name of your file, and name of your Python environment, respectively.
▪ Read about the difference between absolute and relative paths.
▪ The # symbol in Python means a comment.

  • Install a Python distribution: Miniconda or Anaconda. The main difference is that Anaconda has more packages pre-installed and has a graphical user interface called Navigator.
  • Run the Anaconda Prompt and type the commands in bold:
    • To change the path to your working directory: cd Path/to/directory.
    • To create a Conda environment (libraries, packages, etc. are isolated from other environments) with some useful packages (Anaconda might have some pre-installed): conda create -n myenv pip jupyterlab matplotlib ipympl numpy pandas scipy
    • Then, activate your new environment: conda activate myenv
    • Other packages can be installed in your Python environment using pip or conda install: pip install name-of-the-package or conda install name-of-the-package
    • Run the Jupyter notebook (more about notebooks here): jupyter lab


Data examples

I will give examples of how to analyze ABF files generated with the software pClamp (Molecular Devices) that can be also analyzed with the software Clampfit (see post). You can also use Clampfit to convert other data files such as ATF to ABF or to NumPy arrays (CSV, comma delimited).

For the scripts below, I used a whole-cell patch-clamp recording (in current-clamp mode) of a Parvalbumin interneuron (post) from an acute brain slice (see video) containing the Prefrontal Cortex (post). You can download the ABF file here, and the same file in csv. format here. Then, paste the scripts below into your Jupyter Notebook and play around with them.


pyABF

You need this easy and up-to-date package to read ABF files in Python. The package was created by Scott W Harden and you can read tutorials, documentation, and examples on his website and GitHub. The scripts also include some examples to nicely visualize traces. To learn more about plotting in Python, visit the pyABF and matplotlib websites.
Install pyABF (info): pip install pyabf


IPFX by The Allen Institute for Brain Science

This package was created to analyze intrinsic cell features from electrophysiological data collected by The Allen Institute (NWB files). As far as I know, only a couple of features work also on NumPy arrays (but it would be great if the Allen Institute makes all the features available in the future). You can find documentation and tutorials on the IPFX GitHub and website,

Install the package (info). A new environment is required and it seems that IPFX does not work with the latest versions of Python yet.

  • Create a new environment with packages (see above): conda create -n ipfx
  • Activate the environment and type:  conda install python=3.7 anaconda
  • Then, install the IPFX package via pip: pip install ipfx Note: The installation may take several minutes and errors may show during the installation.


Load packages and define paths for ABF files
This is the common part for the rest of the IPFX features so you just need to write it at the beginning of your notebook. Note: Replace path/to/file, and filename, with your actual directory and filename.

# Import the packages 
from ipfx.feature_extractor import (SpikeFeatureExtractor,
                                    SpikeTrainFeatureExtractor)
import pyabf
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# Set the file path
data = 'path/to/file/filename.abf'
abf = pyabf.ABF(data)

print (abf)  # Optional: extract information from the ABF file  


Spike features (sfx) from one or several traces
You get more than 30 features for each action potential with this function but some of them are not very useful. See the IPFX documentation for a definition of each parameter.
You can play around with start, end to adapt the range of the analysis to the first action potentials in the trace (like the example below), to a ramp or short current protocol, etc.
The loop function for i in x also works if you only have one trace. If you do not need the loop, remove for sweepNumber in abf.sweepList, and type instead abf.setSweep(sweepNumber=0, channel=0) (or the number and channel you need), and then move all the lines below to the left.

# Set the table from the dataframe below
dataframe = [] 

# Loop function to analyze each voltage trace of the file
for sweepNumber in abf.sweepList: # e.g. abf.sweepList[1:3] to select the range
    abf.setSweep(sweepNumber)
    time = abf.sweepX
    voltage = abf.sweepY
    current = abf.sweepC
    # Define the region (in seconds)
    start, end = 0.24, 0.265
    # Parameters for analysis: 
    # Filter = None to avoid conflicts with your recording settings
    # dv/dt cutoff (mV/ms), min_peak (mV)
    sfx = SpikeFeatureExtractor (start, end, 
                                 filter=None, 
                                 dv_cutoff=20.0, min_peak=0) 
    sfx_results = sfx.process(time, voltage, current)
    dataframe.append(sfx_results)  # To get the mean: df.append(sfx_results.mean())
    
# Table with the features from all the action potentials
table = pd.concat(dataframe)

# Optional: Plot the trace/s
plt.figure(figsize=(6,4))
plt.xlabel ("Time (s)")
plt.ylabel("Voltage (mV)")
for sweepNumber in abf.sweepList:  # Loop to plot all the traces
    abf.setSweep(sweepNumber)
    plt.plot(abf.sweepX, abf.sweepY, alpha=.6, label="sweep %d" % (sweepNumber))
# To highllight one trace
abf.setSweep(3) 
plt.plot(abf.sweepX, abf.sweepY, linewidth=1, color='black')
plt.xlim(0.24, 0.265)

# Export the graph (png) and the table (csv)
plt.savefig('path/to/file/filename.png', dpi=300)
table.to_csv('path/to/file/spike_features.csv')
    
# Display the graph and the table
plt.show()
table

# Remove the below # to show only selected columns. E.g:
#columns = ['threshold_i', 'threshold_v', 'width', 'upstroke', 'downstroke']
#table[columns]

Output:


Spike train features (stfx) from one trace
You can obtain some useful features related to the firing rate (see post): adaptation, latency, CV of ISI, mean of ISI, and first ISI. You can also change start, end to adapt the analysis to your recording.

# Select the sweep/trace and channel
abf.setSweep(sweepNumber=10, channel=0)

# Define each variable
time = abf.sweepX
voltage = abf.sweepY
current = abf.sweepC

# Optional: Current step
currents = [] # Current value between t1 and t2 (ms) for each step
t1 = int(400*abf.dataPointsPerMs) 
t2 = int(500*abf.dataPointsPerMs)
current_mean = np.average(abf.sweepC[t1:t2])

# Define the region (in sec) and detection parameters
start, end = 0, 1
sfx = SpikeFeatureExtractor (start=start, end=end, 
                             filter=None, 
                             dv_cutoff=20.0, min_peak=0) 
sfx_results = sfx.process(time, voltage, current)
stfx = SpikeTrainFeatureExtractor (start, end)
stfx_results = stfx.process(time, voltage, current, sfx_results)

# Optional: Create a table with the stfx results
table = pd.DataFrame()
length = len(table)
table.loc[length, 'Current_step'] = current_mean
table.loc[length, 'adaptation'] = stfx_results ["adapt"]
table.loc[length, 'Latency'] = stfx_results ["latency"]
table.loc[length, 'ISI_CV'] = stfx_results ["isi_cv"]
table.loc[length, 'ISI_mean'] = stfx_results ["mean_isi"]
table.loc[length, 'ISI_first'] = stfx_results ["first_isi"]
table.loc[length, 'Firing_rate_Hz'] = stfx_results ["avg_rate"]

# Optional: Plotting the voltage trace and the current step
fig = plt.figure(figsize=(7, 5)) #Figure size
# Voltage traces
ax1 = fig.add_subplot(211)
ax1.plot(abf.sweepX, abf.sweepY)
# Show detection of action potential peak
ax1.plot(sfx_results["peak_t"], sfx_results["peak_v"], 'r.')
# Show detection of action potential threshold
ax1.plot(sfx_results["threshold_t"], sfx_results["threshold_v"], 'k.')
ax1.set_ylabel("Voltage (mV)")
# Current stimulus
ax2 = fig.add_subplot(212, sharex=ax1) 
ax2.plot(abf.sweepX, abf.sweepC, color='r')
ax2.set_ylabel("Current (pA)")
ax2.set_xlabel("Time (s)")
# Zoom in and out
ax2.axes.set_xlim(0, 1) # Shared x-axis (range in sec)

# Optional: Export the table (csv file) and the figure
fig.savefig('Path/to/file/filename.png', dpi=300)
table.to_csv('Path/to/file/filename.csv')

# View the graph and the table
plt.show()
table

Output:


Spike train features (stfx) from several traces
To analyze several traces, you can use the ForLoop function in Python. It is also possible to plot the input-output curve (no. of action potentials vs current stimulus).

# Set the dataframe for the table
table = pd.DataFrame()

# Loop function to analyze each voltage and current trace of the file
for sweep in abf.sweepList:
    abf.setSweep(sweep)
    time = abf.sweepX
    voltage = abf.sweepY
    current = abf.sweepC
    
    currents = [] # Current value between t1 and t2 (ms) for each step
    t1 = int(400*abf.dataPointsPerMs) 
    t2 = int(500*abf.dataPointsPerMs)
    current_mean = np.average(abf.sweepC[t1:t2])
    
    start, end = 0.1, 1.1
    sfx = SpikeFeatureExtractor(start=start, end=end, filter=None)
    sfx_results = sfx.process(time, voltage, current)
    stfx = SpikeTrainFeatureExtractor (start=start, end=end)
    stfx_results = stfx.process(time, voltage, current, sfx_results)
    
    # Create a table with the stfx results
    length = len(table)
    table.loc[length, 'current_step'] = current_mean
    table.loc[length, 'avg_rate'] = stfx_results["avg_rate"]
    if stfx_results ['avg_rate'] > 0:
        table.loc[length, 'adaptation'] = stfx_results ["adapt"]
        table.loc[length, 'Latency_s'] = stfx_results ["latency"]
        table.loc[length, 'ISI_CV'] = stfx_results ["isi_cv"]
        table.loc[length, 'ISI_mean_ms'] = stfx_results ["mean_isi"]*1000
        table.loc[length, 'ISI_first_ms'] = stfx_results ["first_isi"]*1000
     
# Optional Plotting: example with three subplots
fig = plt.figure(figsize=(15, 5))

# Input-output curve
ax1 = fig.add_subplot(1, 2, 1)
currents = []
for sweep in abf.sweepList:
    abf.setSweep(sweep)
    currents.append(np.average(abf.sweepC[t1:t2]))
ax1.plot(currents, table.loc[:,'avg_rate'])
ax1.set_xlabel('Current step (pA)')
ax1.set_ylabel('Action potentials')

# All traces
ax2 = fig.add_subplot(2, 2, 2)
for sweep in abf.sweepList:
    abf.setSweep(sweep)
    ax2.plot(abf.sweepX, abf.sweepY, alpha=.6, label="Sweep %d" % (sweepNumber))
ax2.set_ylabel('Membrane voltage (mV)')
# ax2.legend() # Optional
# To highllight one trace
abf.setSweep(3) 
ax2.plot(abf.sweepX, abf.sweepY, linewidth=1, color='black')
ax2.axes.set_xlim(0, 1.1) # Range of the x-axis (seconds)

# All current steps
ax3 = fig.add_subplot(2, 2, 4, sharex=ax2) 
for sweep in abf.sweepList:
    abf.setSweep(sweep)
    current = abf.sweepC
    ax3.plot(abf.sweepX, current)
ax3.set_ylabel("Current (pA)")
ax3.set_xlabel("Time (s)")
# To highllight one current trace
abf.setSweep(3) 
ax3.plot(abf.sweepX, abf.sweepC, linewidth=1, color='black')
fig.tight_layout()

# Optional: Export the graph (png) and the table (csv)
fig.savefig("Path/to/file/filename.png", dpi=300)
table.to_csv('Path/to/file/filename.csv')

# Show the graph and the table results 
plt.show()
table

Output:


Spike Train Features from data in NumPy arrays
This post does not aim to be the “definitive guide” but you can also use NumPy arrays (see tutorial) as input data. In Clampfit: Edit > Transfer traces: Results window (selecting the traces you want to export) and save the results as “.csv” file. I attach one example for the analysis of spike train features.

# Load the file
data = np.loadtxt(fname="Path/to/file/filename.csv", delimiter=",")

table = pd.DataFrame()

# Define the time, voltage, and current columns
time = data [:, 0]/1000
voltages = data [:, 1:14]
currents = data [:, 14:28]

# Create a loop
for voltage in voltages.T:
    start, end = 0.1, 1.1
    stfx = SpikeTrainFeatureExtractor (start=start, end=end)
    sfx = SpikeFeatureExtractor(start=start, end=end, filter=None)
    spikes_df = sfx.process(t=time, v=voltage, i=current)
    stfx = stfx.process(t=time, v=voltage, i=current, spikes_df=spikes_df)

    # Create the table    
    length = len(table)
    table.loc[length, 'avg_rate'] = stfx["avg_rate"]
    if stfx['avg_rate'] > 0:
        table.loc[length, 'adaptation'] = stfx ["adapt"]
        table.loc[length, 'Latency'] = stfx ["latency"]
        table.loc[length, 'ISI_CV'] = stfx ["isi_cv"]
        table.loc[length, 'ISI_mean'] = stfx ["mean_isi"]
        table.loc[length, 'ISI_first'] = stfx ["first_isi"] 

# Optional: plotting
fig = plt.figure(figsize=(6, 4))
plt.plot (time, voltages, linewidth=1.0, alpha=.6)
trace_highlight = data[:, 4]  # Highlight a trace of interest in the plot
plt.plot (time, trace_highlight,  linewidth=1, color='black')
plt.xlabel ("Time (s)")
plt.ylabel("Voltage (mV)")
plt.xlim(0, 1)

# Export the graph (png) and the table (csv)
fig.savefig("Path/to/file/filename.png", dpi=300)
table.to_csv('Path/to/file/filename.csv')

# Display the graph and the table
plt.show()
table

# Optional: Get the current step values
# t1 = currents[(time > 0.3) & (time < 0.4), :]
# t2 = currents[(time > 0) & (time < 0.1), :]
# current_mean = np.mean(t1, axis=0) - np.mean(t2, axis=0)
# list(current_mean)

Output:


eFEL by the Blue Brain Project

EFEL package was created by the Blue Brain Project to extract electrophysiological features recorded from neurons. This package is more versatile than IPFX and allows you to analyze passive membrane properties (see post), different file formats, and a ton of other features (check the description of features here). You can find documentation, examples, and notebooks on the eFEL GitHub and website.

Install the package (info). It should work in any environment (you can install it in the same environment as the IPFX, see above). Activate the environment and install the package in the Anaconda Prompt: pip install efel

Load packages and define paths
This is the common part for the rest of the eFEL features so you just need to write it at the beginning of your notebook. Note: Change path/to/directory and filename for your actual directory and name of your file, respectively.

# Import the packages 
import efel
import pyabf
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import nanmean

# Define the path to your file
# data = 'path/to/file/filename.abf'
abf = pyabf.ABF(data)
# print (abf)  # Optional: extract information from the ABF file. 

# List of available features in the eFEL package
# efel.getFeatureNames() 


Spike features from one sweep
You can use the code below to analyze only one action potential or the first action potential by limiting the stim_start and stim_end. Edit the output features (see above) in efel.getFeatureValues and add the new features to the table (in table and a new table.loc line).

# Set the sweep and channel
sweep = abf.setSweep(sweepNumber=10, channel=0)

# Define the variables
time = abf.sweepX*1000 # in miliseconds
voltage = abf.sweepY
current = abf.sweepC

# Define the variables and region of analysis
trace = {'T': abf.sweepX*1000, 
         'V': abf.sweepY,
         'stim_start': [100],
         'stim_end': [800]} 
traces = [trace]

# Optional: Current step values
currents = [] # Current value between t1 and t2 (ms) for each step
t1 = int(400*abf.dataPointsPerMs) 
t2 = int(500*abf.dataPointsPerMs)
current_mean = np.average(abf.sweepC[t1:t2])

# Detection parameters
efel.api.setThreshold(0)  # Voltage threshold 
efel.api.setDerivativeThreshold(20) # dV/dt threshold 

# Define the output features
# Or use instead efel.getMeanFeatureValues() to get mean values
feature_values = efel.getFeatureValues(traces,
                                       ['AP_amplitude', 'AP_width',
                                        'AP_begin_voltage', 'AP_rise_rate',
                                        'AP_fall_rate'],
                                       raise_warnings=None)[0] # If true, returns warnings (e.g. no spikes in trace)
# Results without table and graph
# traces_results

# OptionaL Create a table with the results
table = pd.DataFrame(columns=['Sweep', 'Current_step', 'AP_amplitude', 'AP_width', 
                              'AP_begin_voltage', 'AP_rise_rate', 
                              'AP_fall_rate'])

# Optional: Create a table with the results
# [0] returns values of first action potential, [1], for 2nd, etc. 
length = len(table)
table.loc[length, 'Sweep'] = 7 # SweepNumber
table.loc[length, 'Current_step'] = current_mean
table.loc[length, 'AP_amplitude'] = list(feature_values['AP_amplitude'])[0]
table.loc[length, 'AP_width'] = feature_values['AP_width'][0]
table.loc[length, 'AP_begin_voltage'] = feature_values['AP_begin_voltage'][0]
table.loc[length, 'AP_rise_rate'] = feature_values['AP_rise_rate'][0]
table.loc[length, 'AP_fall_rate'] = feature_values['AP_fall_rate'][0]

# Optional: ploting the trace
fig = plt.figure(figsize=(10, 4))
ax1 = fig.add_subplot(121)
ax1.plot(time, voltage)
ax1.set_xlabel('Time (ms)')
ax1.set_ylabel('Membrane voltage (mV)')
ax1.axes.set_xlim(0, 1000)
# Zoom in to one action potential
ax2 = fig.add_subplot(122, sharey=ax1)
ax2.plot(time, voltage)
ax2.set_xlabel('Time (ms)')
ax2.set_ylabel('Membrane voltage (mV)')
ax2.axes.set_xlim(240, 254)
plt.tight_layout()
plt.show()

# Optional: export the graph and the table
fig.savefig('Path/to/file/filename.png', dpi=300)
table.to_csv('Path/to/file/filename.csv')

# Display the table and the graph
plt.show()
table

Output:


Firing properties from all the traces
You can adapt the previous block for several traces by using the ForLoop function in Python. The eFEL package has also a cool set of features to analyze bursts (see above efel.getFeatureNames()) that may be even useful for cell-attached recordings (Perkins, 2006). You can change the burst parameter by adding efel.api.setFeatureDouble('burst_factor', [1.5]).

# Define the result outputs for the table
table = pd.DataFrame(columns=['Current_step', 'Spikecount',
                              'adaptation', 'Latency_ms', 
                              'ISI_CV', 'ISI_mean_ms'])  

# Loop function
# Define region of analysis
for sweep in abf.sweepList: # e.g. abf.sweepList[1:3] to select the range
    abf.setSweep(sweep)
    # Defines a trace
    trace = {'T': abf.sweepX*1000, 
             'V': abf.sweepY,
             'stim_start': [100],
             'stim_end': [800]} 
    traces = [trace]
    
    # Define the parameters for detection
    efel.api.setThreshold(0) # Voltage threshold for detection
    efel.api.setDerivativeThreshold(20) # dV/dt threshold for detection
    
    # Define the output results
    feature_values = efel.getFeatureValues(traces,
                                           ['Spikecount','adaptation_index',
                                            'time_to_first_spike', 
                                            'ISI_CV', 'ISI_values'],
                                           raise_warnings=None)[0] # If true, returns warnings

    # Optional: add a column with the current steps
    current = abf.sweepC
    currents = [] # Current value between t1 and t2 (ms) for each step
    t1 = int(400*abf.dataPointsPerMs) 
    t2 = int(500*abf.dataPointsPerMs)
    current_mean = np.average(abf.sweepC[t1:t2])
    currents.append(current_mean)
    

    # Create table from the results
    # Use [0] to extract the values from the Spikecount array
    length = len(table)
    table.loc[length, 'Current_step'] = current_mean
    table.loc[length, 'Spikecount'] = feature_values['Spikecount'][0] 
    # Some features requires AP > 0 or more
    if feature_values['Spikecount'] is not None: 
        table.loc[length, 'Latency_ms'] = feature_values['time_to_first_spike']
        if feature_values['Spikecount'] > 4:  
            table.loc[length, 'adaptation'] = feature_values['adaptation_index'][0]
            table.loc[length, 'ISI_CV'] = feature_values['ISI_CV'][0] 
            table.loc[length, 'ISI_mean'] = feature_values['ISI_values'][0]/1000

# Optional plotting: example with three subplots
# Optional: ploting the trace
fig = plt.figure(figsize=(12, 4))
ax1 = fig.add_subplot(121)
for sweep in abf.sweepList:
    abf.setSweep(sweep)
    ax1.plot(abf.sweepX*1000, abf.sweepY, alpha=0.3)
ax1.set_xlabel('Time (ms)')
ax1.set_ylabel('Membrane voltage (mV)')
ax1.axes.set_xlim(0, 1000)
# Plot an individual trace
ax2 = fig.add_subplot(122, sharey=ax1)
sweep_label = abf.setSweep(10) 
ax2.plot(time, voltage, label='sweep 10')
ax2.set_xlabel('Time (ms)')
ax2.set_ylabel('Membrane voltage (mV)')
ax2.axes.set_xlim(0, 1000)
plt.tight_layout()

# Optional : Export the graph and the table
fig.savefig('Path/to/file/filename.png', dpi=300)
table.to_csv('Path/to/file/filename.csv')
            
# Display the graph and the table
plt.show()
table

Output:


Spike Train Features from data in NumPy arrays
The package pyABF nicely extracts the data from the traces but you can directly use NumPy arrays (see tutorial) as input data. In Clampfit: Edit > Transfer traces: Results window and save the results as “.csv” file. I attach one example for the analysis of spike train features.

data = np.loadtxt(fname="path/to/file/filename.csv", delimiter=",")

# Define the columns with the time, current, and voltage values
time = data[:, 0] # in ms
voltages = data [:, 1:14]
currents = data [:, 14:28]

# Create a table with the values
table = pd.DataFrame(columns=['Spikecount',
                              'adaptation', 'Latency_ms', 
                              'ISI_CV', 'ISI_mean_ms']) 

# A loop for each voltage column
for voltage in voltages.T:
    trace = {'T' : time,
             'V' : voltage,
             'stim_start' : [100],
             'stim_end' : [700]}
    traces = [trace]
    
    # Define the parameters for detection
    efel.api.setThreshold(0) # Voltage threshold for detection
    efel.api.setDerivativeThreshold(20) # dV/dt threshold for detection

    # Select which features you want to calculate. See efel.getFeatureNames()
    feature_values = efel.getFeatureValues(traces, ['Spikecount','adaptation_index',
                                                    'time_to_first_spike','ISI_CV',
                                                    'ISI_values'],
                                                    raise_warnings=None)[0]
             
    # Create a table with the results
    length = len(table)
    table.loc[length, 'Spikecount'] = feature_values['Spikecount'][0] 
    # Some features requires AP > 0 or more
    if feature_values['Spikecount'] is not None: 
        table.loc[length, 'Latency_ms'] = feature_values['time_to_first_spike']
        if feature_values['Spikecount'] > 4:  
            table.loc[length, 'adaptation'] = feature_values['adaptation_index'][0]
            table.loc[length, 'ISI_CV'] = feature_values['ISI_CV'][0] 
            table.loc[length, 'ISI_mean_ms'] = feature_values['ISI_values'][0]/1000

# Optional: Plotting
fig = plt.figure(figsize=(6, 4))
plt.plot (time, voltages, linewidth=1.0, alpha=.6)
trace_highlight = data[:, 8] # Highlight a trace of interest in the plot
plt.plot (time, trace_highlight,  linewidth = 1, color='black')
plt.xlabel ("Time (ms)")
plt.ylabel("Voltage (mV)")
plt.xlim(0, 1000)

# Optional : Export the graph and the table
plt.savefig('Path/to/file/filename.png', dpi=300)
table.to_csv('Path/to/file/filename.csv')

# Display the graph and the table
plt.show()
table

#Optional: Get the current step values
# t1 = currents[(time > 300) & (time < 400), :]
# t2 = currents[(time > 0) & (time < 100), :]
# current_mean = np.mean(t1, axis=0) - np.mean(t2, axis=0)
# list(current_mean)

Output:


Further reading and resources

— Updated on May 7, 2022