Module LFPy
¶
Initialization of LFPy, a module for simulating extracellular potentials.
Group of Computational Neuroscience, Department of Mathematical Sciences and Technology, Norwegian University of Life Sciences.
Copyright (C) 2012 Computational Neuroscience Group, NMBU.
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
Classes: |
|
---|---|
Modules: |
|
class Cell
¶
-
class
LFPy.
Cell
(morphology, v_init=-65.0, passive=True, Ra=150, rm=30000, cm=1.0, e_pas=-65.0, extracellular=True, timeres_NEURON=0.125, timeres_python=0.125, tstartms=0, tstopms=100, nsegs_method='lambda100', lambda_f=100, d_lambda=0.1, max_nsegs_length=None, delete_sections=True, custom_code=None, custom_fun=None, custom_fun_args=None, pt3d=False, celsius=None, verbose=False)[source]¶ Bases:
object
The main cell class used in LFPy.
Arguments:
morphology : [str]: path/to/morphology/file v_init: [-65.]: initial potential passive: [True]/False: passive mechs are initialized if True Ra: [150.]: axial resistance rm: [30000]: membrane resistivity cm: [1.0]: membrane capacitance e_pas: [-65.]: passive mechanism reversal potential extracellular: [True]/False: switch for NEURON's extracellular mechanism timeres_NEURON: [0.1]: internal dt for NEURON simulation timeres_python: [0.1]: overall dt for python simulation tstartms: [0.]: initialization time for simulation <= 0 ms tstopms: [100.]: stop time for simulation > 0 ms nsegs_method: ['lambda100']/'lambda_f'/'fixed_length': nseg rule max_nsegs_length: [None]: max segment length for method 'fixed_length' lambda_f: [100]: AC frequency for method 'lambda_f' d_lambda: [0.1]: parameter for d_lambda rule delete_sections: [True]: delete pre-existing section-references custom_code: [None]: list of model-specific code files ([.py/.hoc]) custom_fun: [None]: list of model-specific functions with args custom_fun_args: [None]: list of args passed to custom_fun functions pt3d: True/[False]: use pt3d-info of the cell geometries switch celsius: [None]: Temperature in celsius. If nothing is specified here or in custom code it is 6.3 C verbose: True/[False]: verbose output switch
Usage of cell class:
import LFPy cellParameters = { 'morphology' : 'path/to/morphology', 'rm' : 30000, 'cm' : 1.0, 'Ra' : 150, 'timeres_NEURON' : 0.1, 'timeres_python' : 0.1, 'tstartms' : -50, 'tstopms' : 50, } cell = LFPy.Cell(**cellParameters) cell.simulate()
-
cellpickler
(filename)[source]¶ Save data in cell to filename, using cPickle. It will however destroy any neuron.h objects upon saving, as c-objects cannot be pickled
Usage:
cell.cellpickler('cell.cpickle')
To load this cell again in another session:
import cPickle f = file('cell.cpickle', 'rb') cell = cPickle.load(f) f.close()
alternatively:
import LFPy cell = LFPy.tools.load('cell.cpickle')
-
chiral_morphology
(axis='x')[source]¶ Mirror the morphology around given axis, (default x-axis), useful to introduce more heterogeneouties in morphology shapes
kwargs:
axis : str 'x' or 'y' or 'z'
-
get_closest_idx
(x=0, y=0, z=0, section='allsec')[source]¶ Get the index number of a segment in specified section which midpoint is closest to the coordinates defined by the user kwargs:
x: float, coordinate y: float, coordinate z: float, coordinate section: str, string matching a section-name
-
get_idx
(section='allsec', z_min=-10000, z_max=10000)[source]¶ Returns neuron idx of segments from sections with names that match the pattern defined in input section on interval [z_min, z_max].
kwargs:
section: str, any entry in cell.allsecnames or just 'allsec'. z_min: float, depth filter z_max: float depth filter
Usage:
idx = cell.get_idx(section='allsec') print idx idx = cell.get_idx(section=['soma', 'dend', 'apic']) print idx
-
get_idx_children
(parent='soma[0]')[source]¶ Get the idx of parent’s children sections, i.e. compartments ids of sections connected to parent-argument
kwargs:
parent: str name-pattern matching a sectionname
-
get_idx_name
(idx=array([0]))[source]¶ Return NEURON convention name of segments with index idx. The returned argument is a list of tuples with corresponding segment idx, section name, and position along the section, like; [(0, ‘neuron.h.soma[0]’, 0.5),]
kwargs:
idx : np.ndarray, dtype int segment indices, must be between 0 and cell.totnsegs
-
get_idx_parent_children
(parent='soma[0]')[source]¶ Get all idx of segments of parent and children sections, i.e. segment idx of sections connected to parent-argument, and also of the parent segments
kwargs:
parent: str name-pattern matching a sectionname
-
get_idx_polygons
(projection=('x', 'z'))[source]¶ for each segment idx in celll create a polygon in the plane determined by the projection kwarg (default (‘x’, ‘z’)), that can be visualized using plt.fill() or mpl.collections.PolyCollection
Returned argument is a list of (np.ndarray, np.ndarray) tuples giving the trajectory of each section
kwargs:
projection : ('x', 'z') tuple of two strings determining projection
The most efficient way of using this would be something like
from matplotlib.collections import PolyCollection import matplotlib.pyplot as plt cell = LFPy.Cell(morphology='PATH/TO/MORPHOLOGY') zips = [] for x, z in cell.get_idx_polygons(projection=('x', 'z')): zips.append(zip(x, z)) polycol = PolyCollection(zips, edgecolors='none', facecolors='gray') fig = plt.figure() ax = fig.add_subplot(111) ax.add_collection(polycol) ax.axis(ax.axis('equal')) plt.show()
-
get_intersegment_distance
(idx0=0, idx1=0)[source]¶ Return the Euclidean distance between midpoints of two segments with indices idx0 and idx1. Will return a float in unit of micrometers.
-
get_intersegment_vector
(idx0=0, idx1=0)[source]¶ Return the distance between midpoints of two segments with index idx0 and idx1. The argument returned is a vector [x, y, z], where x = self.xmid[idx1] - self.xmid[idx0] etc.
kwargs:
idx0 : int idx1 : int
-
get_pt3d_polygons
(projection=('x', 'z'))[source]¶ for each section create a polygon in the plane determined by keyword argument projection=(‘x’, ‘z’), that can be visualized using e.g., plt.fill()
Returned argument is a list of (x, z) tuples giving the trajectory of each section that can be plotted using PolyCollection
from matplotlib.collections import PolyCollection import matplotlib.pyplot as plt cell = LFPy.Cell(morphology='PATH/TO/MORPHOLOGY') zips = [] for x, z in cell.get_idx_polygons(projection=('x', 'z')): zips.append(zip(x, z)) polycol = PolyCollection(zips, edgecolors='none', facecolors='gray') fig = plt.figure() ax = fig.add_subplot(111) ax.add_collection(polycol) ax.axis(ax.axis('equal')) plt.show()
-
get_rand_idx_area_norm
(section='allsec', nidx=1, z_min=-10000, z_max=10000)[source]¶ Return nidx segment indices in section with random probability normalized to the membrane area of segment on interval [z_min, z_max]
kwargs:
section: str, string matching a section-name nidx: int, number of random indices z_min: float, depth filter z_max: float depth filter
-
get_rand_prob_area_norm
(section='allsec', z_min=-10000, z_max=10000)[source]¶ Return the probability (0-1) for synaptic coupling on segments in section sum(prob)=1 over all segments in section. Prob. determined by area.
kwargs:
section: str, string matching a section-name z_min: float, depth filter z_max: float depth filter
-
get_rand_prob_area_norm_from_idx
(idx=array([0]), z_min=-10000, z_max=10000)[source]¶ Return the normalized probability (0-1) for synaptic coupling on segments in idx-array. Normalised probability determined by area of segments.
kwargs:
idx : np.ndarray, dtype=int. array of segment indices z_min: float, depth filter z_max: float depth filter
-
insert_v_ext
(v_ext, t_ext)[source]¶ playback of some extracellular potential v_ext on each cell.totnseg compartments. Assumes that the “extracellular”-mechanism is inserted on each compartment.
Can be used to study ephaptic effects and similar
The inputs will be copied and attached to the cell object as cell.v_ext, cell.t_ext, and converted to (list of) neuron.h.Vector types, to allow playback into each compartment e_extracellular reference.
Can not be deleted prior to running cell.simulate()
Args:
v_ext : cell.totnsegs x t_ext.size np.array, unit mV t_ext : np.array, time vector of v_ext
Simple usage:
import LFPy import numpy as np import matplotlib.pyplot as plt #create cell cell = LFPy.Cell(morphology='morphologies/example_morphology.hoc') #time vector and extracellular field for every segment: t_ext = np.arange(cell.tstopms / cell.timeres_python+ 1) * cell.timeres_python v_ext = np.random.rand(cell.totnsegs, t_ext.size)-0.5 #insert potentials and record response: cell.insert_v_ext(v_ext, t_ext) cell.simulate(rec_imem=True, rec_vmem=True) fig = plt.figure() ax1 = fig.add_subplot(311) ax2 = fig.add_subplot(312) ax3 = fig.add_subplot(313) eim = ax1.matshow(np.array(cell.v_ext), cmap='spectral') cb1 = fig.colorbar(eim, ax=ax1) cb1.set_label('v_ext') ax1.axis(ax1.axis('tight')) iim = ax2.matshow(cell.imem, cmap='spectral') cb2 = fig.colorbar(iim, ax=ax2) cb2.set_label('imem') ax2.axis(ax2.axis('tight')) vim = ax3.matshow(cell.vmem, cmap='spectral') ax3.axis(ax3.axis('tight')) cb3 = fig.colorbar(vim, ax=ax3) cb3.set_label('vmem') ax3.set_xlabel('tstep') plt.show()
-
set_point_process
(idx, pptype, record_current=False, **kwargs)[source]¶ Insert pptype-electrode type pointprocess on segment numbered idx on cell object, with keyword arguments according to types SEClamp, VClamp, IClamp, SinIClamp, ChirpIClamp.
Arguments:
idx : int pptype : str record_current : bool kwargs : arguments passed on from class StimIntElectrode
-
set_pos
(xpos=0, ypos=0, zpos=0)[source]¶ Move the cell geometry so that midpoint of soma section is in (xpos, ypos, zpos). If no soma pos, use the first segment
-
set_rotation
(x=None, y=None, z=None)[source]¶ Rotate geometry of cell object around the x-, y-, z-axis in that order. Input should be angles in radians.
using rotation matrices, takes dict with rot. angles, where x, y, z are the rotation angles around respective axes. All rotation angles are optional.
Usage:
cell = LFPy.Cell(**kwargs) rotation = {'x' : 1.233, 'y' : 0.236, 'z' : np.pi} cell.set_rotation(**rotation)
-
set_synapse
(idx, syntype, record_current=False, record_potential=False, weight=None, **kwargs)[source]¶ Insert syntype (e.g. ExpSyn) synapse on segment with index idx,
Arguments:
idx : int syntype : str record_current : bool record_potential : bool weight : float kwargs : arguments passed on from class Synapse
-
simulate
(electrode=None, rec_imem=False, rec_vmem=False, rec_ipas=False, rec_icap=False, rec_isyn=False, rec_vmemsyn=False, rec_istim=False, rec_variables=[], variable_dt=False, atol=0.001, to_memory=True, to_file=False, file_name=None, dotprodcoeffs=None)[source]¶ This is the main function running the simulation of the NEURON model. Start NEURON simulation and record variables specified by arguments.
Arguments:
electrode: Either an LFPy.RecExtElectrode object or a list of such. If supplied, LFPs will be calculated at every time step and accessible as electrode.LFP. If a list of objects is given, accessible as electrode[0].LFP etc. rec_imem: If true, segment membrane currents will be recorded If no electrode argument is given, it is necessary to set rec_imem=True in order to calculate LFP later on. Units of (nA). rec_vmem: record segment membrane voltages (mV) rec_ipas: record passive segment membrane currents (nA) rec_icap: record capacitive segment membrane currents (nA) rec_isyn: record synaptic currents of from Synapse class (nA) rec_vmemsyn: record membrane voltage of segments with Synapse(mV) rec_istim: record currents of StimIntraElectrode (nA) rec_variables: list of variables to record, i.e arg=['cai', ] variable_dt: boolean, using variable timestep in NEURON atol: absolute tolerance used with NEURON variable timestep to_memory: only valid with electrode, store lfp in -> electrode.LFP to_file: only valid with electrode, save LFPs in hdf5 file format file_name: name of hdf5 file, '.h5' is appended if it doesnt exist dotprodcoeffs : list of N x Nseg np.ndarray. These arrays will at every timestep be multiplied by the membrane currents. Presumably useful for memory efficient csd or lfp calcs
-
class TemplateCell
¶
-
class
LFPy.
TemplateCell
(templatefile='LFPyCellTemplate.hoc', templatename='LFPyCellTemplate', templateargs=None, **kwargs)[source]¶ Bases:
LFPy.cell.Cell
This class allow using cell templates with some limitations
Arguments:
morphology : [str]: path to morphology file; templatefile : [str]: Cell template definition(s) templatename : [str]: Cell template-name used for this cell object templateargs : [str]: Arguments provided to template-definition v_init: [-65.]: initial potential passive: [True]/False: passive mechs are initialized if True Ra: [150.]: axial resistance rm: [30000]: membrane resistivity cm: [1.0]: membrane capacitance e_pas: [-65.]: passive mechanism reversal potential extracellular: [True]/False: switch for NEURON's extracellular mechanism timeres_NEURON: [0.1]: internal dt for NEURON simulation timeres_python: [0.1]: overall dt for python simulation tstartms: [0.]: initialization time for simulation <= 0 ms tstopms: [100.]: stop time for simulation > 0 ms nsegs_method: ['lambda100']/'lambda_f'/'fixed_length': nseg rule max_nsegs_length: [None]: max segment length for method 'fixed_length' lambda_f: [100]: AC frequency for method 'lambda_f' delete_sections: [True]: delete pre-existing section-references custom_code: [None]: list of model-specific code files ([.py/.hoc]) custom_fun: [None]: list of model-specific functions with args custom_fun_args: [None]: list of args passed to custom_fun functions pt3d: True/[False]: use pt3d-info of the cell geometries switch verbose: True/[False]: verbose output switch
Usage of TemplateCell class:
import LFPy cellParameters = { 'morphology' : 'path/to/morphology', 'templatefile' : 'path/to/template-file (.hoc) 'templatename' : 'templatename' 'templateargs' : None 'rm' : 30000, 'cm' : 1.0, 'Ra' : 150, 'timeres_NEURON' : 0.1, 'timeres_python' : 0.1, 'tstartms' : -50, 'tstopms' : 50, } cell = LFPy.TemplateCell(**cellParameters) cell.simulate()
class PointProcess
¶
-
class
LFPy.
PointProcess
(cell, idx, color='k', marker='o', record_current=False, **kwargs)[source]¶ Superclass on top of Synapse, StimIntElectrode, just to import and set some shared variables.
Arguments:
cell : LFPy.Cell object idx : index of segment color : color in plot (optional) marker : marker in plot (optional) record_current : Must be set True for recording of pointprocess currents kwargs : pointprocess specific variables passed on to cell/neuron
class Synapse
¶
-
class
LFPy.
Synapse
(cell, idx, syntype, color='r', marker='o', record_current=False, **kwargs)[source]¶ Bases:
LFPy.pointprocess.PointProcess
The synapse class, pointprocesses that spawn membrane currents. See http://www.neuron.yale.edu/neuron/static/docs/help/neuron/neuron/mech.html#pointprocesses for details, or corresponding mod-files.
This class is meant to be used with synaptic mechanisms, giving rise to currents that will be part of the membrane currents.
Usage:
#!/usr/bin/env python import LFPy import pylab as pl pl.interactive(1) cellParameters = { 'morphology' : 'morphologies/L5_Mainen96_LFPy.hoc', 'tstopms' : 50, } cell = LFPy.Cell(**cellParameters) synapseParameters = { 'idx' : cell.get_closest_idx(x=0, y=0, z=800), 'e' : 0, # reversal potential 'syntype' : 'ExpSyn', # synapse type 'tau' : 2, # syn. time constant 'weight' : 0.01, # syn. weight 'record_current' : True # syn. current record } synapse = LFPy.Synapse(cell, **synapseParameters) synapse.set_spike_times(pl.array([10, 15, 20, 25])) cell.simulate(rec_isyn=True) pl.subplot(211) pl.plot(cell.tvec, synapse.i) pl.title('Synapse current (nA)') pl.subplot(212) pl.plot(cell.tvec, cell.somav) pl.title('Somatic potential (mV)')
-
set_spike_times
(sptimes=array([], dtype=float64))[source]¶ Set the spike times explicitly using numpy arrays
-
set_spike_times_w_netstim
(noise=1.0, start=0.0, number=1000.0, interval=10.0, seed=1234.0)[source]¶ Generate a train of pre-synaptic stimulus times by setting up the neuron NetStim object associated with this synapse
kwargs:
noise : float in [0, 1] Fractional randomness, from deterministic to intervals that drawn from negexp distribution (Poisson spiketimes). start : float ms, (most likely) start time of first spike number : int (average) number of spikes interval : float ms, (mean) time between spikes seed : float random seed value
-
class StimIntElectrode
¶
-
class
LFPy.
StimIntElectrode
(cell, idx, pptype='SEClamp', color='p', marker='*', record_current=False, **kwargs)[source]¶ Bases:
LFPy.pointprocess.PointProcess
Class for NEURON point processes, ie VClamp, SEClamp and ICLamp, SinIClamp, ChirpIClamp with arguments. Electrode currents go here. Membrane currents will no longer sum to zero if these mechanisms are used.
Refer to NEURON documentation @ neuron.yale.edu for kwargs
Usage example:
#/usr/bin/python import LFPy import pylab as pl pl.interactive(1) #define a list of different electrode implementations from NEURON pointprocesses = [ { 'idx' : 0, 'record_current' : True, 'pptype' : 'IClamp', 'amp' : 1, 'dur' : 20, 'delay' : 10, }, { 'idx' : 0, 'record_current' : True, 'pptype' : 'VClamp', 'amp[0]' : -65, 'dur[0]' : 10, 'amp[1]' : 0, 'dur[1]' : 20, 'amp[2]' : -65, 'dur[2]' : 10, }, { 'idx' : 0, 'record_current' : True, 'pptype' : 'SEClamp', 'dur1' : 10, 'amp1' : -65, 'dur2' : 20, 'amp2' : 0, 'dur3' : 10, 'amp3' : -65, }, ] #create a cell instance for each electrode for pointprocess in pointprocesses: cell = LFPy.Cell(morphology='morphologies/L5_Mainen96_LFPy.hoc') stimulus = LFPy.StimIntElectrode(cell, **pointprocess) cell.simulate(rec_istim=True) pl.subplot(211) pl.plot(cell.tvec, stimulus.i, label=pointprocess['pptype']) pl.legend(loc='best') pl.title('Stimulus currents (nA)') pl.subplot(212) pl.plot(cell.tvec, cell.somav, label=pointprocess['pptype']) pl.legend(loc='best') pl.title('Somatic potential (mV)')
class RecExtElectrodeSetup
¶
-
class
LFPy.
RecExtElectrodeSetup
(cell=None, sigma=0.3, x=array([0]), y=array([0]), z=array([0]), N=None, r=None, n=None, shape=None, r_z=None, perCellLFP=False, method='linesource', color='g', marker='o', from_file=False, cellfile=None, verbose=False, seedvalue=None, **kwargs)[source]¶ RecExtElectrode superclass. If (optional) cell argument is given then the it is imported, otherwise the cell argument has to be passed later on to calc_lfp. The argument cell can be an LFPy.cell.Cell or LFPy.templatecell.TemplateCell object Keyword arguments determine properties of later LFP-calculations
Arguments:
cell : object, LFPy.cell.Cell or LFPy.templatecell.TemplateCell sigma : float, extracellular conductivity x, y, z : np.ndarray, coordinates or arrays of coordinates. Must be same length N : np.ndarray, Normal vector [x, y, z] of contact surface, default None r : float, radius of contact surface, default None n : int, if N is not None and r > 0, the number of points to use for each contact point in order to calculate average shape : str, 'circle'/'square' (default 'circle') defines the contact point shape If 'circle' r is the radius, if 'square' r is the side length method : str, ['linesource']/'pointsource'/'som_as_point' switch color : str, color of electrode contact points in plots marker : str, marker of electrode contact points in plots from_file : Bool, if True, load cell object from file cellfile : str, path to cell pickle verbose : Bool, Flag for verbose output seedvalue : int, rand seed when finding random position on contact with r >0
class RecExtElectrode
¶
-
class
LFPy.
RecExtElectrode
(cell=None, sigma=0.3, x=array([0]), y=array([0]), z=array([0]), N=None, r=None, n=0, shape=None, r_z=None, perCellLFP=False, method='linesource', color='g', marker='o', from_file=False, cellfile=None, verbose=False, seedvalue=None, **kwargs)[source]¶ Bases:
LFPy.recextelectrode.RecExtElectrodeSetup
RecExtElectrode class with inheritance from LFPy.RecExtElectrodeSetup able to actually calculate local field potentials from LFPy.Cell objects.
Usage:
import numpy as np import import matplotlib.pyplot as plt import LFPy N = np.empty((16, 3)) for i in xrange(N.shape[0]): N[i,] = [1, 0, 0] #normal vec. of contacts electrodeParameters = { #parameters for RecExtElectrode class 'sigma' : 0.3, #Extracellular potential 'x' : np.zeros(16)+25, #Coordinates of electrode contacts 'y' : np.zeros(16), 'z' : np.linspace(-500,1000,16), 'n' : 20, 'r' : 10, 'N' : N, } cellParameters = { 'morphology' : 'L5_Mainen96_LFPy.hoc', # morphology file 'rm' : 30000, # membrane resistivity 'cm' : 1.0, # membrane capacitance 'Ra' : 150, # axial resistivity 'timeres_NEURON' : 2**-4, # dt for NEURON sim. 'timeres_python' : 2**-4, # dt for python output 'tstartms' : -50, # start t of simulation 'tstopms' : 50, # end t of simulation } cell = LFPy.Cell(**cellParameters) synapseParameters = { 'idx' : cell.get_closest_idx(x=0, y=0, z=800), # compartment 'e' : 0, # reversal potential 'syntype' : 'ExpSyn', # synapse type 'tau' : 2, # syn. time constant 'weight' : 0.01, # syn. weight 'record_current' : True # syn. current record } synapse = LFPy.PointProcessSynapse(cell, **synapseParameters) synapse.set_spike_times(cell, np.array([10, 15, 20, 25])) cell.simulate() electrode = LFPy.RecExtElectrode(cell, **electrodeParameters) electrode.calc_lfp() plt.matshow(electrode.LFP)
submodule lfpcalc
¶
Copyright (C) 2012 Computational Neuroscience Group, NMBU.
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
-
LFPy.lfpcalc.
calc_lfp_choose
()¶ Determine which method to use, line-source for soma default
kwargs:
cell: LFPy.Cell or LFPy.TemplateCell instance x : double, extracellular position, x-axis y : double, extracellular position, y-axis z : double, extracellular position, z-axis sigma : double, extracellular conductivity r_limit : [None]/float/np.ndarray: minimum distance to source current timestep : [None]/int, calculate LFP at this timestep t_indices : [None]/np.ndarray, calculate LFP at specific timesteps method=['linesource']/'pointsource'/'som_as_point' switch for choosing underlying methods
-
LFPy.lfpcalc.
calc_lfp_linesource
()¶ Calculate electric field potential using the line-source method, all compartments treated as line sources, even soma.
kwargs:
cell: LFPy.Cell or LFPy.TemplateCell instance x : double, extracellular position, x-axis y : double, extracellular position, y-axis z : double, extracellular position, z-axis sigma : double, extracellular conductivity r_limit : [None]/float/np.ndarray: minimum distance to source current timestep : [None]/int, calculate LFP at this timestep t_indices : [None]/np.ndarray, calculate LFP at specific timesteps
-
LFPy.lfpcalc.
calc_lfp_pointsource
()¶ Calculate local field potentials using the point-source equation on all compartments
kwargs:
cell: LFPy.Cell or LFPy.TemplateCell instance x : double, extracellular position, x-axis y : double, extracellular position, y-axis z : double, extracellular position, z-axis sigma : double, extracellular conductivity r_limit : [None]/float/np.ndarray: minimum distance to source current timestep : [None]/int, calculate LFP at this timestep t_indices : [None]/np.ndarray, calculate LFP at specific timesteps
-
LFPy.lfpcalc.
calc_lfp_som_as_point
()¶ Calculate electric field potential using the line-source method, soma is treated as point/sphere source
kwargs:
cell: LFPy.Cell or LFPy.TemplateCell instance x : double, extracellular position, x-axis y : double, extracellular position, y-axis z : double, extracellular position, z-axis sigma : double, extracellular conductivity r_limit : [None]/float/np.ndarray: minimum distance to source current timestep : [None]/int, calculate LFP at this timestep t_indices : [None]/np.ndarray, calculate LFP at specific timesteps
submodule tools
¶
Copyright (C) 2012 Computational Neuroscience Group, NMBU.
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
-
LFPy.tools.
noise_brown
(ncols, nrows=1, weight=1, filter=None, filterargs=None)[source]¶ Return 1/f^2 noise of shape(nrows, ncols obtained by taking the cumulative sum of gaussian white noise, with rms weight.
If filter is not None, this function will apply the filter coefficients obtained by:
>>> b, a = filter(**filterargs) >>> signal = scipy.signal.lfilter(b, a, signal)
submodule inputgenerators
¶
Copyright (C) 2012 Computational Neuroscience Group, NMBU.
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
-
LFPy.inputgenerators.
get_normal_input_times
(n, mu, sigma, tstart, tstop)[source]¶ Generates n normal-distributed prosesses with mean mu and deviation sigma
-
LFPy.inputgenerators.
get_normal_spike_times
(nsyn, mu, sigma, tstart, tstop)[source]¶ Generate nsyn normal-distributed processes with mean mu and deviation sigma
-
LFPy.inputgenerators.
get_rand_spike_times
(synpos, nspikes, tstart, tstop)[source]¶ Return synpos times nspikes random spike times on the interval [tstart, tstop]
-
LFPy.inputgenerators.
stationary_gamma
(tstart, tstop, k=2, theta=10, tmin=-1000.0, tmax=1000000.0)[source]¶ Generate spiketimes with interspike interval statistics according to gamma-distribution with ‘shape’ k and ‘scale’ theta between tstart and tstop. Spiketimes from tmin up to tmax is calculated, times between 0 and tstop are returned
submodule run_simulation
¶
Copyright (C) 2012 Computational Neuroscience Group, NMBU.
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.