# -*- coding: utf-8 -*-
"""
Created on Tue May 10 11:07:00 2016
@author: Conor D. cox
@author: Conor D. Cox
Notes:
this requires brian2 and python 2.7, though it might run on python 3 the print functions are all wrong.
WARNING WARNING WARNING
This requires a number of GB of ram equal to the number of seconds you are simulating. If your computer doesn't have that
it will crash this can be ameliorated by changing the recording function down around 436ish. You can turn
those off and it will require less RAM.
The output of this, other than dozens of figures is the output array. You can save it as you see fit
This code representes figure 7
email me for any problems:
cdcox1@gmail.com
"""
# -*- coding: utf-8 -*-
from __future__ import division
import numpy as np
from brian2 import *
import copy
import gc as gc11
from sklearn.cluster import DBSCAN
from sklearn import metrics
from sklearn.datasets.samples_generator import make_blobs
from sklearn.preprocessing import StandardScaler
import matplotlib.pyplot as plt
def cluster_algorithim(nnzoom,i,innnnn):
'''this takes the signal from a region and performs DBSCAN clustering on it
it requires you to scale your time input to a reasonable spatial input,
in this case time was scaled by a factor of 10. Eps represents approximate
cluster size.'''
X=np.copy(nnzoom)
X[:,0]=X[:,0]/10
db = DBSCAN(eps=10, min_samples=5).fit(X)
core_samples_mask = np.zeros_like(db.labels_, dtype=bool)
core_samples_mask[db.core_sample_indices_] = True
labels = db.labels_
n_clusters_ = len(set(labels)) - (1 if -1 in labels else 0)
# Black removed and is used for noise instead.
unique_labels = set(labels)
colors = plt.cm.prism(np.linspace(0, 1, len(unique_labels)))
for k, col in zip(unique_labels, colors):
if k == -1:
# Black used for noise.
col = 'k'
class_member_mask = (labels == k)
xy = X[class_member_mask & core_samples_mask]
plt.plot(xy[:, 0], xy[:, 1], 'o', markerfacecolor=col,
markeredgecolor='k', markersize=14)
xy = X[class_member_mask & ~core_samples_mask]
plt.plot(xy[:, 0], xy[:, 1], 'o', markerfacecolor=col,
markeredgecolor='k', markersize=6)
plt.plot(X[:,0],X[:,1],'+k')
plt.title('Estimated number of clusters: %d' % n_clusters_)
plt.show()
plt.savefig('nnncca3clustconnect'+str(i)+str(innnnn)+'.png')
med_clus_size=median(np.histogram(labels,bins=len(unique_labels)+1,range=(-1,len(unique_labels)))[0][1:-1])
num_clust=len(unique_labels)
return med_clus_size,num_clust
def make_targetting_array():
targetting_count_array=[]
j=1
for i in range(1100):
if i%11!=0:
targetting_count_array.append(j)
j+=1
else:
targetting_count_array.append(0)
targetting_count_array=np.array(targetting_count_array)
return targetting_count_array
def randomize_and_edge_check(real_i,i,targetting_count_array,number_of_neurons):
check=np.random.normal(scale=1000)
target=np.round(check/10)
temp_tar=real_i+target
if ((i+target)<0) or ((temp_tar)>=len(targetting_count_array)):
checker=999999999
return checker,check,12345678
checker=targetting_count_array[temp_tar]
if checker>number_of_neurons:
checker=999999999
return checker,check,12345678
return checker,check,temp_tar
def randomize_and_edge_checkint(real_i,i,targetting_count_array):
check=np.random.normal(scale=100)
target=np.round(check/10)
temp_tar=real_i+target
if ((real_i+target)<0) or ((temp_tar)>=len(targetting_count_array)):
checker=999999999
return checker,check,12345678
checker=targetting_count_array[temp_tar]
if checker>number_of_ca1_neurons:
checker=999999999
return checker,check,12345678
return checker,check,temp_tar
def excite_to_excite_wire_it_up(number_of_neurons,exkij,exkijsd,inkj,inkijsd):
connectivity=[]
connectivityin=[]
delayex=[]
delayin=[]
targetting_count_array=make_targetting_array()
for i in range(1,number_of_neurons+1):
in_count=0
numb_in_connect=np.int(np.round(np.random.normal(loc=inkj,scale=inkijsd)))
if exkij!=0:
numb_connect=np.int(np.round(np.random.normal(loc=exkij,scale=exkijsd)))
else:
numb_connect=0
real_i=np.where(targetting_count_array==i)[0][0]
for j in range(1,numb_connect):
checker,check,temp_tar=randomize_and_edge_check(real_i,i,targetting_count_array,number_of_neurons)
if checker==999999999:
continue
while temp_tar%11==0 or checker==i:
if in_count<numb_in_connect and checker!=i and temp_tar!=0:
connectivityin.append([i,np.int(temp_tar/11)])
delayin.append(np.abs(np.round(check*um/speed,decimals=4)))
in_count+=1
checker,check,temp_tar=randomize_and_edge_check(real_i,i,targetting_count_array,number_of_neurons)
if checker==999999999:
continue
if checker==999999999:
continue
connectivity.append([i,checker])
delayex.append(np.abs(np.round(check*um/speed,decimals=4)))
while in_count<numb_in_connect:
checker,check,temp_tar=randomize_and_edge_check(real_i,i,targetting_count_array,number_of_neurons)
if checker==999999999:
continue
if temp_tar%11==0 and temp_tar!=0:
connectivityin.append([i,np.int(temp_tar/11)])
delayin.append(np.abs(np.round(check*um/speed,decimals=4)))
in_count+=1
connectivity=np.array(connectivity)
connectivityin=np.array(connectivityin)
return connectivity,connectivityin,delayex,delayin
def ca1in_wire_it_up(number_of_ca1_neurons,exkij,exkijsd,inkj,inkijsd):
speed2=100*um/ms
connectivity=[]
connectivityin=[]
delayex=[]
delayin=[]
targetting_count_array=make_targetting_array()
for i in range(1,int(number_of_ca1_neurons/10)+1):
in_count=0
numb_in_in_connect=np.int(np.round(np.random.normal(loc=inkj,scale=inkijsd)))
numb_2_ex_connect=np.int(np.round(np.random.normal(loc=exkij,scale=exkijsd)))
real_i=np.where(targetting_count_array==i)[0][0]
for j in range(1,numb_2_ex_connect):
checker,check,temp_tar=randomize_and_edge_checkint(i*11,i,targetting_count_array)
if checker==999999999:
continue
while temp_tar%11==0 or checker==i:
if in_count<numb_in_in_connect and checker!=i and temp_tar!=0:
connectivityin.append([i,np.int(temp_tar/11)])
delayin.append(np.abs(np.round(check*um/speed2,decimals=4)))
in_count+=1
checker,check,temp_tar=randomize_and_edge_checkint(i*11,i,targetting_count_array)
if checker==999999999:
continue
if checker==999999999:
continue
connectivity.append([i,checker])
delayex.append(np.abs(np.round(check*um/speed2,decimals=4)))
while in_count<numb_in_in_connect:
checker,check,temp_tar=randomize_and_edge_checkint(i*11,i,targetting_count_array)
if checker==999999999:
continue
if temp_tar%11==0 and temp_tar!=0:
connectivityin.append([i,np.int(temp_tar/11)])
delayin.append(np.abs(np.round(check*um/speed2,decimals=4)))
in_count+=1
connectivity=np.array(connectivity)
connectivityin=np.array(connectivityin)
return connectivity,connectivityin,delayex,delayin
def in2ex_wire_it_up(number_of_ca3_neurons,ca3_inpykij,ca3_inpykij_sd):
#wires up CA3 in2ex
targetting_count_array=make_targetting_array()
connectivity=[]
delay=[]
for i in range(1,101):
numb_connect=np.int(np.round(np.random.normal(loc=ca3_inpykij,scale=ca3_inpykij_sd)))
cell_pos=i*11
for j in range(1,numb_connect):
step=np.random.uniform(-1.5*100,1.5*100)
cell_step=step/10
target_cell=np.int(np.round(cell_pos+cell_step))
if target_cell<0 or (target_cell>=len(targetting_count_array)):
continue
real_target=targetting_count_array[target_cell]
while real_target==0:
step=np.random.uniform(-1.5*100,1.5*100)
cell_step=step/10
target_cell=np.int(np.round(cell_pos+cell_step))
if target_cell<0 or (target_cell>=len(targetting_count_array)):
continue
real_target=targetting_count_array[target_cell]
connectivity.append([i,real_target])
delay.append((np.abs(np.round(step*um/speed,decimals=4))))
connectivity=np.array(connectivity)
return connectivity,delay
def randomize_and_edge_checkSC(real_i,i,targetting_count_array,number_of_neurons):
check=np.random.normal(scale=1200)
target=np.round(check/10)
temp_tar=real_i+target
if ((i+target)<0) or ((temp_tar)>=len(targetting_count_array)):
checker=999999999
return checker,check,12345678
checker=targetting_count_array[temp_tar]
if checker>number_of_neurons:
checker=999999999
return checker,check,12345678
return checker,check,temp_tar
def schafer_collaterals(number_of_neurons,exkij,exkijsd,inkj,inkijsd):
connectivity=[]
connectivityin=[]
delayex=[]
delayin=[]
targetting_count_array=make_targetting_array()
for i in range(1,number_of_neurons+1):
if exkij!=0:
numb_connect=np.int(np.round(np.random.normal(loc=exkij,scale=exkijsd)))
else:
numb_connect=0
real_i=np.where(targetting_count_array==i)[0][0]
for j in range(1,numb_connect):
checker,check,temp_tar=randomize_and_edge_checkSC(real_i,i,targetting_count_array,number_of_neurons)
if checker==999999999:
continue
if temp_tar%11==0 and temp_tar!=0:
dupe_rate=13
for k in range(dupe_rate):
connectivityin.append([i,np.int(temp_tar/11)])
delayin.append(np.abs(np.round(check*um/speed,decimals=4)))
else:
dupe_rate=np.abs(np.round(13+np.random.normal(scale=13)))
for k in range(int(dupe_rate)):
connectivity.append([i,checker])
delayex.append(np.abs(np.round(check*um/speed,decimals=4)))
connectivity=np.array(connectivity)
connectivityin=np.array(connectivityin)
delayex=delayex+1000*um/speed
delayin=delayin+1000*um/speed
return connectivity,connectivityin,delayex,delayin
def generate_dentate_gyrus(number_of_cells,frequency_of_spiking,length_of_simulation):
frequency=frequency_of_spiking/1000
firing_array=np.random.rand(number_of_cells,np.int(length_of_simulation/ms))
firing_array=firing_array<frequency
indicies,times=np.where(firing_array)
dentate = SpikeGeneratorGroup(number_of_cells, indicies, times*ms)
return dentate
def main(output_array,innnnn,clust_array):
Pinksy_rinzel_eqs='''
dVs/dt=(-gLs*(Vs-VL)-gNa*(Minfs**2)*hs*(Vs-VNa)-gKdr*ns*(Vs-VK)+(gc/pp)*(Vd-Vs)+(Ip0-Iinssyn)/pp)/Cm : volt
dVd/dt=(-gLd*(Vd-VL)-ICad-gKahp*qd*(Vd-VK)-gKC*cd*chid*(Vd-VK)+(gc*(Vs-Vd))/(1.0-pp)-Isyn/(1.0-pp))/Cm : volt
dCad/dt= -0.13*ICad/uamp/ms*scaler-0.075*Cad/ms : 1
dhs/dt= alphahs-(alphahs+betahs)*hs : 1
dns/dt= alphans-(alphans+betans)*ns : 1
dsd/dt= alphasd-(alphasd+betasd)*sd : 1
dcd/dt= alphacd-(alphacd+betacd)*cd : 1
dqd/dt= alphaqd-(alphaqd+betaqd)*qd : 1
ICad= gCa*sd*sd*(Vd-VCa) : amp
alphams= 0.32*(-46.9-Vs/mV)/(exp((-46.9-Vs/mV)/4.0)-1.0)/ms : Hz
betams= 0.28*(Vs/mV+19.9)/(exp((Vs/mV+19.9)/5.0)-1.0)/ms : Hz
Minfs= alphams/(alphams+betams) : 1
alphans= 0.016*(-24.9-Vs/mV)/(exp((-24.9-Vs/mV)/5.0)-1.0)/ms : Hz
betans= 0.25*exp(-1.0-0.025*Vs/mV)/ms : Hz
alphahs= 0.128*exp((-43.0-Vs/mV)/18.0)/ms : Hz
betahs= 4.0/(1.0+exp((-20.0-Vs/mV)/5.0))/ms : Hz
alphasd= 1.6/(1.0+exp(-0.072*(Vd/mV-5.0)))/ms : Hz
betasd= 0.02*(Vd/mV+8.9)/(exp((Vd/mV+8.9)/5.0)-1.0)/ms : Hz
alphacd=((Vd/mV<=-10)*exp((Vd/mV+50.0)/11-(Vd/mV+53.5)/27)/18.975+(Vd/mV>-10)*2.0*exp((-53.5-Vd/mV)/27.0))/ms : Hz
betacd= ((Vd/mV<=-10)*(2.0*exp((-53.5-Vd/mV)/27.0)-alphacd*ms)+(Vd/mV>-10)*0)/ms : Hz
alphaqd= clip(0.00002*Cad,0,0.01)/ms : Hz
betaqd= 0.001/ms : Hz
chid= clip(Cad/250.0,0,1.0) : 1
INMDA=gNMDA*Si*(1+0.28*exp(-0.062*(Vs/mV)))**(-1)*(Vd-Vsyn) : amp
Hxs=0<=((Vs/mV+50)*1) : 1
Hxw=0<=((Vs/mV+40)*1) : 1
dSi/dt=(Hxs-Si/150)/second : 1
Isyn=Issyn+INMDA+Iapp : amp
Issyn=gbarsyn*clip(ssyn,0,7000)*(Vd-Vsyn): amp
dssyn/dt=-ssyn/tausyn: 1
dinssyn/dt=-inssyn/tauisyn: 1
Iinssyn=gbarsyn*clip(inssyn,0,7000)*(Vs-VK): amp
dIapp/dt=0*amp/ms : amp
'''
CA1_pyramidal = NeuronGroup(number_of_ca1_neurons+1, Pinksy_rinzel_eqs,threshold='Vs>-20*mV',refractory='Vs>-60*mV',method='euler')
CA1_pyramidal.Vs='(randn()*60*.05-60)*mV'
CA1_pyramidal.Vd='(randn()*65*.05-65)*mV'
CA1_pyramidal.hs=.999
CA1_pyramidal.ns=.0001
CA1_pyramidal.sd=.009
CA1_pyramidal.cd=.007
CA1_pyramidal.qd=.01
CA1_pyramidal.Cad=.20
CA1_pyramidal.Si=0
CA1_pyramidal.ssyn=0
CA1_pyramidal.inssyn=0
buzsaki_eqs = '''
dv/dt = (-gNai*m**3*h*(v-VNa)-gKi*n**4*(v-VK)-gL*(v-VL2)-Issyn-Iinssyn-Iapp2)/Cmi : volt
m = alpha_m/(alpha_m+beta_m) : 1
alpha_m = -0.1/mV*(v+35*mV)/(exp(-0.1/mV*(v+35*mV))-1)/ms : Hz
beta_m = 4*exp(-(v+60*mV)/(18*mV))/ms : Hz
dh/dt = 5*(alpha_h*(1-h)-beta_h*h) : 1
alpha_h = 0.07*exp(-(v+58*mV)/(20*mV))/ms : Hz
beta_h = 1./(exp(-0.1/mV*(v+28*mV))+1)/ms : Hz
dn/dt = 5*(alpha_n*(1-n)-beta_n*n) : 1
alpha_n = -0.01/mV*(v+34*mV)/(exp(-0.1/mV*(v+34*mV))-1)/ms : Hz
beta_n = 0.125*exp(-(v+44*mV)/(80*mV))/ms : Hz
Issyn=gbarsyn*(clip(ssyn,0,7000))*(v-Vsyn): amp
dssyn/dt=-ssyn/tausyn: 1
dinssyn/dt=-inssyn/taguisyn: 1
Iinssyn=gbarsyn*clip(inssyn,0,7000)*(v-VK):amp
'''
Ca1_interneuron=NeuronGroup(np.round(number_of_ca1_neurons/10+1),buzsaki_eqs,threshold='v>-20*mV',refractory='v>-60*mV',method='euler')
Ca1_interneuron.v='(randn()*60*.05-60)*mV'
Ca1_interneuron.h=.999
Ca1_interneuron.n=.0001
Ca1_interneuron.ssyn=0
ca1_pypy_kij=0
ca1_pypy_kij_sd=ca1_pypy_kij/20
ca1_pyin_kij=20
ca1_pyin_kij_sd=ca1_pyin_kij/20
connectivity2ex,connectivity2in,delayex,delayin=excite_to_excite_wire_it_up(number_of_ca1_neurons,ca1_pypy_kij,ca1_pypy_kij_sd,ca1_pyin_kij,ca1_pyin_kij_sd)
ca1_inpykij=400
ca1_inpykij_sd=ca1_inpykij/20
ca1in2inkij=100
ca1in2iniij_sd=ca1in2inkij/20
connectivityin2ex,connectivityin2in,delayinex,delayinin=ca1in_wire_it_up(number_of_ca1_neurons,ca1_inpykij,ca1_inpykij_sd,ca1in2inkij,ca1in2iniij_sd)
Ca1_in_in=Synapses(Ca1_interneuron,Ca1_interneuron,pre='inssyn+=5')
Ca1_in_in.connect(i=list(connectivityin2in[0:,0]),j=list(connectivityin2in[0:,1]))
Ca1_in_in.delay=delayinin[0:]
Ca1_py_in=Synapses(CA1_pyramidal,Ca1_interneuron,pre='ssyn+=2.5')
Ca1_py_in.connect(i=list(connectivity2in[0:,0]),j=list(connectivity2in[0:,1]))
Ca1_py_in.delay=delayin[0:]
Ca1_in_py=Synapses(Ca1_interneuron,CA1_pyramidal,pre='inssyn+=15')
Ca1_in_py.connect(i=list(connectivityin2ex[0:,0]),j=list(connectivityin2ex[0:,1]))
Ca1_in_py.delay=delayinex[0:]
number_of_ca3_neurons=1000
CA3_pyramidal = NeuronGroup(number_of_ca3_neurons+1, Pinksy_rinzel_eqs,threshold='Vs>-20*mV',refractory='Vs>-60*mV',method='euler')
CA3_pyramidal.Vs='(randn()*60*.05-60)*mV'
CA3_pyramidal.Vd='(randn()*65*.05-65)*mV'
CA3_pyramidal.hs=.999
CA3_pyramidal.ns=.0001
CA3_pyramidal.sd=.009
CA3_pyramidal.cd=.007
CA3_pyramidal.qd=.01
CA3_pyramidal.Cad=.20
CA3_pyramidal.Si=0
CA3_pyramidal.ssyn=0
CA3_pyramidal.inssyn=0
CA3_pyramidal.Iapp=0*nA
Ca3_interneuron=NeuronGroup(np.round(number_of_ca3_neurons/10+1),buzsaki_eqs,threshold='v>-20*mV',refractory='v>-60*mV',method='euler')
Ca3_interneuron.v=-60*mV
Ca3_interneuron.h=.999
Ca3_interneuron.n=.0001
Ca3_interneuron.ssyn=0
ca3_pypy_kij=55
ca3_pypy_kij_sd=ca3_pypy_kij/20
ca3_pyin_kij=5
ca3_pyin_kij_sd=ca3_pyin_kij/20
connectivity2ex,connectivity2in,delayex,delayin=excite_to_excite_wire_it_up(number_of_ca3_neurons,ca3_pypy_kij,ca3_pypy_kij_sd,ca3_pyin_kij,ca3_pyin_kij_sd)
ca3_inpykij=68
ca3_inpykij_sd=ca3_inpykij/20
connectivityin2ex,delayinin=in2ex_wire_it_up(number_of_ca3_neurons,ca3_inpykij,ca3_inpykij_sd)
Ca3_py_in=Synapses(CA3_pyramidal,Ca3_interneuron,pre='ssyn+=3')
Ca3_py_in.connect(i=list(connectivity2in[0:,0]),j=list(connectivity2in[0:,1]))
Ca3_py_in.delay=delayin[0:]
print 'ca3'
Ca3_in_py=Synapses(Ca3_interneuron,CA3_pyramidal,pre='inssyn+=50')
Ca3_in_py.connect(i=list(connectivityin2ex[0:,0]),j=list(connectivityin2ex[0:,1]))
exkijsc=130
exkijsdsc=exkijsc/20
scconnectivityex,scconnectivityin,scdelayex,scdelayin=schafer_collaterals(number_of_ca1_neurons,exkijsc,exkijsdsc,0,0)
SC_py_py=Synapses(CA3_pyramidal,CA1_pyramidal,pre='ssyn+=15')
SC_py_py.connect(i=list(scconnectivityex[0:,0]),j=list(scconnectivityex[0:,1]))
SC_py_py.delay=scdelayex[0:]
SC_py_in=Synapses(CA3_pyramidal,Ca1_interneuron,pre='ssyn+=15')
SC_py_in.connect(i=list(scconnectivityin[0:,0]),j=list(scconnectivityin[0:,1]))
SC_py_in.delay=scdelayin[0:]
strength_of_input=1
Ca3str=str(np.round(strength_of_input*15,1))
Ca3_py_py=Synapses(CA3_pyramidal,CA3_pyramidal,pre='ssyn+='+Ca3str)
Ca3_py_py.connect(i=list(connectivity2ex[0:,0]),j=list(connectivity2ex[0:,1]))
Ca3_py_py.delay=delayex[0:]
simulation_len=10000*ms
net=Network(Ca3_interneuron,CA1_pyramidal,Ca1_interneuron,CA3_pyramidal,Ca1_py_in,Ca1_in_in,Ca1_in_py,Ca3_py_in,Ca3_in_py,SC_py_py,SC_py_in)
net.add(SpikeMonitor(CA1_pyramidal))
net.add(SpikeMonitor(CA3_pyramidal))
net.add(SpikeMonitor(Ca1_interneuron))
net.add(SpikeMonitor(Ca3_interneuron))
net.add(StateMonitor(CA1_pyramidal, True ,record=True,dt=.5*ms))
net.add(StateMonitor(Ca1_interneuron,True,record=True,dt=.5*ms))
net.add(StateMonitor(CA3_pyramidal, True ,record=True,dt=.5*ms))
net.add(StateMonitor(Ca3_interneuron,True,record=True,dt=.5*ms))
out_spikes_ca1=[]
out_spikes_ca3=[]
med_clust=[]
clust_size=[]
rcParams['figure.figsize'] = 5, 5
#for i in np.arange(6,7,2)*.01:
for i in [15,16,17,18,19,20,21,22,23,24]:
i=i*.01
print str(i)+'nA'
net['neurongroup_2'].Iapp=-1*i*nA
net.store()
net.store()
print len(net.objects)
net.run(simulation_len, report='text')
plt.cla()
plt.clf()
plot(net['spikemonitor'].t/ms,net['spikemonitor'].i,'.k')
plt.axis('off')
print(max(net['spikemonitor'].t/ms))
plt.show()
plt.axis((0,simulation_len/ms,0,1000))
plt.savefig('nnnlowstartdgcellsconnect'+str(i)+str(innnnn)+'.png')
plt.cla()
plt.clf()
plt.clf()
plot(net['spikemonitor_1'].t/ms,net['spikemonitor_1'].i,'.k')
plt.axis('off')
print(max(net['spikemonitor_1'].t/ms))
plt.axis((0,simulation_len/ms,0,1000))
plt.show()
plt.savefig('nnnca3lowstartdgcellsconnect'+str(i)+str(innnnn)+'.png')
plt.cla()
plt.clf()
plt.plot(net['statemonitor'][0].t,(np.sum(net['statemonitor'].Vd[500:550],0)+np.sum(net['statemonitor_1'].v[50:55],0))/55/mV)
plt.axis((0,5,-70,-10))
plt.savefig('nnnlowstartrecordingdgcellsconnect'+str(i)+str(innnnn)+'.eps')
plt.cla()
plt.clf()
#plt.plot(net['statemonitor'][0].t,(np.sum(net['statemonitor_2'].Vd[500:550],0)+np.sum(net['statemonitor_3'].v[50:55],0))/55/mV)
#plt.axis((0,5,-70,-10))
#plt.savefig('nnnlowstartrecordinca3cellsconnect'+str(i)+str(innnnn)+'.eps')
spikes_ca1=np.where(np.histogram(net['spikemonitor'].t,bins=200)[0]>500)[0]
filt_spikes=[]
for sn,spike in enumerate(spikes_ca1):
if spikes_ca1[sn]-spikes_ca1[sn-1]<3:
continue
else:
filt_spikes.append(spike)
spikes_ca3=np.where(np.histogram(net['spikemonitor_1'].t,bins=200)[0]>100)[0]
filt_spikes_ca3=[]
for sn,spike in enumerate(spikes_ca3):
if spikes_ca3[sn]-spikes_ca3[sn-1]<3:
continue
else:
filt_spikes_ca3.append(spike)
out_spikes_ca1.append(len(filt_spikes))
out_spikes_ca3.append(len(filt_spikes_ca3))
nnzoom=np.array([net['spikemonitor_1'].t/ms,net['spikemonitor_1'].i])
nnzoom=nnzoom.T
med_clus_size,num_clust=cluster_algorithim(nnzoom,i,innnnn)
net.restore()
med_clust.append(med_clus_size)
clust_size.append(num_clust)
output_array.append([out_spikes_ca1,out_spikes_ca3])
clust_array.append([med_clust,clust_size])
del net
return output_array,clust_array
speed=500*um/ms
number_of_ca1_neurons=1000
kij={'ca3pypy':55,'ca3pyin':5}
defaultclock.dt = 0.05*ms
'''base values requireed for PR cells'''
Isyn=0*uamp
areapyr=50000*um**2
test=areapyr
Ip0=-0.3*nA
gLs=0.1 *msiemens*cm**-2*test
gLd=0.1 *msiemens*cm**-2*test
gNa=30 *msiemens*cm**-2*test
gKdr=15 *msiemens*cm**-2*test
gCa=7 *msiemens*cm**-2*test
gKahp=0.8*msiemens*cm**-2*test
gKC=15 *msiemens*cm**-2*test
VNa=60*mV
VCa=80*mV
VK=-75*mV
VL=-60*mV
Vsyn=0*mV
vtsyn=0*mV
gc=1.75*usiemens
pp=0.5
Cm=3 *uF*cm**-2*test
gNMDA=0*msiemens*cm**-2*test
gampa=0*msiemens*cm**-2*test
gbarsyn=1*nS
tausyn= 2*ms
tauisyn=7*ms
scaler=(10**8*um**2/(test))
'''base values required for buzsaki cells'''
interneuron_size=20000*um**2
Cmi = 1*uF/cm**2*interneuron_size
gL = 0.1*msiemens/cm**2*interneuron_size
VL2 = -65*mV
gNai = 35*msiemens/cm**2*interneuron_size
gKi = 9*msiemens/cm**2*interneuron_size
Esyn=-75*mV
taguisyn=2*ms
Iapp2=0*nA
'''setting up empty data arrays for data storage'''
output_array=[]
clust_array=[]
'''this allows the simulaiton to be run multiple times to capture a range of randomized connectivity values'''
for i in range(10):
gc11.collect()
output_array,clust_array=main(output_array,i,clust_array)