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

 Download zip file   Auto-launch 
Help downloading and running models
Accession:146949
"... We developed a computational model based primarily on a unified set of brain activity mapping studies of mouse M1. The simulation consisted of 775 spiking neurons of 10 cell types with detailed population-to-population connectivity. Static analysis of connectivity with graph-theoretic tools revealed that the corticostriatal population showed strong centrality, suggesting that would provide a network hub. ... By demonstrating the effectiveness of combined static and dynamic analysis, our results show how static brain maps can be related to the results of brain activity mapping."
Reference:
1 . 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 [PubMed]
Citations  Citation Browser
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 L6 pyramidal corticothalamic cell; Neocortex M1 L2/6 pyramidal intratelencephalic cell; Neocortex fast spiking (FS) 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;
Model Concept(s): Oscillations; Laminar Connectivity;
Implementer(s): Lytton, William [billl at neurosim.downstate.edu]; Neymotin, Sam [samn at neurosim.downstate.edu]; Shepherd, Gordon MG [g-shepherd at northwestern.edu]; Chadderdon, George [gchadder3 at gmail.com]; Kerr, Cliff [cliffk at neurosim.downstate.edu];
Search NeuronDB for information about:  Neocortex V1 L6 pyramidal corticothalamic cell; Neocortex M1 L2/6 pyramidal intratelencephalic cell; GabaA; AMPA; NMDA; Gaba; Glutamate;
/
src
README
infot.mod *
intf6.mod *
intfsw.mod *
matrix.mod
misc.mod *
nstim.mod *
staley.mod *
stats.mod *
vecst.mod *
boxes.hoc *
col.hoc
declist.hoc *
decmat.hoc *
decnqs.hoc *
decvec.hoc *
default.hoc *
drline.hoc *
filtutils.hoc *
gcelldata.hoc
gmgs102.nqs
grvec.hoc *
infot.hoc *
init.hoc
intfsw.hoc *
labels.hoc *
load.py
local.hoc *
main.hoc
misc.h *
miscfuncs.py
network.hoc
neuroplot.py *
nload.hoc
nqs.hoc *
nqsnet.hoc
nrnoc.hoc *
params.hoc
run.hoc
samutils.hoc *
saveoutput.hoc
saveweights.hoc
setup.hoc *
simctrl.hoc *
spkts.hoc *
staley.hoc *
stats.hoc *
stdgui.hoc *
syncode.hoc *
updown.hoc *
wdmaps2.nqs
xgetargs.hoc *
                            
// $Id: run.hoc,v 1.53 2010/09/19 19:29:07 samn Exp $

print "Loading run.hoc..."

objref wf1,wf2,wrec
{sprint(tstr,"o[%d]",numcols)}
{declare("vit",tstr,"nqLFP",tstr,"SPKSZ",25e6*tstop/20e3)}

{method("CN",1)} // CK: use Crank-Nicholson integration; 2nd argument specifies timestep; 1 ms seems fine surprisingly
{cvode.atol(1e-3)}
{cvode.condition_order(1)} // irrelevant to acells?
{declare("sepflds",1)} // whether to record separate fields for each layer @@@ CK: Yes, please do!

//** draw lines between cell subpopulations
proc rasterlines () { localobj o
  g=graphItem
  for ct=0,CTYPi-1 if(col.numc[ct]) {
    drline(0,col.ixe[ct],mytstop,col.ixe[ct],g,2,6)
    // o=mdl2view(g,0,ix[ii]+numc[ii]/10)
    o=mdl2view(g,0.9,col.ix[ct]+col.numc[ct]/4)
    g.label(0.9,o.x[3],CTYP.o(ct).s)
  }
  g.flush
}
//** calls rasterlines, making sure g is set first
proc grlines () { 
  {g=Graph[$1] rasterlines()}
}

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()
}


//** init - called @ start of run
proc init () { 
  initMisc1()
  vseed_stats(392426)
  finitialize() 
  cvode.re_init()
}

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

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

//** initMisc1 -- @@@ CK: replaced with version Sam sent me on 2011feb11
proc initMisc1 () {
  col.intf.global_init()
  jrtm_INTF6=100 // Time step between updates
  for i=0,numcols-1 col[i].cstim.initrands()
}

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,x,ct,cnt localobj tl,vm
  tl=new List()
  for cdx=0,numcols-1 {
    {nqsdel(nqLFP[cdx])}
    if(sepflds) nqLFP[cdx]=new NQS("E2","E4","E5","E6","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(E2)=0 vm.x(E4)=1 vm.x(E5R)=vm.x(E5B)=2 vm.x(E6)=3}
  col.intf.initwrec(tl)
  for x=0,numcols-1 {
    for case(&ct,E2,E4,E5R,E5B,E6) {
      for ii=col[x].ix[ct],col[x].ixe[ct] {
        if(sepflds) {
          col[x].ce.o(ii).wrc(vm.x(ct)+x*5)
          col[x].ce.o(ii).wrc(4+x*5)
        } 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()
objref p
p=printlist
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)        }
        }
      }
    }
    npacsz=20
    if (recv && ce.count>0) {
      if(numcols <= 1 || (co.xpos==int(colW/2) && co.ypos==int(colH/2))) for x=0,CTYPi-1 {      
        if (co.numc[x]>0) max=0 else max=co.numc[x]-1
        for jj=0,max {
          XO=ce.object(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)
      }
    }
  }
  wrecon()
}

obfunc md5 () { localobj aq
  batch_flag=1
  if (numarg()>0) tstr=$s1 else tstr="aa"
  aq=new NQS("t","i")   aq.scpflag=1 aq.setcols(p.o(0).tvec,p.o(0).vec)
  aq.sort("i") aq.sort("t") aq.listvecs()
  prveclist(tstr,aq.vl)
  sprint(tstr,"md5sum %s",tstr) system(tstr)
  batch_flag=0
  return aq
}

//** rub() -- multi-run with saving to veclist
// saves vspks(1), SM(2), SU(2), IN(2)
proc rub () { local ii localobj so
  so=new String()
  clrveclist()
  for ii=0,9 {
    sprint(so.s,"%d",ii)
    shuffle(vspks)
    savevec(vspks)
    time()
    savevec(p.object(0).vec,p.object(0).tvec)
    savevec(p.object(1).vec,p.object(1).tvec)
    savevec(p.object(2).vec,p.object(2).tvec)
  }
}

//** spri(#) puts up graph of just one class
// 0->SM 1->SU 2->IN
// map 0,1,2 onto 1,3,5 ii*2+1
proc spri () { local ii,a,jj
  a=allocvecs(1)  ge(0)
  for ii=0,9 { jj=7*ii+$1*2+1
    mso[a].copy(veclist.object(jj))
    mso[a].add(1.5*ii*(mso[a].max-mso[a].min)-mso[a].min)
    mso[a].mark(g,veclist.object(jj+1),"O",4,cg(ii))
  }
}
//** sprj(#) puts up SM,SU,IN for 1 run (cp a() above)
proc sprj () { local ii,jj
  jj=7*$1+1
  for ii=0,2 veclist.object(jj+ii*2).mark(g,veclist.object(jj+ii*2+1),"O",2,cg(ii))
}
 
flddur=celdur=1000
splshhsz=0.4
objref slicepictypes
slicepictypes=new Vector()
{slicepictypes.append(0,3,2)}
stopoq_INTF6=1

proc rer () { 
  intf.resetall
  prl()
  shock()
  srun()
}

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
  for ctt(&ct) if(numc[ct] && div[ct][ct][0]) turnoff(ct,ct)
}

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

load_file("spkts.hoc")
objref snq,fnq,anq,cvnq
//get CVPNQS
proc getcvnq () {
  if(cvnq!=nil)nqsdel(cvnq)
  if(numarg()==0) cvnq=CVPNQS(snq,1,1,0)
  if(numarg()==1) cvnq=CVPNQS(snq,$1,1,0)
  if(numarg()==2) cvnq=CVPNQS(snq,$1,$2,0)
  if(numarg()==3) cvnq=CVPNQS(snq,$1,$2,$3)
}
//get PActNQS
proc getanq () {
  if(anq!=nil)nqsdel(anq)
  if(numarg()==0) {
    anq=PActNQS(snq)
  } else if(numarg()==1) {
    anq=PActNQS(snq,$1) 
  } else if(numarg()==2) {
    anq=PActNQS(snq,$1,$2) 
  }
}
//get SpikeNQS
proc getsnq () {
  if(snq!=nil)nqsdel(snq)
  snq=SpikeNQS(printlist.o(0),0)
}
//get SpikeNQS,FreqNQS
proc getsfnq () {
  getsnq()
  if(fnq!=nil)nqsdel(fnq)
  fnq=FreqNQS(snq,20,1,0)
}
//display & print average vals of snq,fnq
proc dispfnq () { local ct
  fnq.verbose=0
  for ctt(&ct) {
    fnq.select("Type",ct,"StartT","<=",tmax_INTF6)
    printf("%s mean F = %g Hz\n",CTYP.o(ct).s,fnq.getcol("Freq").mean)
    fnq.gr("Freq","EndT",0,clr,1)
    clr+=1
  }
}
//draw a black box showing stim duration
proc dispsgrdur () {
  if(g==nil)gg()
  drline(0,0,0,$1,g,1,4)
  drline(sgrdur,0,sgrdur,$1,g,1,4)
  drline(0,$1,sgrdur,$1,g,1,4)
  drline(0,0,sgrdur,0,g,1,4)
}
//draw $2 cell type in $3 color from PActNQS $o1
proc dispPActNQ () { local ct,color localobj anq,str  
  anq=$o1 ct=$2 color=$3 str=new String() str.s="act"
  if(numarg()>3)str.s=$s4
  anq.select("ct",ct)
  gvmarkflag=0
  anq.gr(str.s,"ts",0,color,9)
  gvmarkflag=1
  anq.gr(str.s,"ts",0,color,10)
  gvmarkflag=0
}
//get nqs with prumat info
obfunc getprnq () { local prty,poty localobj prnq
  prnq=new NQS("prty","poty","pr")
  for ctt(&prty) for ctt(&poty) if(prumat[prty][poty]!=0) prnq.append(prty,poty,prumat[prty][poty])
  return prnq
}
//get nqs with sprmat info
obfunc getsprnq () { local prty,poty localobj sprnq
  sprnq=new NQS("prty","poty","spr")
  for ctt(&prty) for ctt(&poty) if(sprmat[prty][poty]!=0) sprnq.append(prty,poty,sprmat[prty][poty])
  return sprnq
}

//rerun sim from nqs -- use :
// objref nqn
//nqn=new NQS("/u/samn/vcortex/data/08nov7_dodviptbatch_NORMUNIFNEGEXPDVFIXED_S0_399_eg2p15_update.nqs")
//rernqs(nqn,rowid)
proc rernqs () { local s,tm,ix,w localobj iq,wmnq
  iq=$o1 ix=$2 // ix will be row num in the NQS
  w=iq.fi("wght")
  s=iq.get("s",ix).x setdviPT=iq.get("setdviPT",ix).x tm=iq.get("tmax",ix).x 
  EGain=iq.get("EGain",ix).x IGain=iq.get("IGain",ix).x sgrdur=iq.get("sgrdur",ix).x
  if(iq.fi("vsgrpp")!=-1)vsgrpp=iq.get("vsgrpp",ix).o
  if(iq.fi("vsgrsidx")!=-1)vsgrsidx=iq.get("vsgrsidx",ix).o
  if(iq.fi("wght")!=-1){
    wmnq=iq.get("wght",ix).o setwt(wmnq)
  } else {
    setwt()
  }
  setdvi(s)
  stim()
  print "\nSim should run out to tmax=",tm
  time()
}

//* function calls
prl(1,1)