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

 Download zip file   Auto-launch 
Help downloading and running models
Accession:183014
We developed a 3-layer sensorimotor cortical network of consisting of 704 spiking model-neurons, including excitatory, fast-spiking and low-threshold spiking interneurons. Neurons were interconnected with AMPA/NMDA, and GABAA synapses. We trained our model using spike-timing-dependent reinforcement learning to control a virtual musculoskeletal human arm, with realistic anatomical and biomechanical properties, to reach a target. Virtual arm position was used to simultaneously control a robot arm via a network interface.
References:
1 . 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 [PubMed]
2 . 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
Model Information (Click on a link to find other models with that property)
Model Type: Realistic Network;
Brain Region(s)/Organism:
Cell Type(s): Neocortex M1 pyramidal pyramidal tract L5B cell; Neocortex M1 pyramidal intratelencephalic L2-5 cell; Neocortex M1 interneuron basket PV cell; Neocortex fast spiking (FS) interneuron; Neostriatum fast spiking interneuron; Neocortex spiking regular (RS) neuron; Neocortex spiking low threshold (LTS) neuron;
Channel(s):
Gap Junctions:
Receptor(s): GabaA; AMPA; NMDA;
Gene(s):
Transmitter(s): Gaba; Glutamate;
Simulation Environment: NEURON; Python (web link to model);
Model Concept(s): Synaptic Plasticity; Learning; Reinforcement Learning; STDP; Reward-modulated STDP; Sensory processing; Motor control;
Implementer(s): Neymotin, Sam [samn at neurosim.downstate.edu]; Dura, Salvador [ salvadordura at gmail.com];
Search NeuronDB for information about:  Neocortex M1 pyramidal intratelencephalic L2-5 cell; Neocortex M1 pyramidal pyramidal tract L5B cell; Neocortex M1 interneuron basket PV cell; GabaA; AMPA; NMDA; Gaba; Glutamate;
/
arm2dms_modeldb
mod
msarm
stimdata
README.html
analyse_funcs.py
analysis.py
armGraphs.py
arminterface_pipe.py
basestdp.hoc
bicolormap.py
boxes.hoc *
bpf.h *
col.hoc
colors.hoc *
declist.hoc *
decmat.hoc *
decnqs.hoc *
decvec.hoc *
default.hoc *
drline.hoc *
filtutils.hoc *
grvec.hoc
hinton.hoc *
hocinterface.py
infot.hoc *
init.hoc
intfsw.hoc *
labels.hoc
load.hoc
load.py
local.hoc *
main.hoc
main_demo.hoc
main_neurostim.hoc
misc.h *
misc.py *
msarm.hoc
network.hoc
neuroplot.py *
neurostim.hoc
nload.hoc
nqs.hoc *
nqsnet.hoc *
nrnoc.hoc
params.hoc
perturb.hoc
python.hoc
pywrap.hoc *
run.hoc
runbatch_neurostim.py
runsim_neurostim
samutils.hoc *
saveoutput.hoc
saveoutput2.hoc
setup.hoc *
sim.hoc
sim.py
sim_demo.py
simctrl.hoc *
stats.hoc *
stim.hoc
syncode.hoc *
units.hoc *
vector.py
xgetargs.hoc *
                            
// $Id: run.hoc,v 1.101 2012/06/25 04:16:27 samn Exp $

declare("wf1","o[1]","wf2","o[1]","wrec","o[1]")
sprint(tstr,"o[%d]",numcols)
declare("vit",tstr,"nqLFP",tstr,"SPKSZ",25e6*tstop/20e3,"vitdp","o[1]","maxrec",1)
declare("recLFP",0) // whether to record LFP
declare("recSPKS",1) // whether to record spikes
declare("recV",0) // whether to record cell voltages

declare("nqa","o[1]") // NQS for recording arm joint angles/positions in time

method("local")
cvode.atol(1e-3) 
cvode.condition_order(1) // irrelevant to acells?
declare("sepflds",0) // whether to record separate fields for each layer
declare("initrands",1) // whether to initialize rands @ start
//declare("useNSLOC",1) //  replace drspk cells with NSLOC objects - declared in network.hoc

// rv2 is part of grvec.hoc, as a stub, so it can be
// customized to do something to vector before graphing it
proc rv2 () { }

//** draw lines between cell subpopulations
proc rasterlines () { local i,j,ct,c localobj o
  {g=graphItem c=CDX}
  for CDX=0,numcols-1 for ct=0,CTYPi-1 if(col[CDX].numc[ct]) {
    if(ct==DP) { 
      j=cedp.o(0).id
      i=cedp.o(cedp.count-1).id
    } else {
      i=col[CDX].ce.o(col[CDX].ixe[ct]).gid
      j=col[CDX].ce.o(col[CDX].ix[ct]).gid
    }
    drline(0,i,tstop,i,g,2,6)
    o=mdl2view(g,0.9,j+col[CDX].numc[ct]/4)
    g.label(0.9,o.x[3],CTYP.o(ct).s)
  }
  {CDX=c g.flush}
}
//** calls rasterlines, making sure g is set first
proc grlines () { local id
  if(numarg()>0) id=$1 else id=0
  {g=Graph[id] rasterlines()}
}

//** a
proc a () { local sh,sv
  if(g==nil)gg()
  if (!isobj(aa,"Graph")) aa=g else g=aa
  if (aa.view_count==0) aa=g
  sv=gnum gnum=ojtnum(g)
  graphItem=g
  sh=0 grv_.super=1 g.erase_all
  grv_.gveraseflag=0  grv_.gvmarkflag=grv_.super=1 // gnum=0
  gv(0,1+sh,2) // gv(1,3+sh,2) gv(2,2+sh,2)
  grv_.gvmarkflag=grv_.super=0
  gnum=sv
  rasterlines()
}

//** b
proc b () { rewt() time() a() }

//** initrr - for doing a rerun - not used right now
proc initrr () { 
  rewt()
  col.intf.global_init()
  NStim[0].global_init()
  vseed_stats(392426)
}

//** setMemb - nothing here
proc setMemb () {}

//** initMisc1
proc initMisc1 () { local i localobj xo
  vseed_stats(392426)
  col.intf.global_init()
  if(initrands) for i=0,numcols-1 col[i].cstim.initrands()
  //for ltr(xo,cedp) xo.soma.v = vrefrac_DRSPK //make sure DPD voltages initialized to same level 
}

objref ww // global ww for post processing
{wwht_INTF6=1 wwwid_INTF6=100}
//** wrecon - setup LFP recording, one LFP for each COLUMN
proc wrecon () { local cdx,ii,jj,x,ct,cnt localobj tl,vm
  tl=new List()
  for cdx=0,numcols-1 {
    nqsdel(nqLFP[cdx])
    if(sepflds) nqLFP[cdx]=new NQS("ES","EM","LFP") else {nqLFP[cdx]=new NQS(1) nqLFP[cdx].s[0].s="LFP"}
    {nqLFP[cdx].v.resize(tstop/vdt_INTF6) nqLFP[cdx].pad()}
    for ii=0,nqLFP[cdx].m-1 tl.append(nqLFP[cdx].v[ii])
  }
  {vm=new Vector(CTYPi) vm.x(ES)=0 vm.x(EM)=1}
  col.intf.initwrec(tl)
  for x=0,numcols-1 {
    //for case(&ct,ES,EM,&jj) {
    for case(&ct,EM,&jj) {
      for ii=col[x].ix[ct],col[x].ixe[ct] {
        if(sepflds) {
          col[x].ce.o(ii).wrc(vm.x(jj)+x*3)
          col[x].ce.o(ii).wrc(x*3)
        } else col[x].ce.o(ii).wrc(x)
      }
    }
  }
}

//** wrecoff - turn off LFP recording
proc wrecoff () { local ct
  for ctt(&ct) for ixt(ct) XO.wrec(0)
}

//** finishMisc - called @ end of run()
func finishMisc () { local ii localobj co
  for ltr(XO,printlist) if (isassigned(XO.o)) if (XO.o.fflag) XO.o.fini
  // if(0) panobj.pvplist(ofile,params,100) //dont save printlist for now
  col.intf.global_fini
  for ltr(co,lcol) co.intf.spkstats2(1)
  print "TMAX: ",tmax_INTF6
  return 1
}

//** snapsv() save after printlist items min-max to fixed dt
proc snapsv () { local a,vdt,min,max localobj v1,o
  grv_.bst(3,3)
  vdt=0.2 
  a=allocvecs(v1)
  v1.resize(tstop/vdt)
  for ltr(o,printlist) { 
    if (o.code!=3) continue
    v1.snap(o.vec,o.tvec,vdt)
    o.vec.copy(v1)
    o.pstep=vdt
    o.tvflag=0
  }
  if (numarg()==0) grv_.pvall()
  dealloc(a)
}

proc exeruncall () { for ltr(XO,printlist) if (XO.code==3) XO.tvflag=1 }
proc pvout2 () { snapsv(1) }

//printlist=new List()
if(printlist==nil)printlist=new List()
proc prlclr () { localobj ce,col,intf
  for ltr(XO,printlist) {
    if (isassigned(XO.o)) if (XO.o.fflag) XO.o.recclr
  }
  //  for ltr(XO,ce) XO.wrc(-1)
  for ltr(col,lcol) {
    {ce=col.ce intf=ce.o(0)}
    for ii=0,ce.count-1 ce.o(ii).wrc(-1)
    intf.wwfree(0)
  }
  printlist.remove_all
}

//** prl(recv,recs[,lvextra]) - setup recording in printlist
//$1  = whether to record any cell voltages, default off
//$2  = whether to record spike times, default on
//$o3 = extra cells to record. lv.o(0)=cell,lv.o(1)=param to record,etc.,optional
proc prl () { local a,x,cidx,ii,jj,offst,y,recv,recs,max localobj xo,lvextra,co,ce
  if(numarg()>0) recv=$1 else recv=0
  if(numarg()>1) recs=$2 else recs=1
  if(numarg()>2) lvextra=$o3 else lvextra=nil
  {offst=0 prlclr()}
  for ltr(co,lcol,&ii) {
    {ce=co.ce intf=co.intf}
    if (recs) { sprint(tstr,"%s_SPKS",co.name)
      if (intf.flag("jcn")) { // for use with jitcon()
        printlist.append((vit[ii]=new vitem(tstr,SPKSZ,1)))
        intf.jitrec(vit[ii].vec,vit[ii].tvec)
      } else {
        intf.jitrec() // clear jit recording
        for ltr(xo,ce,&y) {
          if (y==0) vit[ii]=new_printlist_nc(xo, xo.id, tstr) else {
            new_printlist_nc2(vit[ii], xo, xo.id)        }
        }
      }
    }
    if(co.numc[DP]>0) {
      printlist.append(vitdp=new vitem("DPSPKS",100))
      if(useNSLOC == 1) {
          for ltr(xo,cedp,&y) {
              new_printlist_nc(xo, y+co.ix[DP])
          }  
      } else {
          for ltr(xo,cedp,&y) {
              new_printlist_nc(xo.drspk, y+co.ix[DP])
          }
      }
    }
    npacsz=20
    if (recv && ce.count>0) {
      for x=0,CTYPi-1 if(co.numc[x]>0 && x!=DP) {
        max=maxrec-1
        for jj=0,max {
          XO=ce.o(co.ix[x]+jj)
          XO.recclr
          new_printlist_ac(XO,"V",  CTYP.o(x).s,XO.id)
          // new_printlist_ac(XO,"VGB",  CTYP.o(x).s,XO.id)
          // new_printlist_ac(XO,"VGA",  CTYP.o(x).s,XO.id)
          // new_printlist_ac(XO,"AHP",  CTYP.o(x).s,XO.id)
          //      printlist.o(printlist.count-1).code=3 // use code 3 for snapping
          // XO=ce.object(ixe[x]-jj-1)
          // XO.recclr
          // new_printlist_ac(XO,"V",  CTYP.object(x).s,XO.id)
        }
      }
    }
    if(lvextra!=nil){
      for(jj=0;jj<lvextra.count;jj+=2){ XO=lvextra.o(jj)
        new_printlist_ac(XO,lvextra.o(jj+1).s,CTYP.object(ctyp(XO.id)).s,XO.id)
      }
    }
  }
  if(recLFP) wrecon()
}
 
stopoq_INTF6=1

//** turnoff
proc turnoff () { local cel0,cel1,off
  cel0=$1 cel1=$2
  if (argtype(3)==0) off=$3 else off=0
  ind.indgen(ix[cel0],ixe[cel0],1) vec.indgen(ix[cel1],ixe[cel1],1) 
  intf.turnoff(ind,vec,off)
}

//** turn off intralaminar connections
proc intralamoff () { local ct,i
  for i=0,numcols-1 for col[i].ctt(&ct) if(col[i].numc[ct] && col[i].div[ct][ct]) turnoff(ct,ct)
}

//** turn on intralaminar connections
proc intralamon () { local ct,i
  for i=0,numcols-1 for col[i].ctt(&ct) if(col[i].numc[ct] && col[i].div[ct][ct]) turnoff(ct,ct,1)
}

//* mkrecl(col,ty) - setup recording of spikes from a population of cells
// return List with 2 sub-Lists: 1st with NetCons, 2nd with Vectors storing spikes
obfunc mkrecl () { local ct localobj xo,nc,col,ls,myncl,myspkl,myv
  col=$o1 ct=$2
  ls=new List()
  ls.append(myncl=new List())
  ls.append(myspkl=new List())
  for ltr(xo,col.ce) {
    if(xo.type!=ct) continue
    xo.flag("out",1) // make sure NetCon events enabled from this cell
    myncl.append(nc=new NetCon(xo,nil))
    myspkl.append(myv=new Vector())
    nc.record(myv) // record each cell separately
  }
  return ls
}

//* function calls
prl(recV,recSPKS)

Loading data, please wait...