Prosthetic electrostimulation for information flow repair in a neocortical simulation (Kerr 2012)

 Download zip file   Auto-launch 
Help downloading and running models
This model is an extension of a model (<a href="">138379</a>) recently published in Frontiers in Computational Neuroscience. This model consists of 4700 event-driven, rule-based neurons, wired according to anatomical data, and driven by both white-noise synaptic inputs and a sensory signal recorded from a rat thalamus. Its purpose is to explore the effects of cortical damage, along with the repair of this damage via a neuroprosthesis.
1 . Kerr CC, Neymotin SA, Chadderdon GL, Fietkiewicz CT, Francis JT, Lytton WW (2012) Electrostimulation as a prosthesis for repair of information flow in a computer model of neocortex IEEE Transactions on Neural Systems & Rehabilitation Engineering 20(2):153-60 [PubMed]
Model Information (Click on a link to find other models with that property)
Model Type: Realistic Network;
Brain Region(s)/Organism: Neocortex;
Cell Type(s): Neocortex V1 pyramidal corticothalamic L6 cell; Neocortex V1 pyramidal intratelencephalic L2-5 cell; Neocortex V1 interneuron basket PV cell; Neocortex fast spiking (FS) interneuron; Neocortex spiny stellate cell;
Channel(s): I Chloride; I Sodium; I Potassium;
Gap Junctions:
Receptor(s): GabaA; AMPA; NMDA; Gaba;
Transmitter(s): Gaba; Glutamate;
Simulation Environment: NEURON;
Model Concept(s): Activity Patterns; Deep brain stimulation; Information transfer; Brain Rhythms;
Implementer(s): Lytton, William [billl at]; Neymotin, Sam [samn at]; Kerr, Cliff [cliffk at];
Search NeuronDB for information about:  Neocortex V1 pyramidal corticothalamic L6 cell; Neocortex V1 pyramidal intratelencephalic L2-5 cell; Neocortex V1 interneuron basket PV cell; GabaA; AMPA; NMDA; Gaba; I Chloride; I Sodium; I Potassium; Gaba; Glutamate;
infot.mod *
intf6_.mod *
intfsw.mod *
misc.mod *
nstim.mod *
staley.mod *
stats.mod *
vecst.mod *
decmat.hoc *
decnqs.hoc *
default.hoc *
drline.hoc *
infot.hoc *
local.hoc *
misc.h *
ratlfp.dat *
setup.hoc *
simctrl.hoc *
spkts.hoc *
staley.hoc *
stats.hoc *
stdgui.hoc *
syncode.hoc *
updown.hoc *
xgetargs.hoc *

Python module for using other Python modules via the hoc interpreter.
You will need to begin with this:

	objref p
	p = new PythonObject()
	nrnpython("import pyhoc")

Many of the usage examples expect an input consisting of one or two
long vectors representing time series, here called A and B. The following
code generates such vectors. Simply copy and paste this code into the
NEURON interpreter, then copy and paste the example code for a given
module, and the script should run.

	objref A, B, r1, r2 // Initialize variables
	r1 = r2 = new Random() // Create random variables
	r1.normal(0,1) r2.normal(0,1) // Set random variables to normal distribution
	A = new Vector(20000) // Initialize the first time series
	A.indgen(0.01) A.sin(5,0) A.addrand(r1) // Populate it
	B = new Vector(20000) // Initialize the second time series
	B.indgen(0.01) B.sin(5,15) B.addrand(r2) // Populate it

Alternatively, if you're working in an intfcol-based environment,
a more realistic pair of time series (i.e. actual LFP time series)
can be obtained as follows:

	objref A, B

Version: 2011apr28

def bsmart(nqx1,nqx2,ntrls=1,npts=-1,p=12,fs=200,freq=100): # Set defaults for everything except the data x
    This is the wrapper for the BSMART code.
    Usage is similar to
    	x1 = vector representing first time series
        x2 = vector representing second time series
    	ntrls = number of trials in the time series (best set to 1)
    	npts = length of input (if set to -1, is calculated automatically)
    	p = polynomial order for fitting (lower = smoother fit)
    	fs = sampling rate for the time series (in Hz)
    	freq = maximum frequency to be returned (usually fs/2)
    where grangernqs has the following fields:
        F -- vector of frequencies for each of the following
        pp1 -- power spectrum for first time series
        pp2 -- power spectrum for second time series
        cohe -- coherence between the two time series
        Fx2y -- causality from first to second time series
        Fy2x -- causality from second to first time series
        Fxy -- nondirectional causality
    Example usage from NEURON is as follows: 
        objref output
       "Fx2y","F") // Strong causality"Fy2x","F") // No causality

    Version: 2011apr21
## Import packages
    from numpy import array, zeros, size # Shorten useful functions
    from bsmart import pwcausalr
    from neuron import h
## Initialize data vectors
    tmp1=array(nqx1) # Convert NQS table to Numpy arrays
    if npts==-1: npts=size(tmp1,0) # Reset npts if needed
    x=array(zeros((2,npts))) # Store both time series in one matrix
## Do the analysis
    F,pp,cohe,Fx2y,Fy2x,Fxy=pwcausalr(x,int(ntrls),int(npts),int(p),fs,int(freq)); # Do the analysis
## Initialize hoc objects
    h('objref grangernqs') # Initialize NQS object
    h('objref F, pp1, pp2, cohe, Fx2y, Fy2x, Fxy, tmp') # Initialize vectors
    h('F   =new Vector()')
    h('pp1 =new Vector()')
    h('pp2 =new Vector()')
    h('cohe=new Vector()')
    h('Fx2y=new Vector()')
    h('Fy2x=new Vector()')
    h('Fxy=new Vector()')
## Convert from Python to hoc
    h.tmp=F        ; h('F=F.from_python(tmp)')
    h.tmp=pp[0,:]  ; h('pp1=pp1.from_python(tmp)')
    h.tmp=pp[1,:]  ; h('pp2=pp2.from_python(tmp)')
    h.tmp=cohe[0,:]; h('cohe=cohe.from_python(tmp)')
    h.tmp=Fx2y[0,:]; h('Fx2y=Fx2y.from_python(tmp)')
    h.tmp=Fy2x[0,:]; h('Fy2x=Fy2x.from_python(tmp)')
    h.tmp=Fxy[0,:] ; h('Fxy=Fxy.from_python(tmp)')
## Convert from hoc to Python
    h('grangernqs=new NQS("F","pp1","pp2","cohe","Fx2y","Fy2x","Fxy")')
    h('grangernqs.setcol("F",F)') # Save the data to the NQS table
    return grangernqs
def downsample(olddata,oldrate=10000,newrate=200): # Too different from the original code to even call.
    This function downsamples a given vector or matrix.
        newdata = downsampled data
        olddata = data at original sampling rate
        origrate = original sampling rate (default 10 kHz)
        newrate = desired sampling rate (default 200 Hz)
    If olddata has multiple columns, these are assumed to be different time
    series. Thus, an original matrix of N rows by M columns will be downsampled
    to a matrix of N' rows and M columns, where N' = N*origrate/newrate.
    Example usage from NEURON is as follows:
        objref output1, output2
        output1=p.pyhoc.downsample(A) // 
        A.size() // = 20000 -- original vector size
        output1.size() // = 400 -- downsampled by a factor of 50
        output2.size() // = 2000 -- downsampled by a factor of 10
    Version: 2011apr28
## Load packages
    from scipy import array, shape, size, reshape, zeros
    from neuron import h
## Convert data
    ratio=oldrate/float(newrate) # Calculate ratio of sampling rates
    olddata=array(olddata) # Make sure it's an array
    if olddata.ndim==1: olddata=reshape(olddata,(size(olddata,0),1)) # Turn vector into an array
    rows,cols=shape(olddata) # Find out how many rows and columns there are
    newrows=int(rows/ratio); # Calculate how many rows the new file will have
## Perform downsampling
    newdata=zeros((newrows,cols)); # Initialize new array
    for i in range(cols): # Loop over time series
        for j in range(newrows): # Loop over new time points
            tstart=int(j*ratio) # This is the starting time of the average
            tfinish=int((j+1)*ratio) # This is the finishing time of the average
            newdata[j,i]=olddata[tstart:tfinish,i].mean(); # Calculate mean across the original time points
## Convert from PythonObjet to hoc array
    h('objref tmpinput, tmpoutput')
    h('tmpoutput = new Vector()')
    return output

def spklfp():
    This function takes the data structures generated by an
    intfcol-based simulation and uses them to plot every
    quantity of general interest: a spike raster, per-cell
    firing rate histogram, population firing rates, raw LFP
    time series, and LFP spectra.
    Usage is as follows:
    It requires the following to run:
    	- Matplotlib 1.0.1 or later
    	- NQS table storing LFPs called "nqLFP"
    	- NQS table storing spikes called "snq"
    Version: 2011apr28
    print 'Converting data...'
    flattendata() # Convert data
    import spklfp # Run code
    return 0

def spectrogram(ts,fs=200,window=2,maxfreq=50,tsmooth=2,fsmooth=2):
    This function takes a given time series and turns it into a spectrogram
    (i.e. a 3-D plot where one axis is time, one is frequency, and one is
        ts = time series to be spectrogrammed
        fs = sampling rate (in Hz)
        window = length of window for computing spectra (in s)
        maxfreq = maximum frequency to plot (in Hz)
        tsmooth = amount of smoothing to do along time axis
        fsmooth = amount of smoothing to do along frequency axis
    Example usage from NEURON is as follows:
    Version: 2011apr28
    from spectrogram import plotspectrogram
    from pylab import array
    return 0

def viewlfps(ncols=1,trimdata=1,fs=200,tmax=0,fmax=50,fftsmooth=50,mtpar=4,order=12):
    This function provides an interative way of visualizing LFPs for a particular
    simulation -- it allows you to visualize LFP time series or spectra, the latter
    calculated in one of three ways (plain FFT, multitaper spectrum, or auto-
    regressive fitting via BSMART). 
    The non-hoc version allows for the comparison of multiple columns and multiple 
    simulations; however, due to the limitations of the hoc interpreter, only one
    simulation can be viewed at a time in this version. The simulation must have 
    been run in such a way as to generate "nqLFP" (i.e. it must be an intfcol-based
    simulation) , which is then read in by this script.
        ncols = number of columns; M = ncols * number of layers (default 1)
        trimdata = number of seconds' worth of data to trim off each end of the LFP (default 1)
        fs = data sampling rate (default 200 Hz)
        tmax = maximum time to display; set to 0 for all (default 0)
        fmax = maximum frequency to display; set to 0 for all (default 50)
        fftsmooth = the amount of smoothing to do on the FFT (default 50)
        mtpar = the window size for the multitaper method (default 4)
        order = polynomial order for BSMART (default 12)
    Example usage from NEURON is as follows:
    	- NQS table storing LFPs called "nqLFP"
    Version: 2011apr28
    from pylab import loadtxt, shape
    from viewlfps import plotlfps
## Define options
    fs=200.0 # Sampling rate in Hz
    tmax=0 # Maximum time in s
    fmax=50 # Maximum frequency in Hz
    fftsmooth=50 # How much to smooth the raw FFT -- 50-100 is good
    mtpar=3.5 # The parameter for the multitaper method -- 2-4 is good
    order=10 # The polynomial order for BSMART -- 10-30 is good
## Convert data
    print 'Converting data...'
    flattendata() # Convert data
## Import data
    print 'Importing data...'
    filenames=['/tmp/pyhoc-lfp.txt'] # Originally /home/cliffk/bill/ckintf/data/1102/juemo/lfp; alternative '/home/cliffk/bill/ckintf/data/1104/07-chrislfp1.txt'
    killdata=trimdata*fs # How much data to cut off each end
    if npts<=killdata*3: # If killdata is too big for the length of the data, make it smaller
        print 'Warning: trimming data would have result in nothing left!'
    alldata[0]=alldata[0][killdata:-killdata-1,:] # Remove bad data
## Plot LFPs
    print '...done.'
    return 0

def flattendata():
    This function is a hoc script to convert NQS tables generated by an
    intfcol-based simulation to a form readable by Python. Not to be used 
    directly by the user. Based on $ckintf/batch.hoc.
    Version: 2011apr28
    from neuron import h
    from subprocess import call
    h('oldhz=nqLFP.cob.v.size/tstop*1000 // Original sampling rate; *1000 because tstop is in ms')
    h('newhz=200 // The new frequency to sample at, in Hz')
    h('ratio=oldhz/newhz // Calculate the ratio betwen the old and new sampling rates')
    h('npts=tstop/1000*newhz // Number of points in the resampled time seris')
    h('nlayers=nqLFP.m // Number of layers (usually 5 -- 2/3, 4, 5, 6, all)')
    h('objref tempvec // Temporary place to store NQS column as a vector')
    h('objref tempstr // Name of the NQS column being selected')
    h('objref storelfp // Create matrix to store results in')
    h('storelfp = new Matrix(npts, nlayers*numcols) // Combine layers/columns into one dimension')
    h('count=-1 // Set column of storelfp to zero')
    h('for i=0,numcols-1 { for j=0,nlayers-1 { count+=1 tempstr=nqLFP[i].s[j] tempvec=nqLFP[i].getcol(tempstr.s) for k=0,npts-1 {storelfp.x[k][count]=tempvec.mean(k*ratio,(k+1)*ratio-1)}}}')
    h('objref fobj')
    h('fobj = new File("/tmp/pyhoc-lfp.txt")')
    h('storelfp.fprint(fobj,"%10.1f") // Its usually in the thousands so one d.p. should do')
    h('skipsnq=0 // flag to create NQS with spike times, one per column')
    h('initAllMyNQs() // setup of NQS objects with spike/other information')
    h('objref storespikes, tmpt, tmpid, tmptype, tmpcol // Initialize vectors and matrices -- the tmp vectors are for storing parts of the NQS arrays')
    h('totalnumberofspikes=0 // Calculate the total number of spikes generated across all columns')
    h('for i=0,numcols-1 totalnumberofspikes+=snq[i].cob.v.size')
    h('storespikes = new Matrix(totalnumberofspikes, 4) // Four columns: spike time, cell ID, cell type, and spike time')
    h('count=-1 // Initialize row count')
    h('for i=0,numcols-1 { tmpt=snq[i].getcol("t") tmpid=snq[i].getcol("id") tmptype=snq[i].getcol("type") tmpcol=snq[i].getcol("col") for j=0,snq[i].cob.v.size-1 { count+=1 storespikes.x[count][0]=tmpt.x[j] storespikes.x[count][1]=tmpid.x[j] storespikes.x[count][2]=tmptype.x[j] storespikes.x[count][3]=tmpcol.x[j]}}')
    h('objref fobj2')
    h('fobj2 = new File("/tmp/pyhoc-spk.txt")')
    h('storespikes.fprint(fobj2,"%6.0f") // All quantities are integers, so this should be fine')

Kerr CC, Neymotin SA, Chadderdon GL, Fietkiewicz CT, Francis JT, Lytton WW (2012) Electrostimulation as a prosthesis for repair of information flow in a computer model of neocortex IEEE Transactions on Neural Systems & Rehabilitation Engineering 20(2):153-60[PubMed]

References and models cited by this paper

References and models that cite this paper

Adesnik H, Scanziani M (2010) Lateral competition for cortical space by layer-specific horizontal circuits. Nature 464:1155-60 [PubMed]

Binzegger T, Douglas RJ, Martin KA (2004) A quantitative map of the circuit of cat primary visual cortex. J Neurosci 24:8441-53 [PubMed]

Carnevale NT, Hines ML (2006) The NEURON Book

Cui J, Xu L, Bressler SL, Ding M, Liang H (2008) BSMART: a Matlab-C toolbox for analysis of multichannel neural time series. Neural Netw 21:1094-104 [PubMed]

Francis JT, Xu S, Chapin JK (2008) Proprioceptive and cutaneous representations in the rat ventral posterolateral thalamus. J Neurophysiol 99:2291-304 [PubMed]

Gisiger T, Boukadoum M (2011) Mechanisms Gating the Flow of Information in the Cortex: What They Might Look Like and What Their Uses may be. Front Comput Neurosci 5:1-304 [PubMed]

Hines ML, Carnevale NT (2001) NEURON: a tool for neuroscientists. Neuroscientist 7:123-35 [Journal] [PubMed]

   Spatial gridding and temporal accuracy in NEURON (Hines and Carnevale 2001) [Model]

Kamiński M, Ding M, Truccolo WA, Bressler SL (2001) Evaluating causal relations in neural systems: granger causality, directed transfer function and statistical assessment of significance. Biol Cybern 85:145-57 [PubMed]

Lefort S, Tomm C, Floyd Sarria JC, Petersen CC (2009) The excitatory neuronal network of the C2 barrel column in mouse primary somatosensory cortex. Neuron 61:301-16 [PubMed]

Lizier JT, Heinzle J, Horstmann A, Haynes JD, Prokopenko M (2011) Multivariate information-theoretic measures reveal directed information structure and task relevant changes in fMRI connectivity. J Comput Neurosci 30:85-107

Lloyd KG, Davidson L, Hornykiewicz O (1975) The neurochemistry of Parkinson's disease: effect of L-dopa therapy. J Pharmacol Exp Ther 195:453-64 [PubMed]

Lytton WW, Neymotin SA, Hines ML (2008) The virtual slice setup. J Neurosci Methods 171:309-15 [Journal] [PubMed]

   The virtual slice setup (Lytton et al. 2008) [Model]

Lytton WW, Omurtag A (2007) Tonic-clonic transitions in computer simulation. J Clin Neurophysiol 24:175-81 [PubMed]

   Tonic-clonic transitions in a seizure simulation (Lytton and Omurtag 2007) [Model]

Lytton WW, Omurtag A, Neymotin SA, Hines ML (2008) Just in time connectivity for large spiking networks Neural Comput 20(11):2745-56 [Journal] [PubMed]

   JitCon: Just in time connectivity for large spiking networks (Lytton et al. 2008) [Model]

Lytton WW, Stewart M (2005) A rule-based firing model for neural networks Int J Bioelectromagn 7:47-50

Lytton WW, Stewart M (2006) Rule-based firing for network simulations. Neurocomputing 69:1160-1164

Meyer JS, Obara K, Muramatsu K (1993) Diaschisis. Neurol Res 15:362-6 [PubMed]

Neymotin SA, Jacobs KM, Fenton AA, Lytton WW (2011) Synaptic information transfer in computer models of neocortical columns. J Comput Neurosci. 30(1):69-84 [Journal] [PubMed]

   Synaptic information transfer in computer models of neocortical columns (Neymotin et al. 2010) [Model]

Quilodran R, Gariel MA, Markov NT, Falchier A, Vezoli J, Sallet J, Anderson JC, Dehay C, Doug (2008) Strong loops in the neocortex Society for Neuroscience Abstracts 853.4

Rasche D, Rinaldi PC, Young RF, Tronnier VM (2006) Deep brain stimulation for the treatment of various chronic pain syndromes. Neurosurg Focus 21:E8

Richman JS, Moorman JR (2000) Physiological time-series analysis using approximate entropy and sample entropy. Am J Physiol Heart Circ Physiol 278:H2039-49 [PubMed]

Schiff ND, Giacino JT, Kalmar K, Victor JD, Baker K, Gerber M, Fritz B, Eisenberg B, Biondi T (2007) Behavioural improvements with thalamic stimulation after severe traumatic brain injury. Nature 448:600-3 [PubMed]

Schroeder CE, Mehta AD, Foxe JJ (2001) Determinants and mechanisms of attentional modulation of neural processing. Front Biosci 6:D672-84

Shipp S (2005) The importance of being agranular: a comparative account of visual and motor cortex. Philos Trans R Soc Lond B Biol Sci 360:797-814 [PubMed]

Stefani A, Lozano AM, Peppe A, Stanzione P, Galati S, Tropepi D, Pierantozzi M, Brusa L, Scar (2007) Bilateral deep brain stimulation of the pedunculopontine and subthalamic nuclei in severe Parkinson's disease. Brain 130:1596-607 [PubMed]

Stoerig P, Cowey A (1997) Blindsight in man and monkey. Brain 120 ( Pt 3):535-59 [PubMed]

Traub RD, Contreras D, Cunningham MO, Murray H, Lebeau FE, Roopun A, Bibbig A, et al (2005) A single-column thalamocortical network model exhibiting gamma oscillations, sleep spindles and epileptogenic bursts J Neurophysiol 93(4):2194-232 [Journal] [PubMed]

   A single column thalamocortical network model (Traub et al 2005) [Model]
   Collection of simulated data from a thalamocortical network model (Glabska, Chintaluri, Wojcik 2017) [Model]

Van Essen DC, Anderson CH, Felleman DJ (1992) Information processing in the primate visual system: an integrated systems perspective. Science 255:419-23 [PubMed]

Von_monakow C (1914) Die Lokalisation im Grosshirn und der Abbau der Funktion durch kortikale Herde

Chadderdon GL, Mohan A, Suter BA, Neymotin SA, Kerr CC, Francis JT, Shepherd GM, Lytton WW (2014) Motor cortex microcircuit simulation based on brain activity mapping. Neural Comput 26:1239-62 [Journal] [PubMed]

   Motor cortex microcircuit simulation based on brain activity mapping (Chadderdon et al. 2014) [Model]

Dura-Bernal S, Li K, Neymotin SA, Francis JT, Principe JC, Lytton WW (2016) Restoring behavior via inverse neurocontroller in a lesioned cortical spiking model driving a virtual arm. Front. Neurosci. Neuroprosthetics 10:28 [Journal]

   Cortical model with reinforcement learning drives realistic virtual arm (Dura-Bernal et al 2015) [Model]

Dura-Bernal S, Neymotin SA, Kerr CC, Sivagnanam S, Majumdar A, Francis JT, Lytton WW (2017) Evolutionary algorithm optimization of biological learning parameters in a biomimetic neuroprosthesis. IBM Journal of Research and Development (Computational Neuroscience special issue) 61(2/3):6:1-6:14 [Journal]

   Motor system model with reinforcement learning drives virtual arm (Dura-Bernal et al 2017) [Model]

Dura-Bernal S, Zhou X, Neymotin SA, Przekwas A, Francis JT, Lytton WW (2015) Cortical Spiking Network Interfaced with Virtual Musculoskeletal Arm and Robotic Arm. Front Neurorobot 9:13 [Journal] [PubMed]

   Cortical model with reinforcement learning drives realistic virtual arm (Dura-Bernal et al 2015) [Model]

Kerr CC, Van Albada SJ, Neymotin SA, Chadderdon GL, Robinson PA, Lytton WW (2013) Cortical information flow in Parkinson's disease: a composite network-field model. Front Comput Neurosci 7:39:1-14 [Journal] [PubMed]

   Composite spiking network/neural field model of Parkinsons (Kerr et al 2013) [Model]

Neymotin SA, Chadderdon GL, Kerr CC, Francis JT, Lytton WW (2013) Reinforcement learning of 2-joint virtual arm reaching in a computer model of sensorimotor cortex Neural Computation 25(12):3263-93 [Journal] [PubMed]

   Sensorimotor cortex reinforcement learning of 2-joint virtual arm reaching (Neymotin et al. 2013) [Model]

Neymotin SA, Dura-Bernal S, Lakatos P, Sanger TD, Lytton WW (2016) Multitarget Multiscale Simulation for Pharmacological Treatment of Dystonia in Motor Cortex. Front Pharmacol 7:157 [Journal] [PubMed]

   Multitarget pharmacology for Dystonia in M1 (Neymotin et al 2016) [Model]

Neymotin SA, Lee H, Park E, Fenton AA, Lytton WW (2011) Emergence of physiological oscillation frequencies in a computer model of neocortex. Front Comput Neurosci 5:19-75 [Journal] [PubMed]

   Emergence of physiological oscillation frequencies in neocortex simulations (Neymotin et al. 2011) [Model]

Rowan MS, Neymotin SA, Lytton WW (2014) Electrostimulation to reduce synaptic scaling driven progression of Alzheimer's disease. Front Comput Neurosci 8:39 [Journal] [PubMed]

   Electrostimulation to reduce synaptic scaling driven progression of Alzheimers (Rowan et al. 2014) [Model]

(38 refs)