1# encoding: utf-8 2""" 3Common implementation of ID, Population, PopulationView and Assembly classes. 4 5These base classes should be sub-classed by the backend-specific classes. 6 7:copyright: Copyright 2006-2021 by the PyNN team, see AUTHORS. 8:license: CeCILL, see LICENSE for details. 9""" 10 11import numpy as np 12import logging 13import operator 14from itertools import chain 15from functools import reduce 16from collections import defaultdict 17from pyNN import random, recording, errors, standardmodels, core, space, descriptions 18from pyNN.models import BaseCellType 19from pyNN.parameters import ParameterSpace, LazyArray, simplify as simplify_parameter_array 20from pyNN.recording import files 21 22 23deprecated = core.deprecated 24logger = logging.getLogger("PyNN") 25 26 27def is_conductance(target_cell): 28 """ 29 Returns True if the target cell uses conductance-based synapses, False if 30 it uses current-based synapses, and None if the synapse-basis cannot be 31 determined. 32 """ 33 if hasattr(target_cell, 'local') and target_cell.local and hasattr(target_cell, 'celltype'): 34 is_conductance = target_cell.celltype.conductance_based 35 else: 36 is_conductance = None 37 return is_conductance 38 39 40class IDMixin(object): 41 """ 42 Instead of storing ids as integers, we store them as ID objects, 43 which allows a syntax like: 44 p[3,4].tau_m = 20.0 45 where p is a Population object. 46 """ 47 # Simulator ID classes should inherit both from the base type of the ID 48 # (e.g., int or long) and from IDMixin. 49 50 def __getattr__(self, name): 51 if name == "parent": 52 raise Exception("parent is not set") 53 elif name == "set": 54 errmsg = "For individual cells, set values using the parameter name directly, " \ 55 "e.g. population[0].tau_m = 20.0, or use 'set' on a population view, " \ 56 "e.g. population[0:1].set(tau_m=20.0)" 57 raise AttributeError(errmsg) 58 try: 59 val = self.get_parameters()[name] 60 except KeyError: 61 raise errors.NonExistentParameterError(name, 62 self.celltype.__class__.__name__, 63 self.celltype.get_parameter_names()) 64 return val 65 66 def __setattr__(self, name, value): 67 if name == "parent": 68 object.__setattr__(self, name, value) 69 elif self.celltype.has_parameter(name): 70 self.set_parameters(**{name: value}) 71 else: 72 object.__setattr__(self, name, value) 73 74 def set_parameters(self, **parameters): 75 """ 76 Set cell parameters, given as a sequence of parameter=value arguments. 77 """ 78 # if some of the parameters are computed from the values of other 79 # parameters, need to get and translate all parameters 80 if self.local: 81 self.as_view().set(**parameters) 82 else: 83 raise errors.NotLocalError( 84 "Cannot set parameters for a cell that does not exist on this node.") 85 86 def get_parameters(self): 87 """Return a dict of all cell parameters.""" 88 if self.local: 89 parameter_names = self.celltype.get_parameter_names() 90 return dict((k, v) for k, v in zip(parameter_names, self.as_view().get(parameter_names))) 91 else: 92 raise errors.NotLocalError( 93 "Cannot obtain parameters for a cell that does not exist on this node.") 94 95 @property 96 def celltype(self): 97 return self.parent.celltype 98 99 @property 100 def is_standard_cell(self): 101 return isinstance(self.celltype, standardmodels.StandardCellType) 102 103 def _set_position(self, pos): 104 """ 105 Set the cell position in 3D space. 106 107 Cell positions are stored in an array in the parent Population. 108 """ 109 assert isinstance(pos, (tuple, np.ndarray)) 110 assert len(pos) == 3 111 self.parent._set_cell_position(self, pos) 112 113 def _get_position(self): 114 """ 115 Return the cell position in 3D space. 116 117 Cell positions are stored in an array in the parent Population, if any, 118 or within the ID object otherwise. Positions are generated the first 119 time they are requested and then cached. 120 """ 121 return self.parent._get_cell_position(self) 122 123 position = property(_get_position, _set_position) 124 125 @property 126 def local(self): 127 return self.parent.is_local(self) 128 129 def inject(self, current_source): 130 """Inject current from a current source object into the cell.""" 131 current_source.inject_into([self]) 132 133 def get_initial_value(self, variable): 134 """Get the initial value of a state variable of the cell.""" 135 return self.parent._get_cell_initial_value(self, variable) 136 137 def set_initial_value(self, variable, value): 138 """Set the initial value of a state variable of the cell.""" 139 self.parent._set_cell_initial_value(self, variable, value) 140 141 def as_view(self): 142 """Return a PopulationView containing just this cell.""" 143 index = self.parent.id_to_index(self) 144 return self.parent[index:index + 1] 145 146 147class BasePopulation(object): 148 _record_filter = None 149 150 def __getitem__(self, index): 151 """ 152 Return either a single cell (ID object) from the Population, if `index` 153 is an integer, or a subset of the cells (PopulationView object), if 154 `index` is a slice or array. 155 156 Note that __getitem__ is called when using [] access, e.g. 157 p = Population(...) 158 p[2] is equivalent to p.__getitem__(2). 159 p[3:6] is equivalent to p.__getitem__(slice(3, 6)) 160 """ 161 if isinstance(index, (int, np.integer)): 162 return self.all_cells[index] 163 elif isinstance(index, (slice, list, np.ndarray)): 164 return self._get_view(index) 165 elif isinstance(index, tuple): 166 return self._get_view(list(index)) 167 else: 168 raise TypeError( 169 "indices must be integers, slices, lists, arrays or tuples, not %s" % type(index).__name__) 170 171 def __len__(self): 172 """Return the total number of cells in the population (all nodes).""" 173 return self.size 174 175 @property 176 def local_size(self): 177 """Return the number of cells in the population on the local MPI node""" 178 return len(self.local_cells) # would self._mask_local.sum() be faster? 179 180 def __iter__(self): 181 """Iterator over cell ids on the local node.""" 182 return iter(self.local_cells) 183 184 @property 185 def conductance_based(self): 186 """ 187 Indicates whether the post-synaptic response is modelled as a change 188 in conductance or a change in current. 189 """ 190 return self.celltype.conductance_based 191 192 @property 193 def receptor_types(self): 194 return self.celltype.receptor_types 195 196 def is_local(self, id): 197 """ 198 Indicates whether the cell with the given ID exists on the local MPI node. 199 """ 200 assert id.parent is self 201 index = self.id_to_index(id) 202 return self._mask_local[index] 203 204 def all(self): 205 """Iterator over cell ids on all MPI nodes.""" 206 return iter(self.all_cells) 207 208 def __add__(self, other): 209 """ 210 A Population/PopulationView can be added to another Population, 211 PopulationView or Assembly, returning an Assembly. 212 """ 213 assert isinstance(other, BasePopulation) 214 return self._assembly_class(self, other) 215 216 def _get_cell_position(self, id): 217 index = self.id_to_index(id) 218 return self.positions[:, index] 219 220 def _set_cell_position(self, id, pos): 221 index = self.id_to_index(id) 222 self.positions[:, index] = pos 223 224 @property 225 def position_generator(self): # "generator" is a misleading name, has no yield statement 226 def gen(i): 227 return self.positions.T[i] 228 return gen 229 230 def _get_cell_initial_value(self, id, variable): 231 if variable in self.initial_values: 232 assert isinstance(self.initial_values[variable], LazyArray) 233 index = self.id_to_local_index(id) 234 return self.initial_values[variable][index] 235 else: 236 logger.warning( 237 "Variable '{}' is not in initial values, returning 0.0".format(variable)) 238 return 0.0 239 240 def _set_cell_initial_value(self, id, variable, value): 241 assert isinstance(self.initial_values[variable], LazyArray) 242 index = self.id_to_local_index(id) 243 self.initial_values[variable][index] = value 244 245 def nearest(self, position): 246 """Return the neuron closest to the specified position.""" 247 # doesn't always work correctly if a position is equidistant between 248 # two neurons, i.e. 0.5 should be rounded up, but it isn't always. 249 # also doesn't take account of periodic boundary conditions 250 pos = np.array([position] * self.positions.shape[1]).transpose() 251 dist_arr = (self.positions - pos)**2 252 distances = dist_arr.sum(axis=0) 253 nearest = distances.argmin() 254 return self[nearest] 255 256 def sample(self, n, rng=None): 257 """ 258 Randomly sample `n` cells from the Population, and return a 259 PopulationView object. 260 """ 261 assert isinstance(n, int) 262 if not rng: 263 rng = random.NumpyRNG() 264 indices = rng.permutation(np.arange(len(self), dtype=np.int))[0:n] 265 logger.debug("The %d cells selected have indices %s" % (n, indices)) 266 logger.debug("%s.sample(%s)", self.label, n) 267 return self._get_view(indices) 268 269 def get(self, parameter_names, gather=False, simplify=True): 270 """ 271 Get the values of the given parameters for every local cell in the 272 population, or, if gather=True, for all cells in the population. 273 274 Values will be expressed in the standard PyNN units (i.e. millivolts, 275 nanoamps, milliseconds, microsiemens, nanofarads, event per second). 276 """ 277 # if all the cells have the same value for a parameter, should 278 # we return just the number, rather than an array? 279 if isinstance(parameter_names, str): 280 parameter_names = (parameter_names,) 281 return_list = False 282 else: 283 return_list = True 284 if isinstance(self.celltype, standardmodels.StandardCellType): 285 if any(name in self.celltype.computed_parameters() for name in parameter_names): 286 native_names = self.celltype.get_native_names() # need all parameters in order to calculate values 287 else: 288 native_names = self.celltype.get_native_names(*parameter_names) 289 native_parameter_space = self._get_parameters(*native_names) 290 parameter_space = self.celltype.reverse_translate(native_parameter_space) 291 else: 292 parameter_space = self._get_parameters(*parameter_names) 293 # what if parameter space is homogeneous on some nodes but not on others? 294 parameter_space.evaluate(simplify=simplify) 295 # this also causes problems if the population size matches the number of MPI nodes 296 parameters = dict(parameter_space.items()) 297 if gather == True and self._simulator.state.num_processes > 1: 298 # seems inefficient to do it in a loop - should do as single operation 299 for name in parameter_names: 300 values = parameters[name] 301 if isinstance(values, np.ndarray): 302 all_values = {self._simulator.state.mpi_rank: values.tolist()} 303 local_indices = np.arange(self.size,)[self._mask_local].tolist() 304 all_indices = {self._simulator.state.mpi_rank: local_indices} 305 all_values = recording.gather_dict(all_values) 306 all_indices = recording.gather_dict(all_indices) 307 if self._simulator.state.mpi_rank == 0: 308 values = reduce(operator.add, all_values.values()) 309 indices = reduce(operator.add, all_indices.values()) 310 idx = np.argsort(indices) 311 values = np.array(values)[idx] 312 parameters[name] = values 313 try: 314 values = [parameters[name] for name in parameter_names] 315 except KeyError as err: 316 raise errors.NonExistentParameterError("%s. Valid parameters for %s are: %s" % ( 317 err, self.celltype, self.celltype.get_parameter_names())) 318 if return_list: 319 return values 320 else: 321 assert len(parameter_names) == 1 322 return values[0] 323 324 def set(self, **parameters): 325 """ 326 Set one or more parameters for every cell in the population. 327 328 Values passed to set() may be: 329 (1) single values 330 (2) RandomDistribution objects 331 (3) lists/arrays of values of the same size as the population 332 (4) mapping functions, where a mapping function accepts a single 333 argument (the cell index) and returns a single value. 334 335 Here, a "single value" may be either a single number or a list/array of 336 numbers (e.g. for spike times). Values should be expressed in the 337 standard PyNN units (i.e. millivolts, nanoamps, milliseconds, 338 microsiemens, nanofarads, event per second). 339 340 Examples:: 341 342 p.set(tau_m=20.0, v_rest=-65). 343 p.set(spike_times=[0.3, 0.7, 0.9, 1.4]) 344 p.set(cm=rand_distr, tau_m=lambda i: 10 + i/10.0) 345 """ 346 # TODO: add example using of function of (x,y,z) and Population.position_generator 347 if self.local_size > 0: 348 if (isinstance(self.celltype, standardmodels.StandardCellType) 349 and any(name in self.celltype.computed_parameters() for name in parameters) 350 and not isinstance(self.celltype, standardmodels.cells.SpikeSourceArray)): 351 # the last condition above is a bit of hack to avoid calling expand() unecessarily 352 # need to get existing parameter space of models so we can perform calculations 353 native_names = self.celltype.get_native_names() 354 parameter_space = self.celltype.reverse_translate( 355 self._get_parameters(*native_names)) 356 if self.local_size != self.size: 357 parameter_space.expand((self.size,), self._mask_local) 358 parameter_space.update(**parameters) 359 else: 360 parameter_space = ParameterSpace(parameters, 361 self.celltype.get_schema(), 362 (self.size,), 363 self.celltype.__class__) 364 if isinstance(self.celltype, standardmodels.StandardCellType): 365 parameter_space = self.celltype.translate(parameter_space) 366 assert parameter_space.shape == (self.size,), "{} != {}".format( 367 parameter_space.shape, self.size) 368 self._set_parameters(parameter_space) 369 370 @deprecated("set(parametername=value_array)") 371 def tset(self, parametername, value_array): 372 """ 373 'Topographic' set. Set the value of parametername to the values in 374 value_array, which must have the same dimensions as the Population. 375 """ 376 self.set(**{parametername: value_array}) 377 378 @deprecated("set(parametername=rand_distr)") 379 def rset(self, parametername, rand_distr): 380 """ 381 'Random' set. Set the value of parametername to a value taken from 382 rand_distr, which should be a RandomDistribution object. 383 """ 384 # Note that we generate enough random numbers for all cells on all nodes 385 # but use only those relevant to this node. This ensures that the 386 # sequence of random numbers does not depend on the number of nodes, 387 # provided that the same rng with the same seed is used on each node. 388 self.set(**{parametername: rand_distr}) 389 390 def initialize(self, **initial_values): 391 """ 392 Set initial values of state variables, e.g. the membrane potential. 393 394 Values passed to initialize() may be: 395 (1) single numeric values (all neurons set to the same value) 396 (2) RandomDistribution objects 397 (3) lists/arrays of numbers of the same size as the population 398 (4) mapping functions, where a mapping function accepts a single 399 argument (the cell index) and returns a single number. 400 401 Values should be expressed in the standard PyNN units (i.e. millivolts, 402 nanoamps, milliseconds, microsiemens, nanofarads, event per second). 403 404 Examples:: 405 406 p.initialize(v=-70.0) 407 p.initialize(v=rand_distr, gsyn_exc=0.0) 408 p.initialize(v=lambda i: -65 + i/10.0) 409 """ 410 for variable, value in initial_values.items(): 411 logger.debug("In Population '%s', initialising %s to %s" % 412 (self.label, variable, value)) 413 initial_value = LazyArray(value, shape=(self.size,), dtype=float) 414 self._set_initial_value_array(variable, initial_value) 415 self.initial_values[variable] = initial_value 416 417 def find_units(self, variable): 418 """ 419 Returns units of the specified variable or parameter, as a string. 420 Works for all the recordable variables and neuron parameters of all standard models. 421 """ 422 return self.celltype.units[variable] 423 424 def annotate(self, **annotations): 425 self.annotations.update(annotations) 426 427 def can_record(self, variable): 428 """Determine whether `variable` can be recorded from this population.""" 429 return self.celltype.can_record(variable) 430 431 @property 432 def injectable(self): 433 return self.celltype.injectable 434 435 def record(self, variables, to_file=None, sampling_interval=None): 436 """ 437 Record the specified variable or variables for all cells in the 438 Population or view. 439 440 `variables` may be either a single variable name or a list of variable 441 names. For a given celltype class, `celltype.recordable` contains a list of 442 variables that can be recorded for that celltype. 443 444 If specified, `to_file` should be either a filename or a Neo IO instance and `write_data()` 445 will be automatically called when `end()` is called. 446 447 `sampling_interval` should be a value in milliseconds, and an integer 448 multiple of the simulation timestep. 449 """ 450 if variables is None: # reset the list of things to record 451 # note that if record(None) is called on a view of a population 452 # recording will be reset for the entire population, not just the view 453 self.recorder.reset() 454 else: 455 logger.debug("%s.record('%s')", self.label, variables) 456 if self._record_filter is None: 457 self.recorder.record(variables, self.all_cells, sampling_interval) 458 else: 459 self.recorder.record(variables, self._record_filter, sampling_interval) 460 if isinstance(to_file, str): 461 self.recorder.file = to_file 462 self._simulator.state.write_on_end.append((self, variables, self.recorder.file)) 463 464 @deprecated("record('v')") 465 def record_v(self, to_file=True): 466 """ 467 Record the membrane potential for all cells in the Population. 468 """ 469 self.record('v', to_file) 470 471 @deprecated("record(['gsyn_exc', 'gsyn_inh'])") 472 def record_gsyn(self, to_file=True): 473 """ 474 Record synaptic conductances for all cells in the Population. 475 """ 476 self.record(['gsyn_exc', 'gsyn_inh'], to_file) 477 478 def write_data(self, io, variables='all', gather=True, clear=False, annotations=None): 479 """ 480 Write recorded data to file, using one of the file formats supported by 481 Neo. 482 483 `io`: 484 a Neo IO instance 485 `variables`: 486 either a single variable name or a list of variable names. 487 Variables must have been previously recorded, otherwise an 488 Exception will be raised. 489 490 For parallel simulators, if `gather` is True, all data will be gathered 491 to the master node and a single output file created there. Otherwise, a 492 file will be written on each node, containing only data from the cells 493 simulated on that node. 494 495 If `clear` is True, recorded data will be deleted from the `Population`. 496 497 `annotations` should be a dict containing simple data types such as 498 numbers and strings. The contents will be written into the output data 499 file as metadata. 500 """ 501 logger.debug("Population %s is writing %s to %s [gather=%s, clear=%s]" % ( 502 self.label, variables, io, gather, clear)) 503 self.recorder.write(variables, io, gather, self._record_filter, clear=clear, 504 annotations=annotations) 505 506 def get_data(self, variables='all', gather=True, clear=False): 507 """ 508 Return a Neo `Block` containing the data (spikes, state variables) 509 recorded from the Population. 510 511 `variables` - either a single variable name or a list of variable names 512 Variables must have been previously recorded, otherwise an 513 Exception will be raised. 514 515 For parallel simulators, if `gather` is True, all data will be gathered 516 to all nodes and the Neo `Block` will contain data from all nodes. 517 Otherwise, the Neo `Block` will contain only data from the cells 518 simulated on the local node. 519 520 If `clear` is True, recorded data will be deleted from the `Population`. 521 """ 522 return self.recorder.get(variables, gather, self._record_filter, clear) 523 524 @deprecated("write_data(file, 'spikes')") 525 def printSpikes(self, file, gather=True, compatible_output=True): 526 self.write_data(file, 'spikes', gather) 527 528 @deprecated("get_data('spikes')") 529 def getSpikes(self, gather=True, compatible_output=True): 530 return self.get_data('spikes', gather) 531 532 @deprecated("write_data(file, 'v')") 533 def print_v(self, file, gather=True, compatible_output=True): 534 self.write_data(file, 'v', gather) 535 536 @deprecated("get_data('v')") 537 def get_v(self, gather=True, compatible_output=True): 538 return self.get_data('v', gather) 539 540 @deprecated("write_data(file, ['gsyn_exc', 'gsyn_inh'])") 541 def print_gsyn(self, file, gather=True, compatible_output=True): 542 self.write_data(file, ['gsyn_exc', 'gsyn_inh'], gather) 543 544 @deprecated("get_data(['gsyn_exc', 'gsyn_inh'])") 545 def get_gsyn(self, gather=True, compatible_output=True): 546 return self.get_data(['gsyn_exc', 'gsyn_inh'], gather) 547 548 def get_spike_counts(self, gather=True): 549 """ 550 Returns a dict containing the number of spikes for each neuron. 551 552 The dict keys are neuron IDs, not indices. 553 """ 554 # arguably, we should use indices 555 return self.recorder.count('spikes', gather, self._record_filter) 556 557 @deprecated("mean_spike_count()") 558 def meanSpikeCount(self, gather=True): 559 return self.mean_spike_count(gather) 560 561 def mean_spike_count(self, gather=True): 562 """ 563 Returns the mean number of spikes per neuron. 564 """ 565 spike_counts = self.get_spike_counts(gather) 566 total_spikes = sum(spike_counts.values()) 567 if self._simulator.state.mpi_rank == 0 or not gather: # should maybe use allgather, and get the numbers on all nodes 568 if len(spike_counts) > 0: 569 return float(total_spikes) / len(spike_counts) 570 else: 571 return 0 572 else: 573 return np.nan 574 575 def inject(self, current_source): 576 """ 577 Connect a current source to all cells in the Population. 578 """ 579 if not self.celltype.injectable: 580 raise TypeError("Can't inject current into a spike source.") 581 current_source.inject_into(self) 582 583 # name should be consistent with saving/writing data, i.e. save_data() and save_positions() or write_data() and write_positions() 584 def save_positions(self, file): 585 """ 586 Save positions to file. The output format is ``index x y z`` 587 """ 588 if isinstance(file, str): 589 file = recording.files.StandardTextFile(file, mode='w') 590 cells = self.all_cells 591 result = np.empty((len(cells), 4)) 592 result[:, 0] = np.array([self.id_to_index(id) for id in cells]) 593 result[:, 1:4] = self.positions.T 594 if self._simulator.state.mpi_rank == 0: 595 file.write(result, {'population': self.label}) 596 file.close() 597 598 599class Population(BasePopulation): 600 """ 601 A group of neurons all of the same type. "Population" is used as a generic 602 term intended to include layers, columns, nuclei, etc., of cells. 603 604 Arguments: 605 `size`: 606 number of cells in the Population. For backwards-compatibility, 607 `size` may also be a tuple giving the dimensions of a grid, 608 e.g. ``size=(10,10)`` is equivalent to ``size=100`` with ``structure=Grid2D()``. 609 610 `cellclass`: 611 a cell type (a class inheriting from :class:`pyNN.models.BaseCellType`). 612 613 `cellparams`: 614 a dict, or other mapping, containing parameters, which is passed to 615 the neuron model constructor. 616 617 `structure`: 618 a :class:`pyNN.space.Structure` instance, used to specify the 619 positions of neurons in space. 620 621 `initial_values`: 622 a dict, or other mapping, containing initial values for the neuron 623 state variables. 624 625 `label`: 626 a name for the population. One will be auto-generated if this is not 627 supplied. 628 """ 629 _nPop = 0 630 631 def __init__(self, size, cellclass, cellparams=None, structure=None, 632 initial_values={}, label=None): 633 """ 634 Create a population of neurons all of the same type. 635 """ 636 if not hasattr(self, "_simulator"): 637 errmsg = "`common.Population` should not be instantiated directly. " \ 638 "You should import Population from a PyNN backend module, " \ 639 "e.g. pyNN.nest or pyNN.neuron" 640 raise Exception(errmsg) 641 if not isinstance(size, (int, np.integer)): # also allow a single integer, for a 1D population 642 assert isinstance( 643 size, tuple), "`size` must be an integer or a tuple of ints. You have supplied a %s" % type(size) 644 # check the things inside are ints 645 for e in size: 646 assert isinstance( 647 e, int), "`size` must be an integer or a tuple of ints. Element '%s' is not an int" % str(e) 648 649 assert structure is None, "If you specify `size` as a tuple you may not specify structure." 650 if len(size) == 1: 651 structure = space.Line() 652 elif len(size) == 2: 653 nx, ny = size 654 structure = space.Grid2D(nx / float(ny)) 655 elif len(size) == 3: 656 nx, ny, nz = size 657 structure = space.Grid3D(nx / float(ny), nx / float(nz)) 658 else: 659 raise Exception( 660 "A maximum of 3 dimensions is allowed. What do you think this is, string theory?") 661 # NEST doesn't like np.int, so to be safe we cast to Python int 662 size = int(reduce(operator.mul, size)) 663 self.size = size 664 self.label = label or 'population%d' % Population._nPop 665 self._structure = structure or space.Line() 666 self._positions = None 667 self._is_sorted = True 668 if isinstance(cellclass, BaseCellType): 669 self.celltype = cellclass 670 assert cellparams is None # cellparams being retained for backwards compatibility, but use is deprecated 671 elif issubclass(cellclass, BaseCellType): 672 self.celltype = cellclass(**cellparams) 673 # emit deprecation warning 674 else: 675 raise TypeError( 676 "cellclass must be an instance or subclass of BaseCellType, not a %s" % type(cellclass)) 677 self.annotations = {} 678 self.recorder = self._recorder_class(self) 679 # Build the arrays of cell ids 680 # Cells on the local node are represented as ID objects, other cells by integers 681 # All are stored in a single numpy array for easy lookup by address 682 # The local cells are also stored in a list, for easy iteration 683 self._create_cells() 684 self.first_id = self.all_cells[0] 685 self.last_id = self.all_cells[-1] 686 self.initial_values = {} 687 all_initial_values = self.celltype.default_initial_values.copy() 688 all_initial_values.update(initial_values) 689 self.initialize(**all_initial_values) 690 Population._nPop += 1 691 692 def __repr__(self): 693 return "Population(%d, %r, structure=%r, label=%r)" % (self.size, self.celltype, self.structure, self.label) 694 695 @property 696 def local_cells(self): 697 """ 698 An array containing cell ids for the local node. 699 """ 700 return self.all_cells[self._mask_local] 701 702 def id_to_index(self, id): 703 """ 704 Given the ID(s) of cell(s) in the Population, return its (their) index 705 (order in the Population). 706 707 >>> assert p.id_to_index(p[5]) == 5 708 """ 709 if not np.iterable(id): 710 if not self.first_id <= id <= self.last_id: 711 raise ValueError("id should be in the range [%d,%d], actually %d" % ( 712 self.first_id, self.last_id, id)) 713 return int(id - self.first_id) # this assumes ids are consecutive 714 else: 715 if isinstance(id, PopulationView): 716 id = id.all_cells 717 id = np.array(id) 718 if (self.first_id > id.min()) or (self.last_id < id.max()): 719 raise ValueError("ids should be in the range [%d,%d], actually [%d, %d]" % ( 720 self.first_id, self.last_id, id.min(), id.max())) 721 return (id - self.first_id).astype(np.int) # this assumes ids are consecutive 722 723 def id_to_local_index(self, id): 724 """ 725 Given the ID(s) of cell(s) in the Population, return its (their) index 726 (order in the Population), counting only cells on the local MPI node. 727 """ 728 if self._simulator.state.num_processes > 1: 729 return self.local_cells.tolist().index(id) # probably very slow 730 # return np.nonzero(self.local_cells == id)[0][0] # possibly faster? 731 # another idea - get global index, use idx-sum(mask_local[:idx])? 732 else: 733 return self.id_to_index(id) 734 735 def _get_structure(self): 736 """The spatial structure of the Population.""" 737 return self._structure 738 739 def _set_structure(self, structure): 740 assert isinstance(structure, space.BaseStructure) 741 if self._structure is None or structure != self._structure: 742 self._positions = None # setting a new structure invalidates previously calculated positions 743 self._structure = structure 744 structure = property(fget=_get_structure, fset=_set_structure) 745 # arguably structure should be read-only, i.e. it is not possible to change it after Population creation 746 747 def _get_positions(self): 748 """ 749 Try to return self._positions. If it does not exist, create it and then 750 return it. 751 """ 752 if self._positions is None: 753 self._positions = self.structure.generate_positions(self.size) 754 assert self._positions.shape == (3, self.size) 755 return self._positions 756 757 def _set_positions(self, pos_array): 758 assert isinstance(pos_array, np.ndarray) 759 assert pos_array.shape == (3, self.size), "%s != %s" % (pos_array.shape, (3, self.size)) 760 self._positions = pos_array.copy() # take a copy in case pos_array is changed later 761 self._structure = None # explicitly setting positions destroys any previous structure 762 763 positions = property(_get_positions, _set_positions, 764 doc="""A 3xN array (where N is the number of neurons in the Population) 765 giving the x,y,z coordinates of all the neurons (soma, in the 766 case of non-point models).""") 767 768 def describe(self, template='population_default.txt', engine='default'): 769 """ 770 Returns a human-readable description of the population. 771 772 The output may be customized by specifying a different template 773 together with an associated template engine (see :mod:`pyNN.descriptions`). 774 775 If template is None, then a dictionary containing the template context 776 will be returned. 777 """ 778 context = { 779 "label": self.label, 780 "celltype": self.celltype.describe(template=None), 781 "structure": None, 782 "size": self.size, 783 "size_local": len(self.local_cells), 784 "first_id": self.first_id, 785 "last_id": self.last_id, 786 } 787 context.update(self.annotations) 788 if len(self.local_cells) > 0: 789 first_id = self.local_cells[0] 790 context.update({ 791 "local_first_id": first_id, 792 "cell_parameters": {} # first_id.get_parameters(), 793 }) 794 if self.structure: 795 context["structure"] = self.structure.describe(template=None) 796 return descriptions.render(engine, template, context) 797 798 799class PopulationView(BasePopulation): 800 """ 801 A view of a subset of neurons within a Population. 802 803 In most ways, Populations and PopulationViews have the same behaviour, i.e. 804 they can be recorded, connected with Projections, etc. It should be noted 805 that any changes to neurons in a PopulationView will be reflected in the 806 parent Population and vice versa. 807 808 It is possible to have views of views. 809 810 Arguments: 811 selector: 812 a slice or numpy mask array. The mask array should either be a 813 boolean array of the same size as the parent, or an integer array 814 containing cell indices, i.e. if p.size == 5:: 815 816 PopulationView(p, array([False, False, True, False, True])) 817 PopulationView(p, array([2,4])) 818 PopulationView(p, slice(2,5,2)) 819 820 will all create the same view. 821 """ 822 823 def __init__(self, parent, selector, label=None): 824 """ 825 Create a view of a subset of neurons within a parent Population or 826 PopulationView. 827 """ 828 if not hasattr(self, "_simulator"): 829 errmsg = "`common.PopulationView` should not be instantiated directly. " \ 830 "You should import PopulationView from a PyNN backend module, " \ 831 "e.g. pyNN.nest or pyNN.neuron" 832 raise Exception(errmsg) 833 self.parent = parent 834 self.mask = selector # later we can have fancier selectors, for now we just have numpy masks 835 # maybe just redefine __getattr__ instead of the following... 836 self.celltype = self.parent.celltype 837 # If the mask is a slice, IDs will be consecutives without duplication. 838 # If not, then we need to remove duplicated IDs 839 if not isinstance(self.mask, slice): 840 if isinstance(self.mask, list): 841 self.mask = np.array(self.mask) 842 if self.mask.dtype is np.dtype('bool'): 843 if len(self.mask) != len(self.parent): 844 raise Exception("Boolean masks should have the size of Parent Population") 845 self.mask = np.arange(len(self.parent))[self.mask] 846 else: 847 if len(np.unique(self.mask)) != len(self.mask): 848 logging.warning( 849 "PopulationView can contain only once each ID, duplicated IDs are removed") 850 self.mask = np.unique(self.mask) 851 self.mask.sort() # needed by NEST. Maybe emit a warning or exception if mask is not already ordered? 852 self.all_cells = self.parent.all_cells[self.mask] 853 idx = np.argsort(self.all_cells) 854 self._is_sorted = np.all(idx == np.arange(len(self.all_cells))) 855 self.size = len(self.all_cells) 856 self.label = label or "view of '%s' with size %s" % (parent.label, self.size) 857 self._mask_local = self.parent._mask_local[self.mask] 858 self.local_cells = self.all_cells[self._mask_local] 859 # only works if we assume all_cells is sorted, otherwise could use min() 860 self.first_id = np.min(self.all_cells) 861 self.last_id = np.max(self.all_cells) 862 self.annotations = {} 863 self.recorder = self.parent.recorder 864 self._record_filter = self.all_cells 865 866 def __repr__(self): 867 return "PopulationView(parent=%r, selector=%r, label=%r)" % (self.parent, self.mask, self.label) 868 869 @property 870 def initial_values(self): 871 # this is going to be complex - if we keep initial_values as a dict, 872 # need to return a dict-like object that takes account of self.mask 873 raise NotImplementedError 874 875 @property 876 def structure(self): 877 """The spatial structure of the parent Population.""" 878 return self.parent.structure 879 # should we allow setting structure for a PopulationView? Maybe if the 880 # parent has some kind of CompositeStructure? 881 882 @property 883 def positions(self): 884 # make positions N,3 instead of 3,N to avoid all this transposing? 885 return self.parent.positions.T[self.mask].T 886 887 def id_to_index(self, id): 888 """ 889 Given the ID(s) of cell(s) in the PopulationView, return its/their 890 index/indices (order in the PopulationView). 891 892 >>> assert pv.id_to_index(pv[3]) == 3 893 """ 894 if not np.iterable(id): 895 if self._is_sorted: 896 if id not in self.all_cells: 897 raise IndexError("ID %s not present in the View" % id) 898 return np.searchsorted(self.all_cells, id) 899 else: 900 result = np.where(self.all_cells == id)[0] 901 if len(result) == 0: 902 raise IndexError("ID %s not present in the View" % id) 903 else: 904 return result 905 else: 906 if self._is_sorted: 907 return np.searchsorted(self.all_cells, id) 908 else: 909 result = np.array([], dtype=np.int) 910 for item in id: 911 data = np.where(self.all_cells == item)[0] 912 if len(data) == 0: 913 raise IndexError("ID %s not present in the View" % item) 914 elif len(data) > 1: 915 raise Exception("ID %s is duplicated in the View" % item) 916 else: 917 result = np.append(result, data) 918 return result 919 920 @property 921 def grandparent(self): 922 """ 923 Returns the parent Population at the root of the tree (since the 924 immediate parent may itself be a PopulationView). 925 926 The name "grandparent" is of course a little misleading, as it could 927 be just the parent, or the great, great, great, ..., grandparent. 928 """ 929 if hasattr(self.parent, "parent"): 930 return self.parent.grandparent 931 else: 932 return self.parent 933 934 def index_in_grandparent(self, indices): 935 """ 936 Given an array of indices, return the indices in the parent population 937 at the root of the tree. 938 """ 939 indices_in_parent = np.arange(self.parent.size)[self.mask][indices] 940 if hasattr(self.parent, "parent"): 941 return self.parent.index_in_grandparent(indices_in_parent) 942 else: 943 return indices_in_parent 944 945 def index_from_parent_index(self, indices): 946 """ 947 Given an index(indices) in the parent population, return 948 the index(indices) within this view. 949 """ 950 # todo: add check that all indices correspond to cells that are in this view 951 if isinstance(self.mask, slice): 952 start = self.mask.start or 0 953 step = self.mask.step or 1 954 return (indices - start) / step 955 else: 956 if isinstance(indices, int): 957 return np.nonzero(self.mask == indices)[0][0] 958 elif isinstance(indices, np.ndarray): 959 # Lots of ways to do this. Some profiling is in order. 960 # - https://stackoverflow.com/questions/16992713/translate-every-element-in-numpy-array-according-to-key 961 # - https://stackoverflow.com/questions/3403973/fast-replacement-of-values-in-a-numpy-array 962 # - https://stackoverflow.com/questions/13572448/replace-values-of-a-numpy-index-array-with-values-of-a-list 963 parent_indices = self.mask # assert mask is sorted 964 view_indices = np.arange(self.size) 965 index = np.digitize(indices, parent_indices, right=True) 966 return view_indices[index] 967 else: 968 raise ValueError("indices must be an integer or an array of integers") 969 970 def __eq__(self, other): 971 """ 972 Determine whether two views are the same. 973 """ 974 return not self.__ne__(other) 975 976 def __ne__(self, other): 977 """ 978 Determine whether two views are different. 979 """ 980 # We can't use the self.mask, as different masks can select the same cells 981 # (e.g. slices vs arrays), therefore we have to use self.all_cells 982 if isinstance(other, PopulationView): 983 return self.parent != other.parent or not np.array_equal(self.all_cells, other.all_cells) 984 elif isinstance(other, Population): 985 return self.parent != other or not np.array_equal(self.all_cells, other.all_cells) 986 else: 987 return True 988 989 def describe(self, template='populationview_default.txt', engine='default'): 990 """ 991 Returns a human-readable description of the population view. 992 993 The output may be customized by specifying a different template 994 togther with an associated template engine (see ``pyNN.descriptions``). 995 996 If template is None, then a dictionary containing the template context 997 will be returned. 998 """ 999 context = {"label": self.label, 1000 "parent": self.parent.label, 1001 "mask": self.mask, 1002 "size": self.size} 1003 context.update(self.annotations) 1004 return descriptions.render(engine, template, context) 1005 1006 1007class Assembly(object): 1008 """ 1009 A group of neurons, may be heterogeneous, in contrast to a Population where 1010 all the neurons are of the same type. 1011 1012 Arguments: 1013 populations: 1014 Populations or PopulationViews 1015 kwargs: 1016 May contain a keyword argument 'label' 1017 """ 1018 _count = 0 1019 1020 def __init__(self, *populations, **kwargs): 1021 """ 1022 Create an Assembly of Populations and/or PopulationViews. 1023 """ 1024 if not hasattr(self, "_simulator"): 1025 errmsg = "`common.Assembly` should not be instantiated directly. " \ 1026 "You should import Assembly from a PyNN backend module, " \ 1027 "e.g. pyNN.nest or pyNN.neuron" 1028 raise Exception(errmsg) 1029 if kwargs: 1030 assert list(kwargs.keys()) == ['label'] 1031 self.populations = [] 1032 for p in populations: 1033 self._insert(p) 1034 self.label = kwargs.get('label', 'assembly%d' % Assembly._count) 1035 assert isinstance(self.label, str), "label must be a string" 1036 self.annotations = {} 1037 Assembly._count += 1 1038 1039 def __repr__(self): 1040 return "Assembly(*%r, label=%r)" % (self.populations, self.label) 1041 1042 def _insert(self, element): 1043 if not isinstance(element, BasePopulation): 1044 raise TypeError("argument is a %s, not a Population." % type(element).__name__) 1045 if isinstance(element, PopulationView): 1046 if not element.parent in self.populations: 1047 double = False 1048 for p in self.populations: 1049 data = np.concatenate((p.all_cells, element.all_cells)) 1050 if len(np.unique(data)) != len(p.all_cells) + len(element.all_cells): 1051 logging.warning( 1052 'Adding a PopulationView to an Assembly containing elements already present is not posible') 1053 double = True # Should we automatically remove duplicated IDs ? 1054 break 1055 if not double: 1056 self.populations.append(element) 1057 else: 1058 logging.warning( 1059 'Adding a PopulationView to an Assembly when parent Population is there is not possible') 1060 elif isinstance(element, BasePopulation): 1061 if not element in self.populations: 1062 self.populations.append(element) 1063 else: 1064 logging.warning('Adding a Population twice in an Assembly is not possible') 1065 1066 @property 1067 def local_cells(self): 1068 result = self.populations[0].local_cells 1069 for p in self.populations[1:]: 1070 result = np.concatenate((result, p.local_cells)) 1071 return result 1072 1073 @property 1074 def all_cells(self): 1075 result = self.populations[0].all_cells 1076 for p in self.populations[1:]: 1077 result = np.concatenate((result, p.all_cells)) 1078 return result 1079 1080 def all(self): 1081 """Iterator over cell ids on all nodes.""" 1082 return iter(self.all_cells) 1083 1084 @property 1085 def _is_sorted(self): 1086 idx = np.argsort(self.all_cells) 1087 return np.all(idx == np.arange(len(self.all_cells))) 1088 1089 @property 1090 def _homogeneous_synapses(self): 1091 cb = [p.celltype.conductance_based for p in self.populations] 1092 return all(cb) or not any(cb) 1093 1094 @property 1095 def conductance_based(self): 1096 """ 1097 `True` if the post-synaptic response is modelled as a change 1098 in conductance, `False` if a change in current. 1099 """ 1100 return all(p.celltype.conductance_based for p in self.populations) 1101 1102 @property 1103 def receptor_types(self): 1104 """ 1105 Return a list of receptor types that are common to all populations 1106 within the assembly. 1107 """ 1108 rts = self.populations[0].celltype.receptor_types 1109 if len(self.populations) > 1: 1110 rts = set(rts) 1111 for p in self.populations[1:]: 1112 rts = rts.intersection(set(p.celltype.receptor_types)) 1113 return list(rts) 1114 1115 def find_units(self, variable): 1116 """ 1117 Returns units of the specified variable or parameter, as a string. 1118 Works for all the recordable variables and neuron parameters of all standard models. 1119 """ 1120 units = set(p.find_units(variable) for p in self.populations) 1121 if len(units) > 1: 1122 raise ValueError("Inconsistent units") 1123 return units 1124 1125 @property 1126 def _mask_local(self): 1127 result = self.populations[0]._mask_local 1128 for p in self.populations[1:]: 1129 result = np.concatenate((result, p._mask_local)) 1130 return result 1131 1132 @property 1133 def first_id(self): 1134 return np.min(self.all_cells) 1135 1136 @property 1137 def last_id(self): 1138 return np.max(self.all_cells) 1139 1140 def id_to_index(self, id): 1141 """ 1142 Given the ID(s) of cell(s) in the Assembly, return its (their) index 1143 (order in the Assembly):: 1144 1145 >>> assert p.id_to_index(p[5]) == 5 1146 >>> assert p.id_to_index(p.index([1, 2, 3])) == [1, 2, 3] 1147 """ 1148 all_cells = self.all_cells 1149 if not np.iterable(id): 1150 if self._is_sorted: 1151 return np.searchsorted(all_cells, id) 1152 else: 1153 result = np.where(all_cells == id)[0] 1154 if len(result) == 0: 1155 raise IndexError("ID %s not present in the View" % id) 1156 else: 1157 return result 1158 else: 1159 if self._is_sorted: 1160 return np.searchsorted(all_cells, id) 1161 else: 1162 result = np.array([], dtype=np.int) 1163 for item in id: 1164 data = np.where(all_cells == item)[0] 1165 if len(data) == 0: 1166 raise IndexError("ID %s not present in the Assembly" % item) 1167 elif len(data) > 1: 1168 raise Exception("ID %s is duplicated in the Assembly" % item) 1169 else: 1170 result = np.append(result, data) 1171 return result 1172 1173 @property 1174 def positions(self): 1175 result = self.populations[0].positions 1176 for p in self.populations[1:]: 1177 result = np.hstack((result, p.positions)) 1178 return result 1179 1180 @property 1181 def size(self): 1182 return sum(p.size for p in self.populations) 1183 1184 def __iter__(self): 1185 """ 1186 Iterator over cells in all populations within the Assembly, for cells 1187 on the local MPI node. 1188 """ 1189 iterators = [iter(p) for p in self.populations] 1190 return chain(*iterators) 1191 1192 def __len__(self): 1193 """Return the total number of cells in the population (all nodes).""" 1194 return self.size 1195 1196 def __getitem__(self, index): 1197 """ 1198 Where `index` is an integer, return an ID. 1199 Where `index` is a slice, tuple, list or numpy array, return a new Assembly 1200 consisting of appropriate populations and (possibly newly created) 1201 population views. 1202 """ 1203 count = 0 1204 boundaries = [0] 1205 for p in self.populations: 1206 count += p.size 1207 boundaries.append(count) 1208 boundaries = np.array(boundaries, dtype=int) 1209 1210 if isinstance(index, (int, np.integer)): # return an ID 1211 pindex = boundaries[1:].searchsorted(index, side='right') 1212 return self.populations[pindex][index - boundaries[pindex]] 1213 elif isinstance(index, (slice, tuple, list, np.ndarray)): 1214 if isinstance(index, slice) or (hasattr(index, "dtype") and index.dtype == bool): 1215 indices = np.arange(self.size)[index] 1216 else: 1217 indices = np.array(index) 1218 pindices = boundaries[1:].searchsorted(indices, side='right') 1219 views = [self.populations[i][indices[pindices == i] - boundaries[i]] 1220 for i in np.unique(pindices)] 1221 return self.__class__(*views) 1222 else: 1223 raise TypeError("indices must be integers, slices, lists, arrays, not %s" % 1224 type(index).__name__) 1225 1226 def __add__(self, other): 1227 """ 1228 An Assembly may be added to a Population, PopulationView or Assembly 1229 with the '+' operator, returning a new Assembly, e.g.:: 1230 1231 a2 = a1 + p 1232 """ 1233 if isinstance(other, BasePopulation): 1234 return self.__class__(*(self.populations + [other])) 1235 elif isinstance(other, Assembly): 1236 return self.__class__(*(self.populations + other.populations)) 1237 else: 1238 raise TypeError("can only add a Population or another Assembly to an Assembly") 1239 1240 def __iadd__(self, other): 1241 """ 1242 A Population, PopulationView or Assembly may be added to an existing 1243 Assembly using the '+=' operator, e.g.:: 1244 1245 a += p 1246 """ 1247 if isinstance(other, BasePopulation): 1248 self._insert(other) 1249 elif isinstance(other, Assembly): 1250 for p in other.populations: 1251 self._insert(p) 1252 else: 1253 raise TypeError("can only add a Population or another Assembly to an Assembly") 1254 return self 1255 1256 def sample(self, n, rng=None): 1257 """ 1258 Randomly sample `n` cells from the Assembly, and return a Assembly 1259 object. 1260 """ 1261 assert isinstance(n, int) 1262 if not rng: 1263 rng = random.NumpyRNG() 1264 indices = rng.permutation(np.arange(len(self), dtype=np.int))[0:n] 1265 logger.debug("The %d cells recorded have indices %s" % (n, indices)) 1266 logger.debug("%s.sample(%s)", self.label, n) 1267 return self[indices] 1268 1269 def initialize(self, **initial_values): 1270 """ 1271 Set the initial values of the state variables of the neurons in 1272 this assembly. 1273 """ 1274 for p in self.populations: 1275 p.initialize(**initial_values) 1276 1277 def get(self, parameter_names, gather=False, simplify=True): 1278 """ 1279 Get the values of the given parameters for every local cell in the 1280 Assembly, or, if gather=True, for all cells in the Assembly. 1281 """ 1282 if isinstance(parameter_names, str): 1283 parameter_names = (parameter_names,) 1284 return_list = False 1285 else: 1286 return_list = True 1287 1288 parameters = defaultdict(list) 1289 for p in self.populations: 1290 population_values = p.get(parameter_names, gather, simplify=False) 1291 for name, arr in zip(parameter_names, population_values): 1292 parameters[name].append(arr) 1293 for name, value_list in parameters.items(): 1294 parameters[name] = np.hstack(value_list) 1295 if simplify: 1296 parameters[name] = simplify_parameter_array(parameters[name]) 1297 values = [parameters[name] for name in parameter_names] 1298 if return_list: 1299 return values 1300 else: 1301 assert len(parameter_names) == 1 1302 return values[0] 1303 1304 def set(self, **parameters): 1305 """ 1306 Set one or more parameters for every cell in the Assembly. 1307 1308 Values passed to set() may be: 1309 (1) single values 1310 (2) RandomDistribution objects 1311 (3) mapping functions, where a mapping function accepts a single 1312 argument (the cell index) and returns a single value. 1313 1314 Here, a "single value" may be either a single number or a list/array of 1315 numbers (e.g. for spike times). 1316 """ 1317 for p in self.populations: 1318 p.set(**parameters) 1319 1320 @deprecated("set(parametername=rand_distr)") 1321 def rset(self, parametername, rand_distr): 1322 self.set(parametername=rand_distr) 1323 1324 def record(self, variables, to_file=None, sampling_interval=None): 1325 """ 1326 Record the specified variable or variables for all cells in the Assembly. 1327 1328 `variables` may be either a single variable name or a list of variable 1329 names. For a given celltype class, `celltype.recordable` contains a list of 1330 variables that can be recorded for that celltype. 1331 1332 If specified, `to_file` should be either a filename or a Neo IO instance and `write_data()` 1333 will be automatically called when `end()` is called. 1334 """ 1335 for p in self.populations: 1336 p.record(variables, to_file, sampling_interval) 1337 1338 @deprecated("record('v')") 1339 def record_v(self, to_file=True): 1340 """Record the membrane potential from all cells in the Assembly.""" 1341 self.record('v', to_file) 1342 1343 @deprecated("record(['gsyn_exc', 'gsyn_inh'])") 1344 def record_gsyn(self, to_file=True): 1345 """Record synaptic conductances from all cells in the Assembly.""" 1346 self.record(['gsyn_exc', 'gsyn_inh'], to_file) 1347 1348 def get_population(self, label): 1349 """ 1350 Return the Population/PopulationView from within the Assembly that has 1351 the given label. If no such Population exists, raise KeyError. 1352 """ 1353 for p in self.populations: 1354 if label == p.label: 1355 return p 1356 raise KeyError("Assembly does not contain a population with the label %s" % label) 1357 1358 def save_positions(self, file): 1359 """ 1360 Save positions to file. The output format is id x y z 1361 """ 1362 if isinstance(file, str): 1363 file = files.StandardTextFile(file, mode='w') 1364 cells = self.all_cells 1365 result = np.empty((len(cells), 4)) 1366 result[:, 0] = np.array([self.id_to_index(id) for id in cells]) 1367 result[:, 1:4] = self.positions.T 1368 if self._simulator.state.mpi_rank == 0: 1369 file.write(result, {'assembly': self.label}) 1370 file.close() 1371 1372 @property 1373 def position_generator(self): 1374 def gen(i): 1375 return self.positions[:, i] 1376 return gen 1377 1378 def get_data(self, variables='all', gather=True, clear=False, annotations=None): 1379 """ 1380 Return a Neo `Block` containing the data (spikes, state variables) 1381 recorded from the Assembly. 1382 1383 `variables` - either a single variable name or a list of variable names 1384 Variables must have been previously recorded, otherwise an 1385 Exception will be raised. 1386 1387 For parallel simulators, if `gather` is True, all data will be gathered 1388 to all nodes and the Neo `Block` will contain data from all nodes. 1389 Otherwise, the Neo `Block` will contain only data from the cells 1390 simulated on the local node. 1391 1392 If `clear` is True, recorded data will be deleted from the `Assembly`. 1393 """ 1394 name = self.label 1395 description = self.describe() 1396 blocks = [p.get_data(variables, gather, clear) for p in self.populations] 1397 # adjust channel_ids to match assembly channel indices 1398 offset = 0 1399 for block, p in zip(blocks, self.populations): 1400 for segment in block.segments: 1401 for signal_array in segment.analogsignals: 1402 signal_array.array_annotations["channel_index"] += offset 1403 offset += p.size 1404 for i, block in enumerate(blocks): 1405 logger.debug("%d: %s", i, block.name) 1406 for j, segment in enumerate(block.segments): 1407 logger.debug(" %d: %s", j, segment.name) 1408 for arr in segment.analogsignals: 1409 logger.debug(" %s %s", arr.shape, arr.name) 1410 merged_block = blocks[0] 1411 for block in blocks[1:]: 1412 merged_block.merge(block) 1413 merged_block.name = name 1414 merged_block.description = description 1415 if annotations: 1416 merged_block.annotate(**annotations) 1417 return merged_block 1418 1419 @deprecated("get_data('spikes')") 1420 def getSpikes(self, gather=True, compatible_output=True): 1421 return self.get_data('spikes', gather) 1422 1423 @deprecated("get_data('v')") 1424 def get_v(self, gather=True, compatible_output=True): 1425 return self.get_data('v', gather) 1426 1427 @deprecated("get_data(['gsyn_exc', 'gsyn_inh'])") 1428 def get_gsyn(self, gather=True, compatible_output=True): 1429 return self.get_data(['gsyn_exc', 'gsyn_inh'], gather) 1430 1431 def mean_spike_count(self, gather=True): 1432 """ 1433 Returns the mean number of spikes per neuron. 1434 """ 1435 spike_counts = self.get_spike_counts() 1436 total_spikes = sum(spike_counts.values()) 1437 if self._simulator.state.mpi_rank == 0 or not gather: # should maybe use allgather, and get the numbers on all nodes 1438 return float(total_spikes) / len(spike_counts) 1439 else: 1440 return np.nan 1441 1442 def get_spike_counts(self, gather=True): 1443 """ 1444 Returns the number of spikes for each neuron. 1445 """ 1446 try: 1447 spike_counts = self.populations[0].recorder.count( 1448 'spikes', gather, self.populations[0]._record_filter) 1449 except errors.NothingToWriteError: 1450 spike_counts = {} 1451 for p in self.populations[1:]: 1452 try: 1453 spike_counts.update(p.recorder.count('spikes', gather, p._record_filter)) 1454 except errors.NothingToWriteError: 1455 pass 1456 return spike_counts 1457 1458 def write_data(self, io, variables='all', gather=True, clear=False, annotations=None): 1459 """ 1460 Write recorded data to file, using one of the file formats supported by 1461 Neo. 1462 1463 `io`: 1464 a Neo IO instance 1465 `variables`: 1466 either a single variable name or a list of variable names. 1467 Variables must have been previously recorded, otherwise an 1468 Exception will be raised. 1469 1470 For parallel simulators, if `gather` is True, all data will be gathered 1471 to the master node and a single output file created there. Otherwise, a 1472 file will be written on each node, containing only data from the cells 1473 simulated on that node. 1474 1475 If `clear` is True, recorded data will be deleted from the `Population`. 1476 """ 1477 if isinstance(io, str): 1478 io = recording.get_io(io) 1479 if gather is False and self._simulator.state.num_processes > 1: 1480 io.filename += '.%d' % self._simulator.state.mpi_rank 1481 logger.debug("Recorder is writing '%s' to file '%s' with gather=%s" % ( 1482 variables, io.filename, gather)) 1483 data = self.get_data(variables, gather, clear, annotations) 1484 if self._simulator.state.mpi_rank == 0 or gather is False: 1485 logger.debug("Writing data to file %s" % io) 1486 io.write(data) 1487 1488 @deprecated("write_data(file, 'spikes')") 1489 def printSpikes(self, file, gather=True, compatible_output=True): 1490 self.write_data(file, 'spikes', gather) 1491 1492 @deprecated("write_data(file, 'v')") 1493 def print_v(self, file, gather=True, compatible_output=True): 1494 self.write_data(file, 'v', gather) 1495 1496 @deprecated("write_data(['gsyn_exc', 'gsyn_inh'])") 1497 def print_gsyn(self, file, gather=True, compatible_output=True): 1498 self.write_data(file, ['gsyn_exc', 'gsyn_inh'], gather) 1499 1500 def inject(self, current_source): 1501 """ 1502 Connect a current source to all cells in the Assembly. 1503 """ 1504 for p in self.populations: 1505 current_source.inject_into(p) 1506 1507 @property 1508 def injectable(self): 1509 return all(p.injectable for p in self.populations) 1510 1511 def describe(self, template='assembly_default.txt', engine='default'): 1512 """ 1513 Returns a human-readable description of the assembly. 1514 1515 The output may be customized by specifying a different template 1516 togther with an associated template engine (see ``pyNN.descriptions``). 1517 1518 If template is None, then a dictionary containing the template context 1519 will be returned. 1520 """ 1521 context = {"label": self.label, 1522 "populations": [p.describe(template=None) for p in self.populations]} 1523 return descriptions.render(engine, template, context) 1524 1525 def get_annotations(self, annotation_keys, simplify=True): 1526 """ 1527 Get the values of the given annotations for each population in the Assembly. 1528 """ 1529 if isinstance(annotation_keys, str): 1530 annotation_keys = (annotation_keys,) 1531 annotations = defaultdict(list) 1532 1533 for key in annotation_keys: 1534 is_array_annotation = False 1535 for p in self.populations: 1536 annotation = p.annotations[key] 1537 annotations[key].append(annotation) 1538 is_array_annotation = isinstance(annotation, np.ndarray) 1539 if is_array_annotation: 1540 annotations[key] = np.hstack(annotations[key]) 1541 if simplify: 1542 annotations[key] = simplify_parameter_array(np.array(annotations[key])) 1543 return annotations 1544