Skip to content

neuron

Compartment model representing a single neuron in NEURON.

PARAMETER DESCRIPTION
x
    Neuron to generate model for. Has to be in microns!

TYPE: navis.TreeNeuron

res
    Approximate length [um] of segments. This guarantees that
    no section has any segment that is longer than `res` but for
    small branches (i.e. "sections") the segments might be smaller.
    Lower `res` = more detailed simulation.

TYPE: int DEFAULT: 10

Source code in navis/interfaces/neuron/comp.py
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
class CompartmentModel:
    """Compartment model representing a single neuron in NEURON.

    Parameters
    ----------
    x :         navis.TreeNeuron
                Neuron to generate model for. Has to be in microns!
    res :       int
                Approximate length [um] of segments. This guarantees that
                no section has any segment that is longer than `res` but for
                small branches (i.e. "sections") the segments might be smaller.
                Lower `res` = more detailed simulation.

    """

    def __init__(self, x: 'core.TreeNeuron', res=10):
        """Initialize Neuron."""
        utils.eval_param(x, name='x', allowed_types=(core.TreeNeuron, ))

        # Note that we make a copy to make sure that the data underlying the
        # model will not accidentally be changed
        self.skeleton = x.copy()

        # Max section resolution per segment
        self.res = res

        # Some placeholders
        self._sections = []
        self._stimuli = {}
        self._records = {}
        self._synapses = {}

        # Generate the actual model
        self._validate_skeleton()
        self._generate_sections()

    def __repr__(self):
        s = (f'CompartmentModel<id={self.skeleton.label},'
             f'sections={self.n_sections};'
             f'stimuli={self.n_stimuli};'
             f'records={self.n_records}>'
             )
        return s

    @property
    def label(self):
        """Name/label of the neuron."""
        return f'CompartmentModel[{self.skeleton.label}]'

    @property
    def n_records(self):
        """Number of records (across all types) active on this model."""
        return len([r for t in self.records.values() for r in t])

    @property
    def n_sections(self):
        """Number of sections in this model."""
        return len(self.sections)

    @property
    def n_stimuli(self):
        """Number of stimuli active on this model."""
        return len(self.stimuli)

    @property
    def nodes(self) -> pd.DataFrame:
        """Node table of the skeleton."""
        return self.skeleton.nodes

    @property
    def cm(self) -> float:
        """Membran capacity [micro Farads / cm^2] of all sections."""
        return np.array([s.cm for s in self.sections])

    @cm.setter
    def cm(self, value: float):
        """Membran capacity [micro Farads / cm^2] for all sections."""
        for s in self.sections:
            s.cm = value

    @property
    def Ra(self) -> float:
        """Axial resistance [Ohm * cm] of all sections."""
        return np.array([s.Ra for s in self.sections])

    @Ra.setter
    def Ra(self, value: float):
        """Set axial resistance [Ohm * cm] for all sections."""
        for s in self.sections:
            s.Ra = value

    @property
    def records(self) -> dict:
        """Return mapping of node ID(s) to recordings."""
        return self._records

    @property
    def sections(self) -> np.ndarray:
        """List of sections making up this model."""
        return self._sections

    @property
    def stimuli(self) -> dict:
        """Return mapping of node ID(s) to stimuli."""
        return self._stimuli

    @property
    def synapses(self) -> dict:
        """Return mapping of node ID(s) to synapses."""
        return self._synapses

    @property
    def t(self) -> np.ndarray:
        """The global time. Should be the same for all neurons."""
        return main_t

    def _generate_sections(self):
        """Generate sections from the neuron.

        This will automatically be called at initialization and should not be
        called again.

        """
        # First generate sections
        self._sections = []
        nodes = self.skeleton.nodes.set_index('node_id')
        roots = self.skeleton.root
        bp = self.skeleton.branch_points.node_id.values
        G = self.skeleton.graph
        node2sec = {}
        node2pos = {}
        for i, seg in enumerate(self.skeleton.small_segments):
            # Get child -> parent distances in this segment
            dists = np.array([G.edges[(c, p)]['weight']
                              for c, p in zip(seg[:-1], seg[1:])])

            # Invert the sections
            # That's because in navis sections go from tip -> root (i.e.
            # child -> parent) but in neuron section(0) is the base and
            # section(1) is the tip.
            seg = np.asarray(seg)[::-1]
            dists = dists[::-1]

            # Grab the coordinates and radii
            seg_nodes = nodes.loc[seg]
            locs = seg_nodes[['x', 'y', 'z']].values
            radii = seg_nodes.radius.values

            # Generate section
            sec = neuron.h.Section(name=f'segment_{i}')

            # Set 3D points -> this automatically sets length L
            xvec = neuron.h.Vector(locs[:, 0])
            yvec = neuron.h.Vector(locs[:, 1])
            zvec = neuron.h.Vector(locs[:, 2])
            dvec = neuron.h.Vector(radii * 2)
            neuron.h.pt3dadd(xvec, yvec, zvec, dvec, sec=sec)

            # Set number of segments for this section
            # We also will make sure that each section has an odd
            # number of segments
            sec.nseg = 1 + 2 * int(sec.L / (self.res * 2))
            # Keep track of section
            self.sections.append(sec)

            # While we're at it: for each point (except the root of this
            # section) find the relative position within the section

            # Get normalized positions within this segment
            norm_pos = dists.cumsum() / dists.sum()

            # Update positional dictionaries (required for connecting the
            # segments in the next step)
            node2pos.update(dict(zip(seg[1:], norm_pos)))
            node2sec.update({n: i for n in seg[1:]})

            # If this happens to be the segment with the skeleton's root, keep
            # track of it too
            if seg[0] in roots:
                node2pos[seg[0]] = 0
                node2sec[seg[0]] = i

        self._sections = np.array(self.sections)
        self.skeleton.nodes['sec_ix'] = self.skeleton.nodes.node_id.map(node2sec)
        self.skeleton.nodes['sec_pos'] = self.skeleton.nodes.node_id.map(node2pos)

        # Need to grab nodes again after adding `sec_ix` and `sec_pos`
        nodes = self.skeleton.nodes.set_index('node_id')

        # Connect segments
        for i, seg in enumerate(self.skeleton.small_segments):
            # Root is special in that it only needs to be connected if it's also
            # a branch point
            if seg[-1] in roots:
                # Skip if root is not a branch point
                if seg[-1] not in bp:
                    continue
                # If root is also a branch point, it will be part of more than
                # one section but in the positional dicts we will have kept track
                # of only one of them. That's the one we pick as base segment
                if node2sec[seg[-1]] == i:
                    continue

            parent = nodes.loc[seg[-1]]
            parent_sec = self.sections[parent.sec_ix]
            self.sections[i].connect(parent_sec(1))

    def _validate_skeleton(self):
        """Validate skeleton."""
        if self.skeleton.units and not self.skeleton.units.dimensionless:
            not_um = self.skeleton.units.units != config.ureg.Unit('um')
            not_microns = self.skeleton.units.units != config.ureg.Unit('microns')
            if not_um and not_microns:
                logger.warning('Model expects coordinates in microns but '
                               f'neuron has units "{self.skeleton.units}"!')

        if len(self.skeleton.root) > 1:
            logger.warning('Neuron has multiple roots and hence consists of '
                           'multiple disconnected fragments!')

        if 'radius' not in self.skeleton.nodes.columns:
            raise ValueError('Neuron node table must have `radius` column')

        if np.any(self.skeleton.nodes.radius.values <= 0):
            raise ValueError('Neuron node table contains radii <= 0.')

    def add_synaptic_input(self, where, start=5 * ms,
                           spike_no=1, spike_int=10 * ms, spike_noise=0,
                           syn_tau1=.1 * ms, syn_tau2=10 * ms, syn_rev_pot=0,
                           cn_thresh=10, cn_delay=1 * ms, cn_weight=0.05):
        """Add synaptic input to model.

        This uses the Exp2Syn synapse. All targets in `where` are triggered
        by the same NetStim - i.e. they will all receive their spike(s) at the
        same time.

        Parameters
        ----------
        where :         int | list of int
                        Node IDs at which to simulate synaptic input.

        Properties for presynaptic spikes:

        start :         int
                        Onset [ms] of first spike from beginning of simulation.
        spike_no :      int
                        Number of presynaptic spikes to produce.
        spike_int :     int
                        Interval [ms] between consecutive spikes.
        spike_noise :   float [0-1]
                        Fractional randomness in spike timing.

        Synapse properties:

        syn_tau1 :      int
                        Rise time constant [ms].
        syn_tau2 :      int
                        Decay time constant [ms].
        syn_rev_pot :   int
                        Reversal potential (e) [mV].

        Connection properties:

        cn_thresh :     int
                        Presynaptic membrane potential [mV] at which synaptic
                        event is triggered.
        cn_delay :      int
                        Delay [ms] between presynaptic trigger and postsynaptic
                        event.
        cn_weight :     float
                        Weight variable. This bundles a couple of synaptic
                        properties such as e.g. how much transmitter is released
                        or binding affinity at postsynaptic receptors.

        """
        where = utils.make_iterable(where)

        # Make a new stimulator
        stim = neuron.h.NetStim()
        stim.number = spike_no
        stim.start = start
        stim.noise = spike_noise
        stim.interval = spike_int

        # Connect
        self.connect(stim, where, syn_tau1=syn_tau1, syn_tau2=syn_tau2,
                     syn_rev_pot=syn_rev_pot, cn_thresh=cn_thresh,
                     cn_delay=cn_delay, cn_weight=cn_weight)

    def inject_current_pulse(self, where, start=5,
                             duration=1, current=0.1):
        """Add current injection (IClamp) stimulation to model.

        Parameters
        ----------
        where :     int | list of int
                    Node ID(s) at which to stimulate.
        start :     int
                    Onset (delay) [ms] from beginning of simulation.
        duration :  int
                    Duration (dur) [ms] of injection.
        current :   float
                    Amount (i) [nA] of injected current.

        """
        self._add_stimulus('IClamp', where=where, delay=start,
                           dur=duration, amp=current)

    def add_synaptic_current(self, where, start=5, tau=0.1, rev_pot=0,
                             max_syn_cond=0.1):
        """Add synaptic current(s) (AlphaSynapse) to model.

        Parameters
        ----------
        where :         int | list of int
                        Node ID(s) at which to stimulate.
        start :         int
                        Onset [ms] from beginning of simulation.
        tau :           int
                        Decay time constant [ms].
        rev_pot :       int
                        Reverse potential (e) [mV].
        max_syn_cond :  float
                        Max synaptic conductance (gmax) [uS].

        """
        self._add_stimulus('AlphaSynapse', where=where, onset=start,
                           tau=tau, e=rev_pot, gmax=max_syn_cond)

    def _add_stimulus(self, stimulus, where, **kwargs):
        """Add generic stimulus."""
        if not callable(stimulus):
            stimulus = getattr(neuron.h, stimulus)

        where = utils.make_iterable(where)

        nodes = self.nodes.set_index('node_id')
        for node in nodes.loc[where].itertuples():
            sec = self.sections[node.sec_ix](node.sec_pos)
            stim = stimulus(sec)

            for k, v in kwargs.items():
                setattr(stim, k, v)

            self.stimuli[node.Index] = self.stimuli.get(node.Index, []) + [stim]

    def add_voltage_record(self, where, label=None):
        """Add voltage recording to model.

        Parameters
        ----------
        where :     int | list of int
                    Node ID(s) at which to record.
        label :     str, optional
                    If label is given, this recording will be added as
                    `self.records['v'][label]` else  `self.records['v'][node_id]`.

        """
        self._add_record(where, what='v', label=label)

    def add_current_record(self, where, label=None):
        """Add current recording to model.

        This only works if nodes map to sections that have point processes.

        Parameters
        ----------
        where :     int | list of int
                    Node ID(s) at which to record.
        label :     str, optional
                    If label is given, this recording will be added as
                    `self.records['i'][label]` else  `self.records['i'][node_id]`.

        """
        nodes = utils.make_iterable(where)

        # Map nodes to point processes
        secs = self.get_node_segment(nodes)
        where = []
        for n, sec in zip(nodes, secs):
            pp = sec.point_processes()
            if not pp:
                raise TypeError(f'Section for node {n} has no point process '
                                '- unable to add current record')
            elif len(pp) > 1:
                logger.warning(f'Section for node {n} has more than on point '
                               'process. Recording current at first.')
                pp = pp[:1]
            where += pp

        self._add_record(where, what='i', label=label)

    def add_spike_detector(self, where, threshold=20, label=None):
        """Add a spike detector at given node(s).

        Parameters
        ----------
        where :     int | list of int
                    Node ID(s) at which to record.
        threshold : float
                    Threshold in mV for a spike to be counted.
        label :     str, optional
                    If label is given, this recording will be added as
                    `self.records[label]` else  `self.records[node_id]`.

        """
        where = utils.make_iterable(where)

        self.records['spikes'] = self.records.get('spikes', {})
        self._spike_det = getattr(self, '_spike_det', [])
        segments = self.get_node_segment(where)
        sections = self.get_node_section(where)
        for n, sec, seg in zip(where, sections, segments):
            # Generate a NetCon object that has no target
            sp_det = neuron.h.NetCon(seg._ref_v, None, sec=sec)

            # Set threshold
            if threshold:
                sp_det.threshold = threshold

            # Keeping track of this to save it from garbage collector
            self._spike_det.append(sp_det)

            # Create a vector for the spike timings
            vec = neuron.h.Vector()
            # Tell the NetCon object to record into that vector
            sp_det.record(vec)

            if label:
                self.records['spikes'][label] = vec
            else:
                self.records['spikes'][n] = vec

    def _add_record(self, where, what, label=None):
        """Add a recording to given node.

        Parameters
        ----------
        where :     int | list of int | point process | section
                    Node ID(s) (or a section) at which to record.
        what :      str
                    What to record. Can be e.g. `v` or `_ref_v` for Voltage.
        label :     str, optional
                    If label is given, this recording will be added as
                    `self.records[label]` else  `self.records[node_id]`.

        """
        where = utils.make_iterable(where)

        if not isinstance(what, str):
            raise TypeError(f'Required str e.g. "v", got {type(what)}')

        if not what.startswith('_ref_'):
            what = f'_ref_{what}'

        rec_type = what.split('_')[-1]
        if rec_type not in self.records:
            self.records[rec_type] = {}

        # # Get node segments only for nodes
        is_node = ~np.array([is_NEURON_object(w) for w in where])
        node_segs = np.zeros(len(where), dtype=object)
        node_segs[is_node] = self.get_node_segment(where[is_node])

        for i, w in enumerate(where):
            # If this is a neuron object (e.g. segment, section or point
            # process) we assume this does not need mapping
            if is_NEURON_object(w):
                seg = w
            else:
                seg = node_segs[i]

            rec = neuron.h.Vector().record(getattr(seg, what))

            if label:
                self.records[rec_type][label] = rec
            else:
                self.records[rec_type][w] = rec

    def connect(self, pre, where, syn_tau1=.1 * ms, syn_tau2=10 * ms,
                syn_rev_pot=0, cn_thresh=10, cn_delay=1 * ms, cn_weight=0):
        """Connect object to model.

        This uses the Exp2Syn synapse and treats `pre` as the presynaptic
        object.

        Parameters
        ----------
        pre :           NetStim | section
                        The presynaptic object to connect to this neuron.
        where :         int | list of int
                        Node IDs at which to simulate synaptic input.

        Synapse properties:

        syn_tau1 :      int
                        Rise time constant [ms].
        syn_tau2 :      int
                        Decay time constant [ms].
        syn_rev_pot :   int
                        Reversal potential (e) [mV].

        Connection properties:

        cn_thresh :     int
                        Presynaptic membrane potential [mV] at which synaptic
                        event is triggered.
        cn_delay :      int
                        Delay [ms] between presynaptic trigger and postsynaptic
                        event.
        cn_weight :     int
                        Weight variable. This bundles a couple of synaptic
                        properties such as e.g. how much transmitter is released
                        or binding affinity at postsynaptic receptors.

        """
        where = utils.make_iterable(where)

        if not is_NEURON_object(pre):
            raise ValueError(f'Expected NEURON object, got {type(pre)}')

        # Turn section into segment
        if isinstance(pre, neuron.nrn.Section):
            pre = pre()

        # Go over the nodes
        nodes = self.nodes.set_index('node_id')
        for node in nodes.loc[where].itertuples():
            # Generate synapses for the nodes in question
            # Note that we are not reusing existing synapses
            # in case the properties are different
            sec = self.sections[node.sec_ix](node.sec_pos)
            syn = neuron.h.Exp2Syn(sec)
            syn.tau1 = syn_tau1
            syn.tau2 = syn_tau2
            syn.e = syn_rev_pot

            self.synapses[node.Index] = self.synapses.get(node.Index, []) + [syn]

            # Connect spike stimulus and synapse
            if isinstance(pre, neuron.nrn.Segment):
                nc = neuron.h.NetCon(pre._ref_v, syn, sec=pre.sec)
            else:
                nc = neuron.h.NetCon(pre, syn)

            # Set connection parameters
            nc.threshold = cn_thresh
            nc.delay = cn_delay
            nc.weight[0] = cn_weight

            self.stimuli[node.Index] = self.stimuli.get(node.Index, []) + [nc, pre]

    def clear_records(self):
        """Clear records."""
        self._records = {}

    def clear_stimuli(self):
        """Clear stimuli."""
        self._stimuli = {}

    def clear_synapses(self):
        """Clear synapses."""
        self._synapses = {}

    def clear(self):
        """Attempt to remove model from NEURON space.

        This is not guaranteed to work. Check `neuron.h.topology()` to inspect.

        """
        # Basically we have to bring the reference count to zero
        self.clear_records()
        self.clear_stimuli()
        self.clear_synapses()
        for s in self._sections:
            del s
        self._sections = []

    def get_node_section(self, node_ids):
        """Return section(s) for given node(s).

        Parameters
        ----------
        node_ids :  int | list of int
                    Node IDs.

        Returns
        -------
        section(s) :    segment or list of segments
                        Depends on input.

        """
        nodes = self.nodes.set_index('node_id')
        if not utils.is_iterable(node_ids):
            n = nodes.loc[node_ids]
            return self.sections[n.sec_ix]
        else:
            segs = []
            for node in nodes.loc[node_ids].itertuples():
                segs.append(self.sections[node.sec_ix])
            return segs

    def get_node_segment(self, node_ids):
        """Return segment(s) for given node(s).

        Parameters
        ----------
        node_ids :  int | list of int
                    Node IDs.

        Returns
        -------
        segment(s) :    segment or list of segments
                        Depends on input.

        """
        nodes = self.nodes.set_index('node_id')
        if not utils.is_iterable(node_ids):
            n = nodes.loc[node_ids]
            return self.sections[n.sec_ix](n.sec_pos)
        else:
            segs = []
            for node in nodes.loc[node_ids].itertuples():
                segs.append(self.sections[node.sec_ix](node.sec_pos))
            return segs

    def insert(self, mechanism, subset=None, **kwargs):
        """Insert biophysical mechanism for model.

        Parameters
        ----------
        mechanism : str
                    Mechanism to insert - e.g. "hh" for Hodgkin-Huxley kinetics.
        subset :    list of sections | list of int
                    Sections (or indices thereof) to set mechanism for.
                    If `None` will add mechanism to all sections.
        **kwargs
                    Use to set properties for mechanism.

        """
        if isinstance(subset, type(None)):
            sections = self.sections
        else:
            subset = utils.make_iterable(subset)

            if all([is_section(s) for s in subset]):
                sections = subset
            elif all([isinstance(s, Number) for s in subset]):
                sections = self.sections[subset]
            else:
                raise TypeError('`subset` must be None, a list of sections or '
                                'a list of section indices')

        for sec in np.unique(sections):
            _ = sec.insert(mechanism)
            for seg in sec:
                mech = getattr(seg, mechanism)
                for k, v in kwargs.items():
                    setattr(mech, k, v)

    def uninsert(self, mechanism, subset=None):
        """Remove biophysical mechanism from model.

        Parameters
        ----------
        mechanism : str
                    Mechanism to remove - e.g. "hh" for Hodgkin-Huxley kinetics.
        subset :    list of sections | list of int
                    Sections (or indices thereof) to set mechanism for.
                    If `None` will add mechanism to all sections.

        """
        if isinstance(subset, type(None)):
            sections = self.sections
        else:
            subset = utils.make_iterable(subset)

            if all([is_section(s) for s in subset]):
                sections = subset
            elif all([isinstance(s, Number) for s in subset]):
                sections = self.sections[subset]
            else:
                raise TypeError('`subset` must be None, a list of sections or '
                                'a list of section indices')

        for sec in np.unique(sections):
            if hasattr(sec, mechanism):
                _ = sec.uninsert(mechanism)

    def plot_structure(self):
        """Visualize structure in 3D using matplotlib."""
        _ = neuron.h.PlotShape().plot(plt)

    def run_simulation(self, duration=25 * ms, v_init=-65 * mV):
        """Run the simulation."""
        # Add recording of time
        global main_t
        main_t = neuron.h.Vector().record(neuron.h._ref_t)

        # This resets the entire model space not just this neuron!
        neuron.h.finitialize(v_init)
        neuron.h.continuerun(duration)

    def plot_results(self, axes=None):
        """Plot results.

        Parameters
        ----------
        axes :      matplotlib axes
                    Axes to plot onto. Must have one ax for each recording
                    type (mV, spike count, etc) in `self.records`.

        Returns
        -------
        axes

        """
        if isinstance(self.t, type(None)) or not len(self.t):
            logger.warning('Looks like the simulation has not yet been run.')
            return
        if not self.records:
            logger.warning('Nothing to plot: no recordings found.')
            return

        if not axes:
            fig, axes = plt.subplots(len(self.records), sharex=True)

        # Make sure that even a single ax is a list
        if not isinstance(axes, (np.ndarray, list)):
            axes = [axes] * len(self.records)

        for t, ax in zip(self.records, axes):
            for i, (k, v) in enumerate(self.records[t].items()):
                if not len(v):
                    continue
                v = v.as_numpy()
                # For spikes the vector contains the times
                if t == 'spikes':
                    # Calculate spike rate
                    bins = np.linspace(0, max(self.t), 10)
                    hist, _ = np.histogram(v, bins=bins)
                    width = bins[1] - bins[0]
                    rate = hist * (1000 / width)
                    ax.plot(bins[:-1] + (width / 2), rate, label=k)

                    ax.scatter(v, [-i] * len(v), marker='|', s=100)
                else:
                    ax.plot(self.t, v, label=k)

            ax.set_xlabel('time [ms]')
            ax.set_ylabel(f'{t}')

            ax.legend()
        return axes

Axial resistance [Ohm * cm] of all sections.

Membran capacity [micro Farads / cm^2] of all sections.

Name/label of the neuron.

Number of records (across all types) active on this model.

Number of sections in this model.

Number of stimuli active on this model.

Node table of the skeleton.

Return mapping of node ID(s) to recordings.

List of sections making up this model.

Return mapping of node ID(s) to stimuli.

Return mapping of node ID(s) to synapses.

The global time. Should be the same for all neurons.

Initialize Neuron.

Source code in navis/interfaces/neuron/comp.py
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
def __init__(self, x: 'core.TreeNeuron', res=10):
    """Initialize Neuron."""
    utils.eval_param(x, name='x', allowed_types=(core.TreeNeuron, ))

    # Note that we make a copy to make sure that the data underlying the
    # model will not accidentally be changed
    self.skeleton = x.copy()

    # Max section resolution per segment
    self.res = res

    # Some placeholders
    self._sections = []
    self._stimuli = {}
    self._records = {}
    self._synapses = {}

    # Generate the actual model
    self._validate_skeleton()
    self._generate_sections()

Add current recording to model.

This only works if nodes map to sections that have point processes.

PARAMETER DESCRIPTION
where
    Node ID(s) at which to record.

TYPE: int | list of int

label
    If label is given, this recording will be added as
    `self.records['i'][label]` else  `self.records['i'][node_id]`.

TYPE: str DEFAULT: None

Source code in navis/interfaces/neuron/comp.py
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
def add_current_record(self, where, label=None):
    """Add current recording to model.

    This only works if nodes map to sections that have point processes.

    Parameters
    ----------
    where :     int | list of int
                Node ID(s) at which to record.
    label :     str, optional
                If label is given, this recording will be added as
                `self.records['i'][label]` else  `self.records['i'][node_id]`.

    """
    nodes = utils.make_iterable(where)

    # Map nodes to point processes
    secs = self.get_node_segment(nodes)
    where = []
    for n, sec in zip(nodes, secs):
        pp = sec.point_processes()
        if not pp:
            raise TypeError(f'Section for node {n} has no point process '
                            '- unable to add current record')
        elif len(pp) > 1:
            logger.warning(f'Section for node {n} has more than on point '
                           'process. Recording current at first.')
            pp = pp[:1]
        where += pp

    self._add_record(where, what='i', label=label)

Add a spike detector at given node(s).

PARAMETER DESCRIPTION
where
    Node ID(s) at which to record.

TYPE: int | list of int

threshold
    Threshold in mV for a spike to be counted.

TYPE: float DEFAULT: 20

label
    If label is given, this recording will be added as
    `self.records[label]` else  `self.records[node_id]`.

TYPE: str DEFAULT: None

Source code in navis/interfaces/neuron/comp.py
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
def add_spike_detector(self, where, threshold=20, label=None):
    """Add a spike detector at given node(s).

    Parameters
    ----------
    where :     int | list of int
                Node ID(s) at which to record.
    threshold : float
                Threshold in mV for a spike to be counted.
    label :     str, optional
                If label is given, this recording will be added as
                `self.records[label]` else  `self.records[node_id]`.

    """
    where = utils.make_iterable(where)

    self.records['spikes'] = self.records.get('spikes', {})
    self._spike_det = getattr(self, '_spike_det', [])
    segments = self.get_node_segment(where)
    sections = self.get_node_section(where)
    for n, sec, seg in zip(where, sections, segments):
        # Generate a NetCon object that has no target
        sp_det = neuron.h.NetCon(seg._ref_v, None, sec=sec)

        # Set threshold
        if threshold:
            sp_det.threshold = threshold

        # Keeping track of this to save it from garbage collector
        self._spike_det.append(sp_det)

        # Create a vector for the spike timings
        vec = neuron.h.Vector()
        # Tell the NetCon object to record into that vector
        sp_det.record(vec)

        if label:
            self.records['spikes'][label] = vec
        else:
            self.records['spikes'][n] = vec

Add synaptic current(s) (AlphaSynapse) to model.

PARAMETER DESCRIPTION
where
        Node ID(s) at which to stimulate.

TYPE: int | list of int

start
        Onset [ms] from beginning of simulation.

TYPE: int DEFAULT: 5

tau
        Decay time constant [ms].

TYPE: int DEFAULT: 0.1

rev_pot
        Reverse potential (e) [mV].

TYPE: int DEFAULT: 0

max_syn_cond
        Max synaptic conductance (gmax) [uS].

TYPE: float DEFAULT: 0.1

Source code in navis/interfaces/neuron/comp.py
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
def add_synaptic_current(self, where, start=5, tau=0.1, rev_pot=0,
                         max_syn_cond=0.1):
    """Add synaptic current(s) (AlphaSynapse) to model.

    Parameters
    ----------
    where :         int | list of int
                    Node ID(s) at which to stimulate.
    start :         int
                    Onset [ms] from beginning of simulation.
    tau :           int
                    Decay time constant [ms].
    rev_pot :       int
                    Reverse potential (e) [mV].
    max_syn_cond :  float
                    Max synaptic conductance (gmax) [uS].

    """
    self._add_stimulus('AlphaSynapse', where=where, onset=start,
                       tau=tau, e=rev_pot, gmax=max_syn_cond)

Add synaptic input to model.

This uses the Exp2Syn synapse. All targets in where are triggered by the same NetStim - i.e. they will all receive their spike(s) at the same time.

PARAMETER DESCRIPTION
where
        Node IDs at which to simulate synaptic input.

TYPE: int | list of int

Properties

start
        Onset [ms] of first spike from beginning of simulation.

TYPE: int DEFAULT: 5 * ms

spike_no
        Number of presynaptic spikes to produce.

TYPE: int DEFAULT: 1

spike_int
        Interval [ms] between consecutive spikes.

TYPE: int DEFAULT: 10 * ms

spike_noise
        Fractional randomness in spike timing.

TYPE: float [0-1] DEFAULT: 0

Synapse

syn_tau1
        Rise time constant [ms].

TYPE: int DEFAULT: 0.1 * ms

syn_tau2
        Decay time constant [ms].

TYPE: int DEFAULT: 10 * ms

syn_rev_pot
        Reversal potential (e) [mV].

TYPE: int DEFAULT: 0

Connection

cn_thresh
        Presynaptic membrane potential [mV] at which synaptic
        event is triggered.

TYPE: int DEFAULT: 10

cn_delay
        Delay [ms] between presynaptic trigger and postsynaptic
        event.

TYPE: int DEFAULT: 1 * ms

cn_weight
        Weight variable. This bundles a couple of synaptic
        properties such as e.g. how much transmitter is released
        or binding affinity at postsynaptic receptors.

TYPE: float DEFAULT: 0.05

Source code in navis/interfaces/neuron/comp.py
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
def add_synaptic_input(self, where, start=5 * ms,
                       spike_no=1, spike_int=10 * ms, spike_noise=0,
                       syn_tau1=.1 * ms, syn_tau2=10 * ms, syn_rev_pot=0,
                       cn_thresh=10, cn_delay=1 * ms, cn_weight=0.05):
    """Add synaptic input to model.

    This uses the Exp2Syn synapse. All targets in `where` are triggered
    by the same NetStim - i.e. they will all receive their spike(s) at the
    same time.

    Parameters
    ----------
    where :         int | list of int
                    Node IDs at which to simulate synaptic input.

    Properties for presynaptic spikes:

    start :         int
                    Onset [ms] of first spike from beginning of simulation.
    spike_no :      int
                    Number of presynaptic spikes to produce.
    spike_int :     int
                    Interval [ms] between consecutive spikes.
    spike_noise :   float [0-1]
                    Fractional randomness in spike timing.

    Synapse properties:

    syn_tau1 :      int
                    Rise time constant [ms].
    syn_tau2 :      int
                    Decay time constant [ms].
    syn_rev_pot :   int
                    Reversal potential (e) [mV].

    Connection properties:

    cn_thresh :     int
                    Presynaptic membrane potential [mV] at which synaptic
                    event is triggered.
    cn_delay :      int
                    Delay [ms] between presynaptic trigger and postsynaptic
                    event.
    cn_weight :     float
                    Weight variable. This bundles a couple of synaptic
                    properties such as e.g. how much transmitter is released
                    or binding affinity at postsynaptic receptors.

    """
    where = utils.make_iterable(where)

    # Make a new stimulator
    stim = neuron.h.NetStim()
    stim.number = spike_no
    stim.start = start
    stim.noise = spike_noise
    stim.interval = spike_int

    # Connect
    self.connect(stim, where, syn_tau1=syn_tau1, syn_tau2=syn_tau2,
                 syn_rev_pot=syn_rev_pot, cn_thresh=cn_thresh,
                 cn_delay=cn_delay, cn_weight=cn_weight)

Add voltage recording to model.

PARAMETER DESCRIPTION
where
    Node ID(s) at which to record.

TYPE: int | list of int

label
    If label is given, this recording will be added as
    `self.records['v'][label]` else  `self.records['v'][node_id]`.

TYPE: str DEFAULT: None

Source code in navis/interfaces/neuron/comp.py
465
466
467
468
469
470
471
472
473
474
475
476
477
def add_voltage_record(self, where, label=None):
    """Add voltage recording to model.

    Parameters
    ----------
    where :     int | list of int
                Node ID(s) at which to record.
    label :     str, optional
                If label is given, this recording will be added as
                `self.records['v'][label]` else  `self.records['v'][node_id]`.

    """
    self._add_record(where, what='v', label=label)

Attempt to remove model from NEURON space.

This is not guaranteed to work. Check neuron.h.topology() to inspect.

Source code in navis/interfaces/neuron/comp.py
683
684
685
686
687
688
689
690
691
692
693
694
695
def clear(self):
    """Attempt to remove model from NEURON space.

    This is not guaranteed to work. Check `neuron.h.topology()` to inspect.

    """
    # Basically we have to bring the reference count to zero
    self.clear_records()
    self.clear_stimuli()
    self.clear_synapses()
    for s in self._sections:
        del s
    self._sections = []

Clear records.

Source code in navis/interfaces/neuron/comp.py
671
672
673
def clear_records(self):
    """Clear records."""
    self._records = {}

Clear stimuli.

Source code in navis/interfaces/neuron/comp.py
675
676
677
def clear_stimuli(self):
    """Clear stimuli."""
    self._stimuli = {}

Clear synapses.

Source code in navis/interfaces/neuron/comp.py
679
680
681
def clear_synapses(self):
    """Clear synapses."""
    self._synapses = {}

Connect object to model.

This uses the Exp2Syn synapse and treats pre as the presynaptic object.

PARAMETER DESCRIPTION
pre
        The presynaptic object to connect to this neuron.

TYPE: NetStim | section

where
        Node IDs at which to simulate synaptic input.

TYPE: int | list of int

Synapse

syn_tau1
        Rise time constant [ms].

TYPE: int DEFAULT: 0.1 * ms

syn_tau2
        Decay time constant [ms].

TYPE: int DEFAULT: 10 * ms

syn_rev_pot
        Reversal potential (e) [mV].

TYPE: int DEFAULT: 0

Connection

cn_thresh
        Presynaptic membrane potential [mV] at which synaptic
        event is triggered.

TYPE: int DEFAULT: 10

cn_delay
        Delay [ms] between presynaptic trigger and postsynaptic
        event.

TYPE: int DEFAULT: 1 * ms

cn_weight
        Weight variable. This bundles a couple of synaptic
        properties such as e.g. how much transmitter is released
        or binding affinity at postsynaptic receptors.

TYPE: int DEFAULT: 0

Source code in navis/interfaces/neuron/comp.py
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
def connect(self, pre, where, syn_tau1=.1 * ms, syn_tau2=10 * ms,
            syn_rev_pot=0, cn_thresh=10, cn_delay=1 * ms, cn_weight=0):
    """Connect object to model.

    This uses the Exp2Syn synapse and treats `pre` as the presynaptic
    object.

    Parameters
    ----------
    pre :           NetStim | section
                    The presynaptic object to connect to this neuron.
    where :         int | list of int
                    Node IDs at which to simulate synaptic input.

    Synapse properties:

    syn_tau1 :      int
                    Rise time constant [ms].
    syn_tau2 :      int
                    Decay time constant [ms].
    syn_rev_pot :   int
                    Reversal potential (e) [mV].

    Connection properties:

    cn_thresh :     int
                    Presynaptic membrane potential [mV] at which synaptic
                    event is triggered.
    cn_delay :      int
                    Delay [ms] between presynaptic trigger and postsynaptic
                    event.
    cn_weight :     int
                    Weight variable. This bundles a couple of synaptic
                    properties such as e.g. how much transmitter is released
                    or binding affinity at postsynaptic receptors.

    """
    where = utils.make_iterable(where)

    if not is_NEURON_object(pre):
        raise ValueError(f'Expected NEURON object, got {type(pre)}')

    # Turn section into segment
    if isinstance(pre, neuron.nrn.Section):
        pre = pre()

    # Go over the nodes
    nodes = self.nodes.set_index('node_id')
    for node in nodes.loc[where].itertuples():
        # Generate synapses for the nodes in question
        # Note that we are not reusing existing synapses
        # in case the properties are different
        sec = self.sections[node.sec_ix](node.sec_pos)
        syn = neuron.h.Exp2Syn(sec)
        syn.tau1 = syn_tau1
        syn.tau2 = syn_tau2
        syn.e = syn_rev_pot

        self.synapses[node.Index] = self.synapses.get(node.Index, []) + [syn]

        # Connect spike stimulus and synapse
        if isinstance(pre, neuron.nrn.Segment):
            nc = neuron.h.NetCon(pre._ref_v, syn, sec=pre.sec)
        else:
            nc = neuron.h.NetCon(pre, syn)

        # Set connection parameters
        nc.threshold = cn_thresh
        nc.delay = cn_delay
        nc.weight[0] = cn_weight

        self.stimuli[node.Index] = self.stimuli.get(node.Index, []) + [nc, pre]

Return section(s) for given node(s).

PARAMETER DESCRIPTION
node_ids
    Node IDs.

TYPE: int | list of int

RETURNS DESCRIPTION
section(s) : segment or list of segments

Depends on input.

Source code in navis/interfaces/neuron/comp.py
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
def get_node_section(self, node_ids):
    """Return section(s) for given node(s).

    Parameters
    ----------
    node_ids :  int | list of int
                Node IDs.

    Returns
    -------
    section(s) :    segment or list of segments
                    Depends on input.

    """
    nodes = self.nodes.set_index('node_id')
    if not utils.is_iterable(node_ids):
        n = nodes.loc[node_ids]
        return self.sections[n.sec_ix]
    else:
        segs = []
        for node in nodes.loc[node_ids].itertuples():
            segs.append(self.sections[node.sec_ix])
        return segs

Return segment(s) for given node(s).

PARAMETER DESCRIPTION
node_ids
    Node IDs.

TYPE: int | list of int

RETURNS DESCRIPTION
segment(s) : segment or list of segments

Depends on input.

Source code in navis/interfaces/neuron/comp.py
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
def get_node_segment(self, node_ids):
    """Return segment(s) for given node(s).

    Parameters
    ----------
    node_ids :  int | list of int
                Node IDs.

    Returns
    -------
    segment(s) :    segment or list of segments
                    Depends on input.

    """
    nodes = self.nodes.set_index('node_id')
    if not utils.is_iterable(node_ids):
        n = nodes.loc[node_ids]
        return self.sections[n.sec_ix](n.sec_pos)
    else:
        segs = []
        for node in nodes.loc[node_ids].itertuples():
            segs.append(self.sections[node.sec_ix](node.sec_pos))
        return segs

Add current injection (IClamp) stimulation to model.

PARAMETER DESCRIPTION
where
    Node ID(s) at which to stimulate.

TYPE: int | list of int

start
    Onset (delay) [ms] from beginning of simulation.

TYPE: int DEFAULT: 5

duration
    Duration (dur) [ms] of injection.

TYPE: int DEFAULT: 1

current
    Amount (i) [nA] of injected current.

TYPE: float DEFAULT: 0.1

Source code in navis/interfaces/neuron/comp.py
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
def inject_current_pulse(self, where, start=5,
                         duration=1, current=0.1):
    """Add current injection (IClamp) stimulation to model.

    Parameters
    ----------
    where :     int | list of int
                Node ID(s) at which to stimulate.
    start :     int
                Onset (delay) [ms] from beginning of simulation.
    duration :  int
                Duration (dur) [ms] of injection.
    current :   float
                Amount (i) [nA] of injected current.

    """
    self._add_stimulus('IClamp', where=where, delay=start,
                       dur=duration, amp=current)

Insert biophysical mechanism for model.

PARAMETER DESCRIPTION
mechanism
    Mechanism to insert - e.g. "hh" for Hodgkin-Huxley kinetics.

TYPE: str

subset
    Sections (or indices thereof) to set mechanism for.
    If `None` will add mechanism to all sections.

TYPE: list of sections | list of int DEFAULT: None

**kwargs
    Use to set properties for mechanism.

DEFAULT: {}

Source code in navis/interfaces/neuron/comp.py
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
def insert(self, mechanism, subset=None, **kwargs):
    """Insert biophysical mechanism for model.

    Parameters
    ----------
    mechanism : str
                Mechanism to insert - e.g. "hh" for Hodgkin-Huxley kinetics.
    subset :    list of sections | list of int
                Sections (or indices thereof) to set mechanism for.
                If `None` will add mechanism to all sections.
    **kwargs
                Use to set properties for mechanism.

    """
    if isinstance(subset, type(None)):
        sections = self.sections
    else:
        subset = utils.make_iterable(subset)

        if all([is_section(s) for s in subset]):
            sections = subset
        elif all([isinstance(s, Number) for s in subset]):
            sections = self.sections[subset]
        else:
            raise TypeError('`subset` must be None, a list of sections or '
                            'a list of section indices')

    for sec in np.unique(sections):
        _ = sec.insert(mechanism)
        for seg in sec:
            mech = getattr(seg, mechanism)
            for k, v in kwargs.items():
                setattr(mech, k, v)

Plot results.

PARAMETER DESCRIPTION
axes
    Axes to plot onto. Must have one ax for each recording
    type (mV, spike count, etc) in `self.records`.

TYPE: matplotlib axes DEFAULT: None

RETURNS DESCRIPTION
axes
Source code in navis/interfaces/neuron/comp.py
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
def plot_results(self, axes=None):
    """Plot results.

    Parameters
    ----------
    axes :      matplotlib axes
                Axes to plot onto. Must have one ax for each recording
                type (mV, spike count, etc) in `self.records`.

    Returns
    -------
    axes

    """
    if isinstance(self.t, type(None)) or not len(self.t):
        logger.warning('Looks like the simulation has not yet been run.')
        return
    if not self.records:
        logger.warning('Nothing to plot: no recordings found.')
        return

    if not axes:
        fig, axes = plt.subplots(len(self.records), sharex=True)

    # Make sure that even a single ax is a list
    if not isinstance(axes, (np.ndarray, list)):
        axes = [axes] * len(self.records)

    for t, ax in zip(self.records, axes):
        for i, (k, v) in enumerate(self.records[t].items()):
            if not len(v):
                continue
            v = v.as_numpy()
            # For spikes the vector contains the times
            if t == 'spikes':
                # Calculate spike rate
                bins = np.linspace(0, max(self.t), 10)
                hist, _ = np.histogram(v, bins=bins)
                width = bins[1] - bins[0]
                rate = hist * (1000 / width)
                ax.plot(bins[:-1] + (width / 2), rate, label=k)

                ax.scatter(v, [-i] * len(v), marker='|', s=100)
            else:
                ax.plot(self.t, v, label=k)

        ax.set_xlabel('time [ms]')
        ax.set_ylabel(f'{t}')

        ax.legend()
    return axes

Visualize structure in 3D using matplotlib.

Source code in navis/interfaces/neuron/comp.py
808
809
810
def plot_structure(self):
    """Visualize structure in 3D using matplotlib."""
    _ = neuron.h.PlotShape().plot(plt)

Run the simulation.

Source code in navis/interfaces/neuron/comp.py
812
813
814
815
816
817
818
819
820
def run_simulation(self, duration=25 * ms, v_init=-65 * mV):
    """Run the simulation."""
    # Add recording of time
    global main_t
    main_t = neuron.h.Vector().record(neuron.h._ref_t)

    # This resets the entire model space not just this neuron!
    neuron.h.finitialize(v_init)
    neuron.h.continuerun(duration)

Remove biophysical mechanism from model.

PARAMETER DESCRIPTION
mechanism
    Mechanism to remove - e.g. "hh" for Hodgkin-Huxley kinetics.

TYPE: str

subset
    Sections (or indices thereof) to set mechanism for.
    If `None` will add mechanism to all sections.

TYPE: list of sections | list of int DEFAULT: None

Source code in navis/interfaces/neuron/comp.py
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
def uninsert(self, mechanism, subset=None):
    """Remove biophysical mechanism from model.

    Parameters
    ----------
    mechanism : str
                Mechanism to remove - e.g. "hh" for Hodgkin-Huxley kinetics.
    subset :    list of sections | list of int
                Sections (or indices thereof) to set mechanism for.
                If `None` will add mechanism to all sections.

    """
    if isinstance(subset, type(None)):
        sections = self.sections
    else:
        subset = utils.make_iterable(subset)

        if all([is_section(s) for s in subset]):
            sections = subset
        elif all([isinstance(s, Number) for s in subset]):
            sections = self.sections[subset]
        else:
            raise TypeError('`subset` must be None, a list of sections or '
                            'a list of section indices')

    for sec in np.unique(sections):
        if hasattr(sec, mechanism):
            _ = sec.uninsert(mechanism)

Compartment model of an olfactory projection neuron in Drosophila.

This is a CompartmentModel that uses passive membrane properties from Tobin et al. (2017) as presets:

  • specific axial resistivity (Ra) of 266.1 Ohm / cm
  • specific membrane capacitance (cm) of 0.8 mF / cm**2
  • specific leakage conductance (g) of 1/Rm
  • Rm = specific membran resistance of 20800 Ohm cm**2
  • leakage reverse potential of -60 mV
PARAMETER DESCRIPTION
x
    Neuron to generate model for. Has to be in microns!

TYPE: navis.TreeNeuron

res
    Approximate length [um] of segments. This guarantees that
    no section has any segment that is longer than `res` but for
    small branches (i.e. "sections") the segments might be smaller.
    Lower `res` = more detailed simulation.

TYPE: int DEFAULT: 10

Source code in navis/interfaces/neuron/comp.py
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
class DrosophilaPN(CompartmentModel):
    """Compartment model of an olfactory projection neuron in Drosophila.

    This is a `CompartmentModel` that uses passive membrane properties
    from Tobin et al. (2017) as presets:

    - specific axial resistivity (`Ra`) of 266.1 Ohm / cm
    - specific membrane capacitance (`cm`) of 0.8 mF / cm**2
    - specific leakage conductance (`g`) of 1/Rm
    - Rm = specific membran resistance of 20800 Ohm cm**2
    - leakage reverse potential of -60 mV

    Parameters
    ----------
    x :         navis.TreeNeuron
                Neuron to generate model for. Has to be in microns!
    res :       int
                Approximate length [um] of segments. This guarantees that
                no section has any segment that is longer than `res` but for
                small branches (i.e. "sections") the segments might be smaller.
                Lower `res` = more detailed simulation.

    """

    def __init__(self, x, res=10):
        super().__init__(x, res=res)

        self.Ra = 266.1  # specific axial resistivity in Ohm cm
        self.cm = 0.8    # specific membrane capacitance in mF / cm**2

        # Add passive membran properties
        self.insert('pas',
                    g=1/20800,  # specific leakage conductance = 1/Rm; Rm = specific membran resistance in Ohm cm**2
                    e=-60,      # leakage reverse potential
                    )

A Network of Leaky-Integrate-and-Fire (LIF) point processes.

Examples:

>>> import navis.interfaces.neuron as nrn
>>> N = nrn.PointNetwork()
>>>
Source code in navis/interfaces/neuron/network.py
 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
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
class PointNetwork:
    """A Network of Leaky-Integrate-and-Fire (LIF) point processes.

    Examples
    --------
    >>> import navis.interfaces.neuron as nrn
    >>> N = nrn.PointNetwork()
    >>>
    """

    def __init__(self):
        self._neurons = []
        self._neurons_dict = {}
        self._edges = []
        self._stimuli = []
        self._ids = []
        self._labels = []
        self.idx = NetworkIdIndexer(self)

    def __str__(self):
        return f'{type(self).__name__}<neurons={len(self)},edges={len(self.edges)}>'

    def __repr__(self):
        return self.__str__()

    def __contains__(self, id):
        return id in self._ids

    def __getitem__(self, ix):
        """Get point process with given ID."""
        return np.asarray(self.neurons)[ix]

    def __len__(self):
        return len(self._neurons)

    @property
    def edges(self):
        """Edges between nodes of the network.

        Returns
        -------
        list of tuples

                `[(source_ix, target_ix, weight, NetCon object), ...]`

        """
        return self._edges

    @property
    def neurons(self):
        """Neurons in the network.

        Returns
        -------
        list
                List of `PointNeurons`.

        """
        return self._neurons

    @property
    def ids(self):
        """IDs of neurons in the network."""
        return self._ids

    @property
    def labels(self):
        """Labels of neurons in the network."""
        return self._labels

    @classmethod
    def from_edge_list(cls, edges, model='IntFire1', source_col='source',
                       target_col='target', weight_col='weight', **props):
        """Generate network from edge list.

        Parameters
        ----------
        edges :         pd.DataFrame
                        Edge list. Must contain 'source', 'target' and 'weight'
                        columns.
        model :         "IntFire1" | "IntFire2" | "IntFire4"
                        The model to use for the integrate-and-fire point processes.
        source_col :    str
                        Name of the column with the source IDs.
        target_col :    str
                        Name of the column with the target IDs.
        weight_col :    str
                        Name of the column with the weights. The important thing
                        to note here is that weight is expected to be in the 0-1
                        range with 1 effectively guaranteeing that a presynaptic
                        spike triggers a postsynaptic spike.
        **props
                        Keyword arguments are passed through to `add_neurons`.
                        Use to set e.g. labels, threshold or additional
                        model parameters.

        """
        assert isinstance(edges, pd.DataFrame)
        miss = {source_col, target_col, weight_col} - set(edges.columns)
        if miss:
            raise ValueError(f'edge list is missing columns: {miss}')

        # Instantiate
        net = cls()

        # Add neurons
        ids = np.unique(edges[[source_col, target_col]].values)
        net.add_neurons(ids, model=model, **props)

        # Connect neurons
        for (s, t, w) in edges[[source_col, target_col, weight_col]].values:
            net.connect(s, t, w)

        return net

    def add_neurons(self, ids, model='IntFire1', labels=None,
                    skip_existing=False, **props):
        """Add neurons to network.

        Parameters
        ----------
        ids :           list-like
                        Iterable of IDs for which to create neurons.
        model :         "IntFire1" | "IntFire2" | "IntFire4"
                        The model to use for the integrate-and-fire point processes.
        labels :        str | list | dict, optional
                        Labels for neurons. If str will apply the same label to
                        all neurons. If list, must be same length as `ids`.
                        Dictionary must be ID -> label map.
        skip_existing : bool
                        If True, will skip existing IDs.
        **props
                        Additional parameters used when initializing the point
                        processes. Depends on which model type you are using:
                        e.g. `refrac` and `tau` for "IntFire1".

        """
        assert model in ["IntFire1", "IntFire2", "IntFire4"]
        model = getattr(neuron.h, model)

        ids = utils.make_iterable(ids)

        if isinstance(labels, dict):
            labels = [labels.get(i, 'NA') for i in ids]
        elif labels:
            labels = utils.make_iterable(labels)
        else:
            labels = ids

        if len(labels) != len(ids):
            raise ValueError('Must provide a label for each neuron')

        for i, l in zip(ids, labels):
            if i in self:
                if skip_existing:
                    continue
                else:
                    raise ValueError(f'Neuron with id {i} already exists. '
                                     'Try using `skip_existing=True`.')

            # Create neuron
            n = PointNeuron(id=i, model=model, label=l, **props)

            self._neurons.append(n)
            self._ids.append(i)
            self._labels.append(l)
            self._neurons_dict[i] = n

    def add_background_noise(self, ids, frequency, randomness=.5,
                             independent=True):
        """Add background noise to given neurons.

        Parameters
        ----------
        ids :           hashable  | iterable
                        IDs of neurons to add noise to.
        frequency :     int | float | iterable
                        Frequency [Hz] of background noise. If iterable must
                        match length of `ids`.
        randomness :    float [0-1]
                        Randomness of spike timings.
        independent :   bool
                        If True (default), each neuron will get its own
                        independent noise stimulus.

        """
        self.add_stimulus(ids=ids, start=0, stop=9999999999,
                          frequency=frequency, randomness=randomness,
                          independent=independent)

    def add_stimulus(self, ids, start, frequency, stop=None, duration=None,
                     randomness=.5, independent=True, label=None, weight=1):
        """Add stimulus to given neurons.

        Important
        ---------
        Stimuli are implemented via a NetStim object which provides the
        specified stimulation (i.e. start, stop, frequency). This NetStim is
        then connected to the neuron(s) via a NetCon. The response of the
        neuron to the stimulus depends heavily on the `weight` of that
        connection: too low and you won't elicit any spikes, too high and you
        will produce higher frequencies than expected. The "correct" weight
        depends on the model & parameters you use for your point processes, and
        I highly recommend you check if you get the expected stimulation.

        Parameters
        ----------
        ids :           int | iterable
                        IDs of neurons to add stimulus to.
        start :         int
                        Start time [ms] for the stimulus. Note that
                        the exact start and stop time will vary depending on
                        `randomness`.
        stop/duration : int
                        Either stop time or duration of stimulus [ms]. Must
                        provide one or the other but not both.
        frequency :     int | iterable
                        Frequency [Hz] of background noise . If iterable must
                        match length of `ids`. Values <= 0 are silently skipped.
        randomness :    float [0-1]
                        Randomness of spike timings.
        independent :   bool
                        If True (default), each neuron will get its own
                        independent stimulus.
        label :         str, optional
                        A label to identify the stimulus.
        weight :        float
                        Weight for the connection between the stimulator and
                        the neuron. This really should be 1 to make sure each
                        spike in the stimulus elicits a spike in the target.

        """
        ids = utils.make_iterable(ids)
        if not utils.is_iterable(frequency):
            frequency = [frequency] * len(ids)
        elif not independent:
            raise ValueError('Stimuli/noises must be independent when '
                             'providing individual frequencies')

        if len(frequency) != len(ids):
            raise ValueError('Must provide either a single frequency for all '
                             'neurons or a frequency for each neuron.')

        if (duration and stop) or (not duration and not stop):
            raise ValueError('Must provide either duration or stop (but not both).')

        if stop:
            duration = stop - start
        else:
            stop = start + duration

        if duration <= 0:
            raise ValueError(f'Duration must be greater than zero.')

        if not independent:
            ns = neuron.h.NetStim()
            ns.interval = 1000 / frequency[0]
            ns.noise = randomness
            ns.number = int(duration / 1000 * frequency[0])
            ns.start = start

        for i, f in zip(ids, frequency):
            # Skip frequencies lower than 0
            if f <= 0:
                continue

            interval = 1000 / f
            proc = self.idx[i].process

            if independent:
                ns = neuron.h.NetStim()
                ns.interval = interval
                ns.noise = randomness
                ns.number = int(duration / 1000 * f)
                ns.start = start

            nc = neuron.h.NetCon(ns, proc)
            nc.weight[0] = weight
            nc.delay = 0

            self._stimuli.append(Stimulus(start, stop, f, randomness, i, ns, nc, label))

    def clear_stimuli(self):
        """Clear stimuli."""
        self._stimuli = []

    def connect(self, source, target, weight, delay=5):
        """Connect two neurons.

        Parameters
        ----------
        source :    int | str
                    ID of the source.
        target :    int | str
                    ID of the target
        weight :    float
                    Weight of the edge. The important thing to note here is that
                    the weight is expected to be in the 0-1 range with 1
                    effectively guaranteeing that a presynaptic spike triggers
                    a postsynaptic spike.
        delay :     int
                    Delay in ms between a pre- and a postsynaptic spike.

        """
        # Get the point processes corresponding to source and target
        pre = self.idx[source]
        post = self.idx[target]

        # Connect
        nc = neuron.h.NetCon(pre.process, post.process)
        nc.weight[0] = weight
        nc.delay = delay

        # Keep track
        self._edges.append([source, target, weight, nc])

    def plot_raster(self, subset=None, group=False, stimuli=True, ax=None, label=False,
                    backend='auto', **kwargs):
        """Raster plot of spike timings.

        Parameters
        ----------
        subset :        list-like
                        Subset of IDs to plot. You can also use this to
                        determine the order of appearance.
        ax :            matplotlib axis |  plotly figure, optional
                        Axis/figure to plot onto.
        label :         bool
                        Whether to label individual neurons.
        backend :       "auto" | "plotly" | "matplotlib"
                        Which backend to use. If "auto" will use plotly in
                        Jupyter environments and matplotlib elsewhere.

        """
        if not isinstance(subset, type(None)):
            ids = utils.make_iterable(subset)
            if not len(ids):
                raise ValueError('`ids` must not be empty')
        else:
            ids = self._ids

        # Collect spike timings
        x = []
        y = []
        i = 0
        for id in ids:
            pp = self.idx[id]
            x += list(pp.spk_timings)
            y += [i] * len(pp.spk_timings)
            i += 1

        if not x:
            raise ValueError('No spikes detected.')

        if label:
            ld = dict(zip(self._ids, self._labels))
            labels = [ld[i] for i in ids]
        else:
            labels = None

        # Turn into lines
        x = np.vstack((x, x, [None] * len(x))).T.flatten()
        y = np.array(y)
        y = np.vstack((y, y + .9, [None] * len(y))).T.flatten()

        if backend == 'auto':
            if utils.is_jupyter():
                backend = 'plotly'
            else:
                backend = 'matplotlib'

        if backend == 'plotly':
            return _plot_raster_plotly(x, y, ids, fig=ax, labels=labels, **kwargs)
        elif backend == 'matplotlib':
            ax = _plot_raster_mpl(x, y, ids, ax=ax, labels=labels,
                                  stimuli=self._stimuli if stimuli else None,
                                  **kwargs)
            ax.set_xlim(0, neuron.h.t)
            return ax
        else:
            raise ValueError(f'Unknown backend "{backend}"')

    def plot_traces(self, bin_size=100, subset=None, rolling_window=None,
                    group=False, stimuli=True, ax=None, backend='auto', **kwargs):
        """Plot mean firing rate.

        Parameters
        ----------
        bin_size :          int
                            Size [ms] of the bins over which to average spike
                            frequency.
        subset :            list-like
                            Subset of IDs to plot. You can also use this to
                            determine the order of appearance.
        rolling_window :    int
                            Average firing rates over a given rolling window.
        group :             bool | array
                            Set to True to group traces by label, showing mean
                            firing rate and standard error as envelope. Pass
                            an array of labels for each neuron to group by
                            arbitrary labels.
        ax :                matplotlib axis | plotly figure, optional
                            Axis/figure to plot onto.
        stimuli :           bool
                            Whether to plot stimuli.
        backend :           "auto" | "plotly" | "matplotlib"
                            Which backend to use. If "auto" will use plotly in
                            Jupyter environments and matplotlib elsewhere.

        """
        if not isinstance(subset, type(None)):
            ids = utils.make_iterable(subset)
            if not len(ids):
                raise ValueError('`ids` must not be empty')
        else:
            ids = self._ids

        # Collect spike frequencies
        spks = self.get_spike_counts(bin_size=bin_size, subset=subset,
                                     rolling_window=rolling_window)
        freq = spks * 1000 / bin_size

        if self.labels:
            ld = dict(zip(self._ids, self._labels))
            labels = [f'{i} ({ld[i]})' for i in ids]
        else:
            labels = None

        if backend == 'auto':
            if utils.is_jupyter():
                backend = 'plotly'
            else:
                backend = 'matplotlib'

        if isinstance(group, bool) and group:
            sem = freq.groupby(dict(zip(self._ids, self._labels))).sem()
            freq = freq.groupby(dict(zip(self._ids, self._labels))).mean()
            labels = freq.index.values.tolist()
        elif not isinstance(group, bool):
            sem = freq.groupby(group).sem()
            freq = freq.groupby(group).mean()
            labels = freq.index.values.tolist()
        else:
            sem = None

        if backend == 'plotly':
            return _plot_traces_plotly(freq, fig=ax, labels=labels,
                                       stimuli=self._stimuli if stimuli else None,
                                       env=sem, **kwargs)
        elif backend == 'matplotlib':
            return _plot_traces_mpl(freq, ax=ax, env=sem,
                                    stimuli=self._stimuli if stimuli else None,
                                    **kwargs)
        else:
            raise ValueError(f'Unknown backend "{backend}"')

    def set_labels(self, labels):
        """Set labels for neurons.

        Parameters
        ----------
        labels :    dict | list-like
                    If list, must provide a label for every neuron.

        """
        if isinstance(labels, dict):
            for i, n in enumerate(self.neurons):
                n.label = self._labels[i] = labels.get(n.id, n.id)
        elif utils.is_iterable(labels):
            if len(labels) != len(self.neurons):
                raise ValueError(f'Got {len(labels)} labels for {len(self)} neurons.')
            for i, n in enumerate(self.neurons):
                n.label = self._labels[i] = labels[i]
        else:
            raise TypeError(f'`labels` must be dict or list-like, got "{type(labels)}"')

    def run_simulation(self, duration=25 * ms, v_init=-65 * mV):
        """Run the simulation."""

        # This resets the entire model space not just this neuron!
        neuron.h.finitialize(v_init)
        neuron.h.continuerun(duration)

    def get_spike_counts(self, bin_size=50, subset=None, rolling_window=None,
                         group=False):
        """Get matrix of spike counts.

        Parameters
        ----------
        bin_size :          int | None
                            Size [ms] of the bins over which to count spikes.
                            If None, will simply return total counts.
        rolling_window :    int, optional
                            Average spike counts in a rolling window.
        group :             bool
                            If True, will return the spike counts per unique
                            label.

        Returns
        -------
        pd.DataFrame

        """
        if not isinstance(subset, type(None)):
            ids = utils.make_iterable(subset)
        else:
            ids = self._ids

        end_time = neuron.h.t

        if not end_time:
            raise ValueError('Looks like simulation has not yet been run.')

        if not bin_size:
            bin_size = end_time

        bins = np.arange(0, end_time + bin_size, bin_size)
        counts = np.zeros((len(ids), len(bins) - 1))
        # Collect spike counts
        for i, id in enumerate(ids):
            pp = self.idx[id]
            timings = list(pp.spk_timings)
            if timings:
                hist, _ = np.histogram(timings, bins)
                counts[i, :] = hist

        counts = pd.DataFrame(counts, index=ids, columns=bins[1:])

        if group:
            if not self._labels:
                raise ValueError('Unable to group: Network has no labels.')
            counts = counts.groupby(counts.index.map(dict(zip(self.ids, self._labels)))).sum()

        if rolling_window:
            avg = sliding_window_view(counts, rolling_window, axis=1).mean(axis=2)
            counts.iloc[:, :-(rolling_window - 1)] = avg

        return counts

Edges between nodes of the network.

RETURNS DESCRIPTION
list of tuples

[(source_ix, target_ix, weight, NetCon object), ...]

IDs of neurons in the network.

Labels of neurons in the network.

Neurons in the network.

RETURNS DESCRIPTION
list

List of PointNeurons.

Add background noise to given neurons.

PARAMETER DESCRIPTION
ids
        IDs of neurons to add noise to.

TYPE: hashable | iterable

frequency
        Frequency [Hz] of background noise. If iterable must
        match length of `ids`.

TYPE: int | float | iterable

randomness
        Randomness of spike timings.

TYPE: float [0-1] DEFAULT: 0.5

independent
        If True (default), each neuron will get its own
        independent noise stimulus.

TYPE: bool DEFAULT: True

Source code in navis/interfaces/neuron/network.py
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
def add_background_noise(self, ids, frequency, randomness=.5,
                         independent=True):
    """Add background noise to given neurons.

    Parameters
    ----------
    ids :           hashable  | iterable
                    IDs of neurons to add noise to.
    frequency :     int | float | iterable
                    Frequency [Hz] of background noise. If iterable must
                    match length of `ids`.
    randomness :    float [0-1]
                    Randomness of spike timings.
    independent :   bool
                    If True (default), each neuron will get its own
                    independent noise stimulus.

    """
    self.add_stimulus(ids=ids, start=0, stop=9999999999,
                      frequency=frequency, randomness=randomness,
                      independent=independent)

Add neurons to network.

PARAMETER DESCRIPTION
ids
        Iterable of IDs for which to create neurons.

TYPE: list-like

model
        The model to use for the integrate-and-fire point processes.

TYPE: "IntFire1" | "IntFire2" | "IntFire4" DEFAULT: 'IntFire1'

labels
        Labels for neurons. If str will apply the same label to
        all neurons. If list, must be same length as `ids`.
        Dictionary must be ID -> label map.

TYPE: str | list | dict DEFAULT: None

skip_existing
        If True, will skip existing IDs.

TYPE: bool DEFAULT: False

**props
        Additional parameters used when initializing the point
        processes. Depends on which model type you are using:
        e.g. `refrac` and `tau` for "IntFire1".

DEFAULT: {}

Source code in navis/interfaces/neuron/network.py
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
def add_neurons(self, ids, model='IntFire1', labels=None,
                skip_existing=False, **props):
    """Add neurons to network.

    Parameters
    ----------
    ids :           list-like
                    Iterable of IDs for which to create neurons.
    model :         "IntFire1" | "IntFire2" | "IntFire4"
                    The model to use for the integrate-and-fire point processes.
    labels :        str | list | dict, optional
                    Labels for neurons. If str will apply the same label to
                    all neurons. If list, must be same length as `ids`.
                    Dictionary must be ID -> label map.
    skip_existing : bool
                    If True, will skip existing IDs.
    **props
                    Additional parameters used when initializing the point
                    processes. Depends on which model type you are using:
                    e.g. `refrac` and `tau` for "IntFire1".

    """
    assert model in ["IntFire1", "IntFire2", "IntFire4"]
    model = getattr(neuron.h, model)

    ids = utils.make_iterable(ids)

    if isinstance(labels, dict):
        labels = [labels.get(i, 'NA') for i in ids]
    elif labels:
        labels = utils.make_iterable(labels)
    else:
        labels = ids

    if len(labels) != len(ids):
        raise ValueError('Must provide a label for each neuron')

    for i, l in zip(ids, labels):
        if i in self:
            if skip_existing:
                continue
            else:
                raise ValueError(f'Neuron with id {i} already exists. '
                                 'Try using `skip_existing=True`.')

        # Create neuron
        n = PointNeuron(id=i, model=model, label=l, **props)

        self._neurons.append(n)
        self._ids.append(i)
        self._labels.append(l)
        self._neurons_dict[i] = n

Add stimulus to given neurons.

Important

Stimuli are implemented via a NetStim object which provides the specified stimulation (i.e. start, stop, frequency). This NetStim is then connected to the neuron(s) via a NetCon. The response of the neuron to the stimulus depends heavily on the weight of that connection: too low and you won't elicit any spikes, too high and you will produce higher frequencies than expected. The "correct" weight depends on the model & parameters you use for your point processes, and I highly recommend you check if you get the expected stimulation.

PARAMETER DESCRIPTION
ids
        IDs of neurons to add stimulus to.

TYPE: int | iterable

start
        Start time [ms] for the stimulus. Note that
        the exact start and stop time will vary depending on
        `randomness`.

TYPE: int

stop
        Either stop time or duration of stimulus [ms]. Must
        provide one or the other but not both.

DEFAULT: None

frequency
        Frequency [Hz] of background noise . If iterable must
        match length of `ids`. Values <= 0 are silently skipped.

TYPE: int | iterable

randomness
        Randomness of spike timings.

TYPE: float [0-1] DEFAULT: 0.5

independent
        If True (default), each neuron will get its own
        independent stimulus.

TYPE: bool DEFAULT: True

label
        A label to identify the stimulus.

TYPE: str DEFAULT: None

weight
        Weight for the connection between the stimulator and
        the neuron. This really should be 1 to make sure each
        spike in the stimulus elicits a spike in the target.

TYPE: float DEFAULT: 1

Source code in navis/interfaces/neuron/network.py
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
def add_stimulus(self, ids, start, frequency, stop=None, duration=None,
                 randomness=.5, independent=True, label=None, weight=1):
    """Add stimulus to given neurons.

    Important
    ---------
    Stimuli are implemented via a NetStim object which provides the
    specified stimulation (i.e. start, stop, frequency). This NetStim is
    then connected to the neuron(s) via a NetCon. The response of the
    neuron to the stimulus depends heavily on the `weight` of that
    connection: too low and you won't elicit any spikes, too high and you
    will produce higher frequencies than expected. The "correct" weight
    depends on the model & parameters you use for your point processes, and
    I highly recommend you check if you get the expected stimulation.

    Parameters
    ----------
    ids :           int | iterable
                    IDs of neurons to add stimulus to.
    start :         int
                    Start time [ms] for the stimulus. Note that
                    the exact start and stop time will vary depending on
                    `randomness`.
    stop/duration : int
                    Either stop time or duration of stimulus [ms]. Must
                    provide one or the other but not both.
    frequency :     int | iterable
                    Frequency [Hz] of background noise . If iterable must
                    match length of `ids`. Values <= 0 are silently skipped.
    randomness :    float [0-1]
                    Randomness of spike timings.
    independent :   bool
                    If True (default), each neuron will get its own
                    independent stimulus.
    label :         str, optional
                    A label to identify the stimulus.
    weight :        float
                    Weight for the connection between the stimulator and
                    the neuron. This really should be 1 to make sure each
                    spike in the stimulus elicits a spike in the target.

    """
    ids = utils.make_iterable(ids)
    if not utils.is_iterable(frequency):
        frequency = [frequency] * len(ids)
    elif not independent:
        raise ValueError('Stimuli/noises must be independent when '
                         'providing individual frequencies')

    if len(frequency) != len(ids):
        raise ValueError('Must provide either a single frequency for all '
                         'neurons or a frequency for each neuron.')

    if (duration and stop) or (not duration and not stop):
        raise ValueError('Must provide either duration or stop (but not both).')

    if stop:
        duration = stop - start
    else:
        stop = start + duration

    if duration <= 0:
        raise ValueError(f'Duration must be greater than zero.')

    if not independent:
        ns = neuron.h.NetStim()
        ns.interval = 1000 / frequency[0]
        ns.noise = randomness
        ns.number = int(duration / 1000 * frequency[0])
        ns.start = start

    for i, f in zip(ids, frequency):
        # Skip frequencies lower than 0
        if f <= 0:
            continue

        interval = 1000 / f
        proc = self.idx[i].process

        if independent:
            ns = neuron.h.NetStim()
            ns.interval = interval
            ns.noise = randomness
            ns.number = int(duration / 1000 * f)
            ns.start = start

        nc = neuron.h.NetCon(ns, proc)
        nc.weight[0] = weight
        nc.delay = 0

        self._stimuli.append(Stimulus(start, stop, f, randomness, i, ns, nc, label))

Clear stimuli.

Source code in navis/interfaces/neuron/network.py
342
343
344
def clear_stimuli(self):
    """Clear stimuli."""
    self._stimuli = []

Connect two neurons.

PARAMETER DESCRIPTION
source
    ID of the source.

TYPE: int | str

target
    ID of the target

TYPE: int | str

weight
    Weight of the edge. The important thing to note here is that
    the weight is expected to be in the 0-1 range with 1
    effectively guaranteeing that a presynaptic spike triggers
    a postsynaptic spike.

TYPE: float

delay
    Delay in ms between a pre- and a postsynaptic spike.

TYPE: int DEFAULT: 5

Source code in navis/interfaces/neuron/network.py
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
def connect(self, source, target, weight, delay=5):
    """Connect two neurons.

    Parameters
    ----------
    source :    int | str
                ID of the source.
    target :    int | str
                ID of the target
    weight :    float
                Weight of the edge. The important thing to note here is that
                the weight is expected to be in the 0-1 range with 1
                effectively guaranteeing that a presynaptic spike triggers
                a postsynaptic spike.
    delay :     int
                Delay in ms between a pre- and a postsynaptic spike.

    """
    # Get the point processes corresponding to source and target
    pre = self.idx[source]
    post = self.idx[target]

    # Connect
    nc = neuron.h.NetCon(pre.process, post.process)
    nc.weight[0] = weight
    nc.delay = delay

    # Keep track
    self._edges.append([source, target, weight, nc])

Generate network from edge list.

PARAMETER DESCRIPTION
edges
        Edge list. Must contain 'source', 'target' and 'weight'
        columns.

TYPE: pd.DataFrame

model
        The model to use for the integrate-and-fire point processes.

TYPE: "IntFire1" | "IntFire2" | "IntFire4" DEFAULT: 'IntFire1'

source_col
        Name of the column with the source IDs.

TYPE: str DEFAULT: 'source'

target_col
        Name of the column with the target IDs.

TYPE: str DEFAULT: 'target'

weight_col
        Name of the column with the weights. The important thing
        to note here is that weight is expected to be in the 0-1
        range with 1 effectively guaranteeing that a presynaptic
        spike triggers a postsynaptic spike.

TYPE: str DEFAULT: 'weight'

**props
        Keyword arguments are passed through to `add_neurons`.
        Use to set e.g. labels, threshold or additional
        model parameters.

DEFAULT: {}

Source code in navis/interfaces/neuron/network.py
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
@classmethod
def from_edge_list(cls, edges, model='IntFire1', source_col='source',
                   target_col='target', weight_col='weight', **props):
    """Generate network from edge list.

    Parameters
    ----------
    edges :         pd.DataFrame
                    Edge list. Must contain 'source', 'target' and 'weight'
                    columns.
    model :         "IntFire1" | "IntFire2" | "IntFire4"
                    The model to use for the integrate-and-fire point processes.
    source_col :    str
                    Name of the column with the source IDs.
    target_col :    str
                    Name of the column with the target IDs.
    weight_col :    str
                    Name of the column with the weights. The important thing
                    to note here is that weight is expected to be in the 0-1
                    range with 1 effectively guaranteeing that a presynaptic
                    spike triggers a postsynaptic spike.
    **props
                    Keyword arguments are passed through to `add_neurons`.
                    Use to set e.g. labels, threshold or additional
                    model parameters.

    """
    assert isinstance(edges, pd.DataFrame)
    miss = {source_col, target_col, weight_col} - set(edges.columns)
    if miss:
        raise ValueError(f'edge list is missing columns: {miss}')

    # Instantiate
    net = cls()

    # Add neurons
    ids = np.unique(edges[[source_col, target_col]].values)
    net.add_neurons(ids, model=model, **props)

    # Connect neurons
    for (s, t, w) in edges[[source_col, target_col, weight_col]].values:
        net.connect(s, t, w)

    return net

Get matrix of spike counts.

PARAMETER DESCRIPTION
bin_size
            Size [ms] of the bins over which to count spikes.
            If None, will simply return total counts.

TYPE: int | None DEFAULT: 50

rolling_window
            Average spike counts in a rolling window.

TYPE: int DEFAULT: None

group
            If True, will return the spike counts per unique
            label.

TYPE: bool DEFAULT: False

RETURNS DESCRIPTION
pd.DataFrame
Source code in navis/interfaces/neuron/network.py
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
def get_spike_counts(self, bin_size=50, subset=None, rolling_window=None,
                     group=False):
    """Get matrix of spike counts.

    Parameters
    ----------
    bin_size :          int | None
                        Size [ms] of the bins over which to count spikes.
                        If None, will simply return total counts.
    rolling_window :    int, optional
                        Average spike counts in a rolling window.
    group :             bool
                        If True, will return the spike counts per unique
                        label.

    Returns
    -------
    pd.DataFrame

    """
    if not isinstance(subset, type(None)):
        ids = utils.make_iterable(subset)
    else:
        ids = self._ids

    end_time = neuron.h.t

    if not end_time:
        raise ValueError('Looks like simulation has not yet been run.')

    if not bin_size:
        bin_size = end_time

    bins = np.arange(0, end_time + bin_size, bin_size)
    counts = np.zeros((len(ids), len(bins) - 1))
    # Collect spike counts
    for i, id in enumerate(ids):
        pp = self.idx[id]
        timings = list(pp.spk_timings)
        if timings:
            hist, _ = np.histogram(timings, bins)
            counts[i, :] = hist

    counts = pd.DataFrame(counts, index=ids, columns=bins[1:])

    if group:
        if not self._labels:
            raise ValueError('Unable to group: Network has no labels.')
        counts = counts.groupby(counts.index.map(dict(zip(self.ids, self._labels)))).sum()

    if rolling_window:
        avg = sliding_window_view(counts, rolling_window, axis=1).mean(axis=2)
        counts.iloc[:, :-(rolling_window - 1)] = avg

    return counts

Raster plot of spike timings.

PARAMETER DESCRIPTION
subset
        Subset of IDs to plot. You can also use this to
        determine the order of appearance.

TYPE: list-like DEFAULT: None

ax
        Axis/figure to plot onto.

TYPE: matplotlib axis | plotly figure DEFAULT: None

label
        Whether to label individual neurons.

TYPE: bool DEFAULT: False

backend
        Which backend to use. If "auto" will use plotly in
        Jupyter environments and matplotlib elsewhere.

TYPE: "auto" | "plotly" | "matplotlib" DEFAULT: 'auto'

Source code in navis/interfaces/neuron/network.py
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
def plot_raster(self, subset=None, group=False, stimuli=True, ax=None, label=False,
                backend='auto', **kwargs):
    """Raster plot of spike timings.

    Parameters
    ----------
    subset :        list-like
                    Subset of IDs to plot. You can also use this to
                    determine the order of appearance.
    ax :            matplotlib axis |  plotly figure, optional
                    Axis/figure to plot onto.
    label :         bool
                    Whether to label individual neurons.
    backend :       "auto" | "plotly" | "matplotlib"
                    Which backend to use. If "auto" will use plotly in
                    Jupyter environments and matplotlib elsewhere.

    """
    if not isinstance(subset, type(None)):
        ids = utils.make_iterable(subset)
        if not len(ids):
            raise ValueError('`ids` must not be empty')
    else:
        ids = self._ids

    # Collect spike timings
    x = []
    y = []
    i = 0
    for id in ids:
        pp = self.idx[id]
        x += list(pp.spk_timings)
        y += [i] * len(pp.spk_timings)
        i += 1

    if not x:
        raise ValueError('No spikes detected.')

    if label:
        ld = dict(zip(self._ids, self._labels))
        labels = [ld[i] for i in ids]
    else:
        labels = None

    # Turn into lines
    x = np.vstack((x, x, [None] * len(x))).T.flatten()
    y = np.array(y)
    y = np.vstack((y, y + .9, [None] * len(y))).T.flatten()

    if backend == 'auto':
        if utils.is_jupyter():
            backend = 'plotly'
        else:
            backend = 'matplotlib'

    if backend == 'plotly':
        return _plot_raster_plotly(x, y, ids, fig=ax, labels=labels, **kwargs)
    elif backend == 'matplotlib':
        ax = _plot_raster_mpl(x, y, ids, ax=ax, labels=labels,
                              stimuli=self._stimuli if stimuli else None,
                              **kwargs)
        ax.set_xlim(0, neuron.h.t)
        return ax
    else:
        raise ValueError(f'Unknown backend "{backend}"')

Plot mean firing rate.

PARAMETER DESCRIPTION
bin_size
            Size [ms] of the bins over which to average spike
            frequency.

TYPE: int DEFAULT: 100

subset
            Subset of IDs to plot. You can also use this to
            determine the order of appearance.

TYPE: list-like DEFAULT: None

rolling_window
            Average firing rates over a given rolling window.

TYPE: int DEFAULT: None

group
            Set to True to group traces by label, showing mean
            firing rate and standard error as envelope. Pass
            an array of labels for each neuron to group by
            arbitrary labels.

TYPE: bool | array DEFAULT: False

ax
            Axis/figure to plot onto.

TYPE: matplotlib axis | plotly figure DEFAULT: None

stimuli
            Whether to plot stimuli.

TYPE: bool DEFAULT: True

backend
            Which backend to use. If "auto" will use plotly in
            Jupyter environments and matplotlib elsewhere.

TYPE: "auto" | "plotly" | "matplotlib" DEFAULT: 'auto'

Source code in navis/interfaces/neuron/network.py
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
def plot_traces(self, bin_size=100, subset=None, rolling_window=None,
                group=False, stimuli=True, ax=None, backend='auto', **kwargs):
    """Plot mean firing rate.

    Parameters
    ----------
    bin_size :          int
                        Size [ms] of the bins over which to average spike
                        frequency.
    subset :            list-like
                        Subset of IDs to plot. You can also use this to
                        determine the order of appearance.
    rolling_window :    int
                        Average firing rates over a given rolling window.
    group :             bool | array
                        Set to True to group traces by label, showing mean
                        firing rate and standard error as envelope. Pass
                        an array of labels for each neuron to group by
                        arbitrary labels.
    ax :                matplotlib axis | plotly figure, optional
                        Axis/figure to plot onto.
    stimuli :           bool
                        Whether to plot stimuli.
    backend :           "auto" | "plotly" | "matplotlib"
                        Which backend to use. If "auto" will use plotly in
                        Jupyter environments and matplotlib elsewhere.

    """
    if not isinstance(subset, type(None)):
        ids = utils.make_iterable(subset)
        if not len(ids):
            raise ValueError('`ids` must not be empty')
    else:
        ids = self._ids

    # Collect spike frequencies
    spks = self.get_spike_counts(bin_size=bin_size, subset=subset,
                                 rolling_window=rolling_window)
    freq = spks * 1000 / bin_size

    if self.labels:
        ld = dict(zip(self._ids, self._labels))
        labels = [f'{i} ({ld[i]})' for i in ids]
    else:
        labels = None

    if backend == 'auto':
        if utils.is_jupyter():
            backend = 'plotly'
        else:
            backend = 'matplotlib'

    if isinstance(group, bool) and group:
        sem = freq.groupby(dict(zip(self._ids, self._labels))).sem()
        freq = freq.groupby(dict(zip(self._ids, self._labels))).mean()
        labels = freq.index.values.tolist()
    elif not isinstance(group, bool):
        sem = freq.groupby(group).sem()
        freq = freq.groupby(group).mean()
        labels = freq.index.values.tolist()
    else:
        sem = None

    if backend == 'plotly':
        return _plot_traces_plotly(freq, fig=ax, labels=labels,
                                   stimuli=self._stimuli if stimuli else None,
                                   env=sem, **kwargs)
    elif backend == 'matplotlib':
        return _plot_traces_mpl(freq, ax=ax, env=sem,
                                stimuli=self._stimuli if stimuli else None,
                                **kwargs)
    else:
        raise ValueError(f'Unknown backend "{backend}"')

Run the simulation.

Source code in navis/interfaces/neuron/network.py
536
537
538
539
540
541
def run_simulation(self, duration=25 * ms, v_init=-65 * mV):
    """Run the simulation."""

    # This resets the entire model space not just this neuron!
    neuron.h.finitialize(v_init)
    neuron.h.continuerun(duration)

Set labels for neurons.

PARAMETER DESCRIPTION
labels
    If list, must provide a label for every neuron.

TYPE: dict | list-like

Source code in navis/interfaces/neuron/network.py
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
def set_labels(self, labels):
    """Set labels for neurons.

    Parameters
    ----------
    labels :    dict | list-like
                If list, must provide a label for every neuron.

    """
    if isinstance(labels, dict):
        for i, n in enumerate(self.neurons):
            n.label = self._labels[i] = labels.get(n.id, n.id)
    elif utils.is_iterable(labels):
        if len(labels) != len(self.neurons):
            raise ValueError(f'Got {len(labels)} labels for {len(self)} neurons.')
        for i, n in enumerate(self.neurons):
            n.label = self._labels[i] = labels[i]
    else:
        raise TypeError(f'`labels` must be dict or list-like, got "{type(labels)}"')