// CA1 heteroassociative memory network: Storage and recall // CA1 PCs, BCs, AACs, BSCs and OLMs (using moderate cell models) // EC, CA3 (excitatory) and Septal (inhibitory) inputs // Cycle is: Recall-Storage-Recall etc // Parallel code adapted from Hines' ran3par.hoc // VCU & BPG 8-1-09 {load_file("nrngui.hoc")} // load the GUI and standard run libraries {load_file("pyramidal_cell4.hoc")} {load_file("basket_cell17S.hoc")} {load_file("axoaxonic_cell17S.hoc")} {load_file("bistratified_cell13S.hoc")} {load_file("olm_cell2.hoc")} {load_file("stim_cell.hoc")} {load_file("burst_cell.hoc")} {load_file("ranstream.hoc")} // to give each cell its own sequence generator STARTDEL = 50 // msecs THETA = 250 // msecs (4 Hz) GAMMA =25 // msecs (40 Hz) ECCA3DEL = 9 // msecs npcell = 100 nbcell = 2 naacell = 1 nbscell = 1 nolm = 1 nCA3 = 100 nEC = 20 nSEP = 10 ncell = npcell+nbcell+naacell+nbscell+nolm // total number of cells nstim = nCA3+nEC+nSEP // total number of inputs ntot = ncell+nstim // gid ordering: // PCs:0..npcell-1 // BCs:npcell..npcell+nbcell-1 // etc // indices of first cell of each type in list "cells" iPC = 0 iBC = npcell iAAC = npcell+nbcell iBSC = npcell+nbcell+naacell iOLM = npcell+nbcell+naacell+nbscell iCA3 = npcell+nbcell+naacell+nbscell+nolm iEC = npcell+nbcell+naacell+nbscell+nolm+nCA3 iSEP = npcell+nbcell+naacell+nbscell+nolm+nCA3+nEC ////////////////////////////////////////////////////////////// // Steps 2 and 3 are to create the cells and connect the cells ////////////////////////////////////////////////////////////// C_P = 1 // probability of excitatory connections received by each CA1 PC // from CA3 inputs (1 gives full connectivity) NPATT = 0 // number of stored patterns CFRAC = 1 // fraction of active cells in cue iNPPC=1 // index of a non-pattern PC (1st patt in 5 patterns) CPATT=1 ECPATT=1 strdef FCONN, FSTORE // file name of connection weights and patterns MOLT_THETA = 2 // multiplicator of THETA parameter SIMDUR = STARTDEL + (THETA*MOLT_THETA) // simulation duration (msecs) STORE_NEW_WEIGHTS = 1 // enables the weight matrix storage, the matrix will be the input of the Recall phase NSTORE = 1 // indicates the number of patterns to store SPATT= 20 // number of active cells per pattern SETPATT = 1 //indicates the set of patterns to be used; default is 1 (setpatt1) ORT=0 // indicates if orthogonal patterns are presented NORT=1 // /indicates the set of orthogonal patterns to be used; CREB=0 // indicates if CREB case runs ECWGT = 0.001 //EC versuc PC CLWGT = 0.00045 //cue low weight STDPPFAC = 0.6 //potentiation factor STDPDFAC = 0.8 //depression factor inputrand = 1 // seed for randomstream N_PRINT_INFO_CELLS = 7 // number of cells whose membrane potential will be printed sprint(FCONN,"Weights/wgtsN%dS%dP0.dat",npcell,SPATT) if (ORT==0) { sprint(FSTORE,"Weights/setpatt%d/pattsN%dS%dP%d.dat",SETPATT,npcell,SPATT,NSTORE) } else { sprint(FSTORE,"Weights/setpatt%d/pattsN%dS%dP%do%d.dat",SETPATT,npcell,SPATT,NSTORE,NORT) } NDIM=NSTORE // varibales for STDP window STDPDMed=-22 STDPDVar=8 STDPPTAU=10 // variables for CA1 template vinit=-65 somakap =0.007 somakad =0.007 somacaL=0.001 ghsoma = 0.00004 somacaT=0.0005 // output directory strdef running_proc,position,position2,s_sys, path, pathNORT, sfSCOPE if (ORT==0) { sprint(path, "MT_%d_NS_%d_SET_%d_INPR_%d", MOLT_THETA, NSTORE, SETPATT, inputrand) } else { sprint(path, "MT_%d_NS_%d_SET_%d_NORT_%d_INPR_%d", MOLT_THETA, NSTORE, SETPATT, NORT, inputrand) } sprint(running_proc,"bpattrun2") if (CREB==0) { sprint(position,"Results/%s/%s", running_proc, path) sprint(position2,"Results/bpattrun/%s_RECALL",path) } else { sprint(position,"Results/%s/%s_CREB", running_proc, path) sprint(position2,"Results/bpattrun/%s_RECALL_CREB", path) } sprint(s_sys,"mkdir -p %s", position) system(s_sys) sprint(s_sys,"mkdir -p %s", position2) system(s_sys) if (CREB==0) { redsAHP=1 redmAHP=1 tauca=3.6 thna3=-28 thna3dend=-28 shift=0 CHWGT = 0.000738 //cue heigh weight } else { redsAHP=0.36 redmAHP=0.48 tauca=3.6 thna3=-28 thna3dend=-28 shift=0 CHWGT = 0.001197 //cue heigh weight } // Simple connectivity CA3_PC = nCA3 // # of connections received by each PC from CA3 cells (excit) CA3_BC = nCA3 // # of connections received by each BC from CA3 cells (excit) CA3_AAC = nCA3 // # of connections received by each BC from CA3 cells (excit) CA3_BSC = nCA3 // # of connections received by each BC from CA3 cells (excit) EC_PC = nEC // # of connections received by each PC from EC cells (excit) EC_BC = nEC // # of connections received by each BC from EC cells (excit) EC_AAC = nEC // # of connections received by each AAC from EC cells (excit) SEP_BC = nSEP // # of connections received by each basket cell from septum (inhib) SEP_AAC = nSEP // # of connections received by each AAC cell from septum (inhib) SEP_BSC = nSEP // # of connections received by each BSC cell from septum (inhib) SEP_OLM = nSEP // # of connections received by each OLM cell from septum (inhib) PC_PC = 1 // # of connections received by each PC from other PCs (excit) PC_BC = npcell // # of connections received by each basket cell from PCs (excit) PC_BSC = npcell // # of connections received by each bistratified cell from PCs (excit) PC_AAC = npcell // # of connections received by each bistratified cell from PCs (excit) PC_OLM = npcell // # of connections received by each OLM cell from PCs (excit) BC_PC = 2 // # of connections received by each PC from basket cells (inhib) BC_BSC = 2 // # of connections received by each BSC from basket cells (inhib) BC_OLM = 2 // # of connections received by each OLM from basket cells (inhib) BC_BC = 1 // # of connections received by each BC from other BCs (inhib) AAC_PC = 1 // # of connections received by each PC from axoaxonic cells (inhib) BSC_PC = 1 // # of connections received by each PC from bistratified cells (inhib) BSC_BC = 1 // # of connections received by each BC from bistratified cells (inhib) OLM_PC = 1 // # of connections received by each PC from OLM cells (inhib) OLM_BC = 1 // # of connections received by each basket cell from OLM cells (inhib) Pcell2Pcell_weight = 0.001 Pcell2Pcell_delay = 1 Bcell2Pcell_weight =0.03 //0.02 Bcell2Pcell_delay = 1 Pcell2Bcell_weight =0.0005 Pcell2Bcell_delay = 1 Bcell2Bcell_weight = 0.001 Bcell2Bcell_delay = 1 Bcell2BScell_weight =0.02 Bcell2BScell_delay = 1 Bcell2OLMcell_weight = 0.0 Bcell2OLMcell_delay = 1 AAcell2Pcell_weight = 0.05 //0.04 AAcell2Pcell_delay = 1 Pcell2AAcell_weight =0.0005 Pcell2AAcell_delay = 1 BScell2Pcell_weight =0.003 //0.002 BScell2Pcell_delay = 1 BScell2Pcell_GABAB_weight=0.0004 //0.0004 BScell2Pcell_delay = 1 Pcell2BScell_weight =0.0005 Pcell2BScell_delay = 1 BScell2Bcell_weight = 0.01 BScell2Bcell_delay = 1 OLMcell2Pcell_weight =0.05//0.04 OLMcell2Pcell_GABAB_weight =0.0004 OLMcell2Pcell_delay = 1 OLMcell2Bcell_weight = 0.01 OLMcell2Bcell_delay = 1 Pcell2OLMcell_weight =0.00005 Pcell2OLMcell_delay = 1 OLMcell2Bcell_weight = 0.0 OLMcell2Bcell_delay = 1 // Synapse indices // onto CA1 PCs E_EC = 0 // EC AMPA excit to medium SLM (2 of) E_CA3 = 2 // CA3 AMPA excit to medium SR EN_CA3 = 3 // CA3 NMDA excit to medium SR EM_CA3 = 23 // CA3 modifiable (STDP) AMPA excit to medium SR E_PC = 4 // CA1 recurrent AMPA excit to proximal SR I_BC = 5 // ff&fb inhib via BCs to soma I_AAC = 6 // ff&fb inhib via AACs to axon initial segment I_BSC = 11 // ff&fb inhib via BSCs to SR med (12 of: 6 GABAA, 6 GABAB) I_OLM = 7 // fb inhib via OLMs to SLM (4 of: 2 GABAA, 2 GABAB) // onto INs (BC, AAC, BSC) EI_EC = 0 // EC AMPA excit (2 of; not onto BSC) EI_CA3 = 2 // CA3 AMPA excit (4 of) EI_PC = 6 // CA1 PC AMPA excit (2 of) II_SAME = 8 // inhib from neighbouring INs (BC->BC; BSC->BSC) II_OPP = 9 // inhib from other INs (BSC->BC; BC->BSC) II_SEP = 10 // inhib from septum (4 of: 2 GABAA, 2 GABAB) // onto INs (OLM) EO_PC = 0 // CA1 PC AMPA excit (2 of) IO_IN = 2 // inhib from INs and septum (2 of: 1 GABAA, 1 GABAB) // Septal inhibition SEPNUM = 100000 // number of SEP spikes SEPSTART = STARTDEL+(THETA/12) // time of first SEP spike SEPINT = 20 // SEP spike ISI (during burst) SEPNOISE = 0.4 // SEP ISI noise SEPBINT = 2*THETA/3 // SEP interburst interval SEPBLEN = THETA/3 // SEP burst length SEPWGT = 0.02 // SEP weight to BCs and AACs SEPWGTL = 0.0002 //0.0002 // SEP weight to BSCs and OLMs SEPWGTL2 = 0.0002 //0.0002 // SEP weight to BSCs and OLMs SEPDEL = 1 // SEP delay // Background excitation (not used) ENUM = 0 // number of spikes ESTART = 0 // time of first spike EINT = 200 // spike ISI ENOISE = 1 // ISI noise EWGT = 0.001 // excitatory weights (AMPA) ENWGT = 0.002 // excitatory weights (NMDA) EDEL = 1 // delay (msecs) // EC excitation ECNUM = 100000 // number of EC spikes ECSTART = STARTDEL // time of first EC spike ECINT = GAMMA // EC spike ISI ECNOISE = 0.2 // EC ISI noise //ECWGT = 0.0005 // 0.001// EC weight to PCs ECWGT0=0.0 ECDEL = 1 // EC delay ECIWGT =0.00015 // excitatory weights to INs CA3IWGT =0.00015 //0.00015 // excitatory weights to INs CA3IWGT2 =0.00015 //0.00015 // excitatory weights to INs EIDEL = 1 // delay (msecs) // Cue (CA3) excitation CNUM = 100000 // number of cue spikes CSTART = STARTDEL+ECCA3DEL // time of first cue spike CINT = GAMMA // cue spike ISI CNOISE = 0.2 // cue ISI noise //CHWGT = 0.0015 // cue weight //CLWGT = 0.0005 // unlearnt weight (usually 0) CNWGT = CLWGT //0.0005 // 0.0005 // excitatory weights (NMDA) CDEL = 1 // cue delay // STDP configuration //STDPPFAC = 0 // potentiation factor //STDPDFAC = 0.2 // depression factor AMPASUPP = 0.4//0.4 // fraction of AMPA during storage STDPTHRESH = -55 // voltage threshold for STDP STDPSTART = STARTDEL+(THETA/2) // STDP starts at same time as EC input STDPINT = THETA/2 // STDP interburst (recall) interval STDPLEN = THETA/2 // STDP burst (storage) length connect_random_low_start_ = 1 // low seed for mcell_ran4_init() objref cells, nclist, ncslist, ncelist, wslist, wtlist, wsEClist, wtEClist objref ranlist // for RandomStreams on this host objref stims, stimlist, cuelist, EClist // phasic and tonic cell stimulation objref gidvec // to associate gid and position in cells List // Make the network proc mknet() { mcell_ran4_init(connect_random_low_start_) nclist = new List() print "Make cells..." mkcells() // create the cells print "Make inputs..." mkinputs() // create the CA3, EC and septal inputs print "Connect cells..." //print " EC..." // EC to BC connectcells(nbcell, iBC, nEC, iEC, EC_BC, EI_EC, EI_EC+1, EIDEL, ECIWGT) // EC to AAC connectcells(naacell, iAAC, nEC, iEC, EC_AAC, EI_EC, EI_EC+1, EIDEL, ECIWGT) //print " CA3..." // CA3 to BC connectcells(nbcell, iBC, nCA3, iCA3, CA3_BC, EI_CA3, EI_CA3+3, EIDEL, CA3IWGT) // CA3 to AAC connectcells(naacell, iAAC, nCA3, iCA3, CA3_AAC, EI_CA3, EI_CA3+3, EIDEL, CA3IWGT) // CA3 to BSC connectcells(nbscell, iBSC, nCA3, iCA3, CA3_BSC, EI_CA3, EI_CA3+3, EIDEL, CA3IWGT2) //print " septum..." // SEP to BC connectcells(nbcell, iBC, nSEP, iSEP, SEP_BC, II_SEP, II_SEP+1, SEPDEL, SEPWGT) // SEP to AAC connectcells(naacell, iAAC, nSEP, iSEP, SEP_AAC, II_SEP, II_SEP+1, SEPDEL, SEPWGT) // SEP to BSC connectcells(nbscell, iBSC, nSEP, iSEP, SEP_BSC, II_SEP, II_SEP+1, SEPDEL, SEPWGTL) // SEP to OLM connectcells(nolm, iOLM, nSEP, iSEP, SEP_OLM, IO_IN, IO_IN, SEPDEL, SEPWGTL2) //print " PC..." // PC to PC connectcells(npcell, iPC, npcell, iPC, PC_PC, E_PC, E_PC, Pcell2Pcell_delay, Pcell2Pcell_weight) // PC to BC connectcells(nbcell, iBC, npcell, iPC, PC_BC, EI_PC, EI_PC+1, Pcell2Bcell_delay, Pcell2Bcell_weight) // PC to AAC connectcells(naacell, iAAC, npcell, iPC, PC_AAC, EI_PC, EI_PC+1, Pcell2AAcell_delay, Pcell2AAcell_weight) // PC to BSC connectcells(nbscell, iBSC, npcell, iPC, PC_BSC, EI_PC, EI_PC+1, Pcell2BScell_delay, Pcell2BScell_weight) // PC to OLM connectcells(nolm, iOLM, npcell, iPC, PC_OLM, EO_PC, EO_PC+1, Pcell2OLMcell_delay, Pcell2OLMcell_weight) //print " INs..." // BC to PC connectcells(npcell, iPC, nbcell, iBC, BC_PC, I_BC, I_BC, Bcell2Pcell_delay, Bcell2Pcell_weight) // BC to BC connectcells(nbcell, iBC, nbcell, iBC, BC_BC, II_SAME, II_SAME, Bcell2Bcell_delay, Bcell2Bcell_weight) // BC to BSC connectcells(nbscell, iBSC, nbcell, iBC, BC_BSC, II_OPP, II_OPP, Bcell2BScell_delay, Bcell2BScell_weight) // BC to OLM //connectcells(nolm, iOLM, nbcell, iBC, BC_OLM, IO_IN, IO_IN, Bcell2OLMcell_delay, Bcell2OLMcell_weight) // AAC to PC connectcells(npcell, iPC, naacell, iAAC, AAC_PC, I_AAC, I_AAC, AAcell2Pcell_delay, AAcell2Pcell_weight) // BSC to PC connectcells(npcell, iPC, nbscell, iBSC, BSC_PC, I_BSC, I_BSC+5, BScell2Pcell_delay, BScell2Pcell_weight) connectcells(npcell, iPC, nbscell, iBSC, BSC_PC, I_BSC+6, I_BSC+11, BScell2Pcell_delay, BScell2Pcell_GABAB_weight) // BSC to BSC //connectcells(nbscell, iBSC, nbscell, iBSC, BSC_BSC, II_SAME, II_SAME, BScell2BScell_delay, BScell2BScell_weight) // BSC to BC connectcells(nbcell, iBC, nbscell, iBSC, BSC_PC, II_OPP, II_OPP, BScell2Bcell_delay, BScell2Bcell_weight) // OLM to PC connectcells(npcell, iPC, nolm, iOLM, OLM_PC, I_OLM, I_OLM+1, OLMcell2Pcell_delay, OLMcell2Pcell_weight) connectcells(npcell, iPC, nolm, iOLM, OLM_PC, I_OLM+2, I_OLM+3, OLMcell2Pcell_delay, OLMcell2Pcell_GABAB_weight) // OLM to BC //connectcells(nbcell, iBC, nolm, iOLM, OLM_BC, II_OPP, II_OPP, OLMcell2Bcell_delay, OLMcell2Bcell_weight) //print "Connect inputs..." // EC input to PCs connectEC(E_EC,2) connectEC_2(FSTORE, $3, NDIM) // store new pattern // CA3 input to PCs connectCA3(FCONN, $2, EM_CA3, EN_CA3) // with modifiable synapses } // creates the cells and appends them to a List called cells // argument is the number of cells to be created proc mkcells() {local i,j localobj cell, nc, nil cells = new List() ranlist = new List() gidvec = new Vector() // each host gets every nhost'th cell, // starting from the id of the host // and continuing until no more cells are left for i=0, ntot-1 { if (i < iBC) { cell = new PyramidalCell(redsAHP,redmAHP,v_init,tauca,thna3,thna3dend,somakap,somakad,ghsoma,somacaT,somacaL,shift) } else if (i < iAAC) { cell = new BasketCell() // BC } else if (i < iBSC) { cell = new AACell() // AAC } else if (i < iOLM) { cell = new BistratifiedCell() // BSC } else if (i < iCA3) { cell = new OLMCell() // OLM } else if (i < iEC) { cell = new StimCell() // CA3 input } else if (i < iSEP) { cell = new StimCell() // EC input } else { cell = new BurstCell() // Septal input } cells.append(cell) ranlist.append(new RandomStream(inputrand*iSEP+i)) // ranlist.o(i) corresponds to gidvec.append(i) } } // Procedure for CA1 templeate proc caT_insert() { for (x) if (x>0 && x<1) { xdist = distance(x) insert cat if (xdist < 300) { gcatbar_cat(x) = $1*(1+xdist/60) } else { gcatbar_cat(x) = $1*6 } } } proc caL_insert() { cal_distance=200 for (x) if (x>0 && x<1) { xdist = distance(x) insert calH mytau_calH=$2 if (xdist < cal_distance) { gcalbar_calH(x) = $1 *(1-xdist/cal_distance) } else { gcalbar_calH(x) = 0 } } } proc A_h_insert(){ ghend=$1*7 dhalf=280 steep=50 KMA=3 insert h ghdbar_h=0 insert kap gkabar_kap = 0 insert kad gkabar_kad = 0 for (x) if (x>0 && x<1) { xdist=distance(x) ghdbar_h(x)= $1 + (ghend - $1)/(1.0 + exp((dhalf-xdist)/steep)) if (xdist < 100){ gkabar_kap(x) = $2*(1+KMA*xdist/100) vhalfl_h=-73 }else{ vhalfl_h=-81 gkabar_kad(x) = $3*(1+KMA*xdist/100) // print secname(), " ",xdist, gkabar_kad(x) } } } proc acc_dist(){local i for i=0, npcell-1 { access cells.object(i).soma distance(0,x) access cells.object(i).radTprox A_h_insert(ghsoma,somakap,somakad) caT_insert(somacaT) caL_insert(somacaL,tauca) access cells.object(i).radTmed A_h_insert(ghsoma,somakap,somakad) caT_insert(somacaT) caL_insert(somacaL,tauca) access cells.object(i).radTdist A_h_insert(ghsoma,somakap,somakad) caT_insert(somacaT) caL_insert(somacaL,tauca) access cells.object(i).lm_thick1 A_h_insert(ghsoma,somakap,somakad) access cells.object(i).lm_medium1 A_h_insert(ghsoma,somakap,somakad) access cells.object(i).lm_thin1 A_h_insert(ghsoma,somakap,somakad) access cells.object(i).lm_thick2 A_h_insert(ghsoma,somakap,somakad) access cells.object(i).lm_medium2 A_h_insert(ghsoma,somakap,somakad) access cells.object(i).lm_thin2 A_h_insert(ghsoma,somakap,somakad) cells.object(i).current_balance(v_init) } } // sets the CA3, EC and Septal background inputs proc mkinputs() {local i localobj stim, rs for i=0, cells.count-1 { gid = gidvec.x[i] // id of cell if (gid >= iCA3 && gid < ntot-nSEP-nEC) { // appropriate target cell // set background activity for excitatory inputs stim = cells.object(i).stim stim.number = ENUM stim.start = ESTART stim.interval = EINT stim.noise = ENOISE } if (gid >= iEC && gid < ntot-nSEP) { // appropriate target cell // set background activity for excitatory inputs stim = cells.object(i).stim stim.number = ENUM stim.start = ESTART stim.interval = EINT stim.noise = ENOISE } if (gid >= iSEP && gid < ntot) { // appropriate target cell // set background activity for septum stim = cells.object(i).stim rs = ranlist.object(i) stim.number = SEPNUM stim.start = SEPSTART stim.interval = SEPINT stim.noise = SEPNOISE stim.burstint = SEPBINT stim.burstlen = SEPBLEN // Use the gid-specific random generator so random streams are // independent of where and how many stims there are. stim.noiseFromRandom(rs.r) rs.r.negexp(1) rs.start() } } } // Target cells receive "convergence" number of inputs from // the pool of source cells (only one input per source cell at most) // ("convergence" not reached if no. of sources < convergence) // connectcells(number of targets, first target cell, // number of source cells, first source cell, // convergence, first synapse, // last synapse, connection delay, weight) // appends the NetCons to a List called nclist proc connectcells() {local i, j, gid, nsyn, r localobj syn, nc, rs, u // initialize the pseudorandom number generator u = new Vector($3) // for sampling without replacement // print "NUMERO ELEMENTI: ", nclist.count, " per il processore ", pc.id for i=0, cells.count-1 { // loop over possible target cells gid = gidvec.x[i] // id of cell if (gid >= $2 && gid < $1+$2) { // appropriate target cell rs = ranlist.object(i) // RandomStream for cells.object(i) rs.start() rs.r.discunif($4, $4+$3-1) // return source cell index u.fill(0) // u.x[i]==1 means spike source i has already been chosen nsyn = 0 while (nsyn < $5 && nsyn < $3) { r = rs.repick() // no self-connection and only one connection from any source if (r != gidvec.x[i]) if (u.x[r-$4] == 0) { // target synapses for j = $6, $7 { // set up connection from source to target syn = cells.object(i).pre_list.object(j) // nc = pc.gid_connect(r, syn) nc = cells.object(r).connect2target(syn) nclist.append(nc) nc.delay = $8 nc.weight = $9 } u.x[r-$4] = 1 nsyn += 1 } } } } } // connects the EC input layer to PC cells // appends the PC NetCons to a List called ncslist proc connectEC() {local i, j, k, gid localobj syn, src, nc, fp, target ncelist = new List() wsEClist = new Vector() wtEClist = new Vector() for i=0, cells.count-1 { // loop over possible target cells gid = gidvec.x[i] // id of cell if (gid >= iPC && gid < npcell+iPC) { // appropriate target cell print "cell ", gid target = cells.object(i+iPC) // target all synapses for k = $1, $1+$2-1 { syn = target.pre_list.object(k) // excitatory synapse // create pattern stimulus for j=0, nEC-1 { src = cells.object(j+iEC).stim nc = new NetCon(src, syn) nc.delay = ECDEL nc.weight = ECWGT0 ncelist.append(nc) wsEClist.append(j+iEC) wtEClist.append(gid) } } } } } proc connectEC_2() {local i, j,gid, ncue, index localobj cue, cstim, syn, src, fp, target fp = new File($s1) fp.ropen() cue = new Vector(npcell) cue.scanf(fp, $2, $3) // read pattern fp.close() ncue = 0 // find active cells in pattern for i=0, cue.size()-1 { for j=0, 2*nEC-1 { ncelist.object(i*2*SPATT+j).weight = ECWGT0 if (cue.x[i] == 1) { ncelist.object(i*2*SPATT+j).weight = ECWGT } } } } // connects the CA3 input layer to output cells // appends the PC NetCons to a List called ncslist proc connectCA3() {local i, j, cp, gid localobj src, syn, synN, nc, fc, rs, conns, rc cp = $2 // connection probability mcell_ran4_init(connect_random_low_start_) conns = new Vector(nCA3) // connection weights rc = new Vector(nCA3) // random physical connectivity ncslist = new List() wslist = new Vector() wtlist = new Vector() for i=0, cells.count-1 { // loop over possible target cells gid = gidvec.x[i] // id of cell if (gid >= iPC && gid < npcell+iPC) { // appropriate target cell print "cell ", gid syn = cells.object(i).pre_list.object($3) // AMPA synapse with STDP syn.wmax = CHWGT syn.wmin = CLWGT syn.d = STDPDFAC // depression syn.p = STDPPFAC // potentiation syn.ptau = STDPPTAU // potentiation time syn.gscale = AMPASUPP // fraction of AMPA during storage syn.thresh = STDPTHRESH // threshold for postsynaptic voltage detection syn.gbdel = STDPSTART syn.gbint = STDPINT syn.gblen = STDPLEN synN = cells.object(i).pre_list.object($4) // NMDA synapse rs = ranlist.object(i) // the corresponding RandomStream rs.start() rs.r.uniform(0, 1) // return integer in range 0..1 rc.setrand(rs.r) // generate random connectivity // open connections file fc = new File($s1) fc.ropen() conns.scanf(fc, gid+1, nCA3) // read incoming weights for cell gid fc.close() for j=0, nCA3-1 { // only connection if physical connection exists if (rc.x[j] <= cp) { src = cells.object(j+iCA3).stim nc = new NetCon(src, synN) ncslist.append(nc) wslist.append(j+iCA3) wtlist.append(gid) nc.delay = CDEL nc.weight = CNWGT // NMDA weight same for all connections // set up connection from source to AMPA target nc = new NetCon(src, syn) ncslist.append(nc) wslist.append(j+iCA3) wtlist.append(gid) nc.delay = CDEL nc.weight = CLWGT // unlearned weight } } } } } proc connectCA3_2() {local i, j, k, gid, index localobj src, syn, nc, fc, conns conns = new Vector(nCA3) // connection weights // inputs to PCs determined by weight matrix for i=0, cells.count-1 { // loop over possible target cells gid = gidvec.x[i] // id of cell if (gid >= iPC && gid < npcell+iPC) { // appropriate target cell // open connections file fc = new File($s1) fc.ropen() conns.scanf(fc, gid+1, nCA3) // read incoming weights for cell gid fc.close() k=0 for(j=1; j<2*nCA3; j=j+2){ // set up connection from source to target ncslist.object(i*2*nCA3+j).weight[0]=conns.x[k] k=k+1 } // for j } // if gid } // for i } strdef s, s_system strdef running_proc objref m, v1, v2 proc cancstim() {local i, j, old_cstim, new_cstim localobj cstim for i=0, cells.count-1 { gid = gidvec.x[i] // id of cell if (gid >= iCA3 && gid < nCA3+iCA3) { // appropriate target cell cstim = cells.object(i).stim old_cstim=cstim.number cstim.number = 0 new_cstim=cells.object(i).stim.number } if (gid >= iEC && gid < nEC+iEC) { // appropriate target cell cstim = cells.object(i).stim old_cstim=cstim.number cstim.number = 0 new_cstim=cells.object(i).stim.number } } } proc aweights() {local i, j, simul, count localobj matrix, f, f3, f_T simul=$1 f3 = new File() sprint(s,"%s/AMPA_WEIGHT.dat",$s2) f3.wopen(s) //Print AMPA synapses for(j=1; j= iEC && gid < iEC+nEC) { // appropriate target cell // create cue stimulus cstim = cells.object(i).stim rs = ranlist.object(i) cstim.number = ECNUM cstim.start = ECSTART cstim.interval = ECINT cstim.noise = ECNOISE // Use the gid-specific random generator so random streams are // independent of where and how many stims there are. cstim.noiseFromRandom(rs.r) rs.r.normal(0, 1) rs.start() EClist.append(i) necs += 1 } } } objref cue, fp // setup activity pattern in input cue stims proc mkcue() {local i, j, ncue localobj cstim, target, rs cuelist = new Vector() // open patterns file fp = new File($s1) fp.ropen() cue = new Vector(nCA3) cue.scanf(fp, $2, $4) // read pattern // cue.printf() fp.close() ncue = 0 // find active cells in pattern for i=0, cue.size()-1 { if (ncue <= SPATT*$3) { // fraction of active cells in cue if (cue.x[i] == 1) { print "Cue cell ", i cstim = cells.object(i+iCA3).stim for j=0, cells.count-1 { if (gidvec.x[j] == i+iCA3) {break} // find cell index } rs = ranlist.object(j) // create cue stimulus //cstim = target.stim cstim.number = CNUM cstim.start = CSTART cstim.interval = CINT cstim.noise = CNOISE // Use the gid-specific random generator so random streams are // independent of where and how many stims there are. cstim.noiseFromRandom(rs.r) rs.r.normal(0, 1) rs.start() cuelist.append(i) ncue += 1 } } } // print "cue size ", ncue } // remove activity pattern in input cue stims proc erasecue() {local i, j localobj cstim for i=0, cue.size()-1 { cstim = cells.object(cuelist.x[i]+iCA3).stim cstim.number = 0 } } //mkcue(FSTORE, CPATT, CFRAC, NSTORE) // cue from new pattern //mkEC() // Spike recording objref tvec, idvec // will be Vectors that record all spike times (tvec) // and the corresponding id numbers of the cells that spiked (idvec) proc spikerecord() {local i localobj nc, nil tvec = new Vector() idvec = new Vector() for i=0, cells.count-1 { nc = cells.object(i).connect2target(nil) //nc.record(tvec, idvec, nc.srcgid) nc.record(tvec, idvec, i) // the Vector will continue to record spike times // even after the NetCon has been destroyed } } spikerecord() // Record cell voltage traces objref pvsoma, pvsr, pvslm // Vectors that record voltages from pattern PC objref npvsoma, npvsr, npvslm // Vectors that record voltages from non-pattern PC objref vBC, vAAC, vBSC, vOLM // Vectors that record voltages from INs objref list_pvsoma, list_pvsr, list_pvslm proc vrecord2() {local i, gid list_pvsoma = new List() list_pvsr = new List() list_pvslm = new List() for i=0, cells.count-1 { // loop over possible target cells gid = gidvec.x[i] // id of cell if (gid< N_PRINT_INFO_CELLS) { pvsoma = new Vector() pvsr = new Vector() pvslm = new Vector() pvsoma.record(&cells.object(i).soma.v(0.5)) pvsr.record(&cells.object(i).radTmed.v(0.5)) pvslm.record(&cells.object(i).lm_thick1.v(0.5)) list_pvsoma.append(pvsoma) list_pvsr.append(pvsr) list_pvslm.append(pvslm) } if (gid==iNPPC) { npvsoma = new Vector() npvsr = new Vector() npvslm = new Vector() npvsoma.record(&cells.object(i).soma.v(0.5)) npvsr.record(&cells.object(i).radTmed.v(0.5)) npvslm.record(&cells.object(i).lm_thick1.v(0.5)) } if (gid==iBC) { vBC = new Vector() vBC.record(&cells.object(i).soma.v(0.5)) } if (gid==iAAC) { vAAC = new Vector() vAAC.record(&cells.object(i).soma.v(0.5)) } if (gid==iBSC) { vBSC = new Vector() vBSC.record(&cells.object(i).soma.v(0.5)) } if (gid==iOLM) { vOLM = new Vector() vOLM.record(&cells.object(i).soma.v(0.5)) } } } vrecord2() //////////////////////////// // Simulation control //////////////////////////// strdef fstem, s_sys,position,position2 tstop = SIMDUR celsius = 34 proc bpattrun2() { local i for i=1, NSTORE { if (CREB==0) { sprint(position,"Results/bpattrun2/%s", path) sprint(position2,"Results/bpattrun/%s_RECALL",path) } else { sprint(position,"Results/bpattrun2/%s_CREB",path) sprint(position2,"Results/bpattrun/%s_RECALL_CREB", path) } print "Cue pattern ", i sprint(fstem, "%s/HAM_P%dR%d",position,NPATT,i) if(i!=1) { sprint(FCONN,"%s/wgtsN%dS%dP%db_bpattrun2.dat" ,position,npcell,SPATT, i-1) // EC input to PCs connectEC_2(FSTORE, i, NDIM) // store new pattern // CA3 input to PCs connectCA3_2(FCONN) // with modifiable synapses } mkEC() {mkcue(FSTORE, i, CFRAC, NDIM)} // cue from stored pattern stdinit() run() printf("pos %d %s",i,position) aweights(i,position,position2) spikeout() // vout2() //erasecue() cancstim() } } //////////////////////////// // Report simulation results //////////////////////////// objref fo strdef fno proc spikeout() { local i printf("\ntime\t cell\n") // print header once sprint(fno,"%s_spt.dat", fstem) fo = new File(fno) fo.wopen() for i=0, tvec.size-1 { printf("%g\t %d\n", tvec.x[i], idvec.x[i]) fo.printf("%g\t %d\n", tvec.x[i], idvec.x[i]) } fo.close() } proc vout2() { local i, j, gid for j=0, cells.count-1 { // loop over possible target cells gid = gidvec.x[j] // id of cell if (gid< N_PRINT_INFO_CELLS) { sprint(fno,"%s_pvsoma_%d.dat", fstem,gid) fo = new File(fno) fo.wopen() for i=0, pvsoma.size-1 { fo.printf("%g\n", list_pvsoma.o(j).x[i]) } fo.close() sprint(fno,"%s_pvsr_%d.dat", fstem,gid) fo = new File(fno) fo.wopen() for i=0, pvsr.size-1 { fo.printf("%g\n", list_pvsr.o(j).x[i]) } fo.close() sprint(fno,"%s_pvslm_%d.dat", fstem,gid) fo = new File(fno) fo.wopen() for i=0, pvslm.size-1 { fo.printf("%g\n", list_pvslm.o(j).x[i]) } fo.close() } if (gid==iNPPC) { sprint(fno,"%s_npvsoma.dat", fstem) fo = new File(fno) fo.wopen() for i=0, npvsoma.size-1 { fo.printf("%g\n", npvsoma.x[i]) } fo.close() sprint(fno,"%s_npvsr.dat", fstem) fo = new File(fno) fo.wopen() for i=0, npvsr.size-1 { fo.printf("%g\n", npvsr.x[i]) } fo.close() sprint(fno,"%s_npvslm.dat", fstem) fo = new File(fno) fo.wopen() for i=0, npvslm.size-1 { fo.printf("%g\n", npvslm.x[i]) } fo.close() } if (gid==iBC) { sprint(fno,"%s_BC.dat", fstem) fo = new File(fno) fo.wopen() for i=0, vBC.size-1 { fo.printf("%g\n", vBC.x[i]) } fo.close() } if (gid==iAAC) { sprint(fno,"%s_AAC.dat", fstem) fo = new File(fno) fo.wopen() for i=0, vAAC.size-1 { fo.printf("%g\n", vAAC.x[i]) } fo.close() } if (gid==iBSC) { sprint(fno,"%s_BSC.dat", fstem) fo = new File(fno) fo.wopen() for i=0, vBSC.size-1 { fo.printf("%g\n", vBSC.x[i]) } fo.close() } if (gid==iOLM) { sprint(fno,"%s_OLM.dat", fstem) fo = new File(fno) fo.wopen() for i=0, vOLM.size-1 { fo.printf("%g\n", vOLM.x[i]) } fo.close() } } } objref gs proc spikeplot() { local i gs = new Graph() gs.size(0, tstop, -1, ntot) for i=0, tvec.size-1 { gs.mark(tvec.x[i], idvec.x[i], "|", 8) } gs.flush() } // panel for simulation results proc xspikeres() { xpanel("Spike results") xbutton("Plot", "spikeplot()") xpanel() } xspikeres() xopen("HAM_SR1.ses") bpattrun2() //quit()