File: examples-frompapers_Rossant_et_al_2011.txt

package info (click to toggle)
brian 1.4.3-1
  • links: PTS, VCS
  • area: main
  • in suites: sid, stretch
  • size: 23,436 kB
  • sloc: python: 68,707; cpp: 29,040; ansic: 5,182; sh: 111; makefile: 61
file content (142 lines) | stat: -rw-r--r-- 3,952 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
.. currentmodule:: brian

.. index::
   pair: example usage; subplot
   pair: example usage; plot
   pair: example usage; PopulationRateMonitor
   pair: example usage; run
   pair: example usage; SpikeGeneratorGroup
   pair: example usage; xlim
   pair: example usage; show
   pair: example usage; sum
   pair: example usage; title
   pair: example usage; Connection
   pair: example usage; rate
   pair: example usage; reinit_default_clock
   pair: example usage; figure
   pair: example usage; NeuronGroup
   pair: example usage; ylim
   pair: example usage; mean

.. _example-frompapers_Rossant_et_al_2011:

Example: Rossant_et_al_2011 (frompapers)
========================================

Coincidence detection example
=============================
Fig. 4 from:

    Rossant C, Leijon S, Magnusson AK, Brette R (2011).
    "Sensitivity of noisy neurons to coincident inputs".
    Journal of Neuroscience, 31(47).

Two distant or coincident spikes are injected into a noisy balanced
leaky integrate-and-fire neuron. The PSTH of the neuron in response to 
these inputs is calculated along with the extra number of spikes 
in the two cases. This number is higher for the coincident spikes,
showing the sensitivity of a noisy neuron to coincident inputs.

::

    from brian import *
    import matplotlib.patches as patches
    import matplotlib.path as path
    
    def histo(bins, cc, ax):
        # get the corners of the rectangles for the histogram
        left = array(bins[:-1])
        right = array(bins[1:])
        bottom = zeros(len(left))
        top = bottom + cc
        
        # we need a (numrects x numsides x 2) numpy array for the path helper
        # function to build a compound path
        XY = array([[left,left,right,right], [bottom,top,top,bottom]]).T
        
        # get the Path object
        barpath = path.Path.make_compound_path_from_polys(XY)
        
        # make a patch out of it
        patch = patches.PathPatch(barpath, facecolor='blue', edgecolor='gray', alpha=0.8)
        ax.add_patch(patch)
        
        # update the view limits
        ax.set_xlim(left[0], right[-1])
        ax.set_ylim(bottom.min(), top.max())
    
    # neuron parameters
    theta = -55*mV
    vmean = -65*mV
    taum = 5*ms
    taue = 3*ms
    taun = 15*ms
    sigma = 4*mV
    
    # input times
    t1 = 100*ms
    t2 = 120*ms
    
    # simulation duration
    dur = 200*ms
    
    # number of neuron
    N = 10000
    bin = 2*ms
    
    # EPSP size
    int_EPSP=taue
    int_EPSP2=taue*taue/(2*(taum+taue))
    max_EPSP=(taum/taue)**(taum/(taue-taum))
    we = 3.0*mV/max_EPSP
    
    # model equations
    eqs = '''
    V=V0+noise : volt
    dV0/dt=(-V0+psp)/taum : volt
    dpsp/dt=-psp/taue : volt
    dnoise/dt=(vmean-noise)/taun+sigma*(2./taun)**.5*xi : volt
    '''
    threshold = 'V>theta'
    reset = vmean
    
    # initialization of the NeuronGroup
    reinit_default_clock()
    group = NeuronGroup(2*N, model=eqs, reset=reset, threshold=threshold)
    group.V0 = group.psp = 0*volt
    group.noise = vmean + sigma * randn(2*N)
    
    # input spikes
    input_spikes = [(0, t1), (0, t2), (1, t1)]
    input = SpikeGeneratorGroup(2, array(input_spikes))
    
    # connections
    C = Connection(input, group, 'psp')
    C.connect_full(input[0], group[:N], weight=we)
    C.connect_full(input[1], group[N:], weight=2*we)
    
    # monitors
    prM1 = PopulationRateMonitor(group[:N], bin=bin)
    prM2 = PopulationRateMonitor(group[N:], bin=bin)
    
    # launch simulation
    run(dur)
    
    # PSTH plot
    figure(figsize=(10,10))
    prMs = [prM1, prM2]
    for i in [0,1]:
        prM = prMs[i]
        r = prM.rate[:-1]*bin
        m = mean(r[:len(r)/2])
    
        ax = subplot(211+i)
        histo(prM.times, r, ax)
        plot([0,dur],[m,m],'--r')
        title("%.2f extra spikes" % sum(r[t1/bin:(t2+20*ms)/bin]-m))
        xlim(.05, .2)
        ylim(0, .125)
    
    show()