Completed
Push — master ( 27be37...a60125 )
by Andrei
01:07
created

legion_network.__len__()   A

Complexity

Conditions 2

Size

Total Lines 10

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 2
dl 0
loc 10
rs 9.4285
c 0
b 0
f 0
1
"""!
2
3
@brief Neural Network: Local Excitatory Global Inhibitory Oscillatory Network (LEGION)
4
@details Based on article description:
5
         - D.Wang, D.Terman. Image Segmentation Based on Oscillatory Correlation. 1997.
6
         - D.Wang, D.Terman. Locally Excitatory Globally Inhibitory Oscillator Networks. 1995.
7
8
@authors Andrei Novikov ([email protected])
9
@date 2014-2017
10
@copyright GNU Public License
11
12
@cond GNU_PUBLIC_LICENSE
13
    PyClustering is free software: you can redistribute it and/or modify
14
    it under the terms of the GNU General Public License as published by
15
    the Free Software Foundation, either version 3 of the License, or
16
    (at your option) any later version.
17
    
18
    PyClustering is distributed in the hope that it will be useful,
19
    but WITHOUT ANY WARRANTY; without even the implied warranty of
20
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
21
    GNU General Public License for more details.
22
    
23
    You should have received a copy of the GNU General Public License
24
    along with this program.  If not, see <http://www.gnu.org/licenses/>.
25
@endcond
26
27
"""
28
29
import numpy;
0 ignored issues
show
Configuration introduced by
The import numpy could not be resolved.

This can be caused by one of the following:

1. Missing Dependencies

This error could indicate a configuration issue of Pylint. Make sure that your libraries are available by adding the necessary commands.

# .scrutinizer.yml
before_commands:
    - sudo pip install abc # Python2
    - sudo pip3 install abc # Python3
Tip: We are currently not using virtualenv to run pylint, when installing your modules make sure to use the command for the correct version.

2. Missing __init__.py files

This error could also result from missing __init__.py files in your module folders. Make sure that you place one file in each sub-folder.

Loading history...
30
import math;
31
import random;
32
33
import pyclustering.core.legion_wrapper as wrapper;
34
35
from pyclustering.nnet import *;
0 ignored issues
show
Unused Code introduced by
IntEnum was imported with wildcard, but is not used.
Loading history...
Unused Code introduced by
initial_type was imported with wildcard, but is not used.
Loading history...
36
37
from pyclustering.utils import heaviside, allocate_sync_ensembles;
38
39
from scipy.integrate import odeint;
0 ignored issues
show
Configuration introduced by
The import scipy.integrate could not be resolved.

This can be caused by one of the following:

1. Missing Dependencies

This error could indicate a configuration issue of Pylint. Make sure that your libraries are available by adding the necessary commands.

# .scrutinizer.yml
before_commands:
    - sudo pip install abc # Python2
    - sudo pip3 install abc # Python3
Tip: We are currently not using virtualenv to run pylint, when installing your modules make sure to use the command for the correct version.

2. Missing __init__.py files

This error could also result from missing __init__.py files in your module folders. Make sure that you place one file in each sub-folder.

Loading history...
40
41
42
class legion_parameters:
43
    """!
44
    @brief Describes parameters of LEGION.
45
    @details Contained parameters affect on output dynamic of each oscillator of the network.
46
    
47
    @see legion_network
48
    
49
    """  
50
    
51
    def __init__(self):
52
        """!
53
        @brief    Default constructor of parameters for LEGION (local excitatory global inhibitory oscillatory network).
54
        @details  Constructor initializes parameters by default non-zero values that can be
55
                  used for simple simulation.
56
        """
57
        
58
        ## Coefficient that affects intrinsic inhibitor of each oscillator. Should be the same as 'alpha'.
59
        self.eps         = 0.02;
60
        
61
        ## Coefficient is chosen to be on the same order of magnitude as 'eps'. Affects on exponential function that decays on a slow time scale.
62
        self.alpha       = 0.005;
63
        
64
        ## Coefficient that is used to control the ratio of the times that the solution spends in these two phases. For a larger value of g, the solution spends a shorter time in the active phase.
65
        self.gamma       = 6.0;
66
        
67
        ## Coefficient that affects on intrinsic inhibitor of each oscillator. Specifies the steepness of the sigmoid function.
68
        self.betta       = 0.1;
69
        
70
        ## Scale coefficient that is used by potential, should be greater than 0.
71
        self.lamda       = 0.1;
72
        
73
        ## Threshold that should be exceeded by a potential to switch on potential.
74
        self.teta        = 0.9;
75
        
76
        ## Threshold that should be exceeded by a single oscillator to affect its neighbors.
77
        self.teta_x      = -1.5;
78
        
79
        ## Threshold that should be exceeded to activate potential. If potential less than the threshold then potential is relaxed to 0 on time scale 'mu'.
80
        self.teta_p      = 1.5;
81
        
82
        ## Threshold that should be exceeded by any oscillator to activate global inhibitor.
83
        self.teta_xz     = 0.1;
84
        
85
        ## Threshold that should be exceeded to affect on a oscillator by the global inhibitor.
86
        self.teta_zx     = 0.1;
87
        
88
        ## Weight of permanent connections.
89
        self.T           = 2.0;
90
        
91
        ## Defines time scaling of relaxing of oscillator potential.
92
        self.mu          = 0.01;
93
        
94
        ## Weight of global inhibitory connections.
95
        self.Wz          = 1.5;
96
        
97
        ## Total dynamic weights to a single oscillator from neighbors. Sum of weights of dynamic connections to a single oscillator can not be bigger than Wt.
98
        self.Wt          = 8.0;
99
        
100
        ## Rate at which the global inhibitor reacts to the stimulation from the oscillator network.
101
        self.fi          = 3.0;
102
        
103
        ## Multiplier of oscillator noise. Plays important role in desynchronization process.
104
        self.ro          = 0.02;
105
        
106
        ## Value of external stimulus.
107
        self.I           = 0.2;
108
        
109
        ## Defines whether to use potentional of oscillator or not.
110
        self.ENABLE_POTENTIONAL = True;
111
112
113
class legion_dynamic:
114
    """!
115
    @brief Represents output dynamic of LEGION.
116
    
117
    """
118
    
119
    @property
120
    def output(self):
121
        """!
122
        @brief Returns output dynamic of the network.
123
        
124
        """
125
        if (self.__ccore_legion_dynamic_pointer is not None):
126
            return wrapper.legion_dynamic_get_output(self.__ccore_legion_dynamic_pointer);
127
            
128
        return self.__output;
129
    
130
131
    @property
132
    def inhibitor(self):
133
        """!
134
        @brief Returns output dynamic of the global inhibitor of the network.
135
        
136
        """
137
        
138
        if (self.__ccore_legion_dynamic_pointer is not None):
139
            return wrapper.legion_dynamic_get_inhibitory_output(self.__ccore_legion_dynamic_pointer);
140
            
141
        return self.__inhibitor;
142
    
143
    
144
    @property
145
    def time(self):
146
        """!
147
        @brief Returns simulation time.
148
        
149
        """
150
        if (self.__ccore_legion_dynamic_pointer is not None):
151
            return wrapper.legion_dynamic_get_time(self.__ccore_legion_dynamic_pointer);
152
        
153
        return list(range(len(self)));
154
    
155
    
156
    def __init__(self, output, inhibitor, time, ccore = None):
157
        """!
158
        @brief Constructor of legion dynamic.
159
        
160
        @param[in] output (list): Output dynamic of the network represented by excitatory values of oscillators.
161
        @param[in] inhibitor (list): Output dynamic of the global inhibitor of the network.
162
        @param[in] time (list): Simulation time.
163
        @param[in] ccore (POINTER): Pointer to CCORE legion_dynamic. If it is specified then others arguments can be omitted.
164
        
165
        """
166
        
167
        self.__output = output;
168
        self.__inhibitor = inhibitor;
169
        self._time = time;
170
        
171
        self.__ccore_legion_dynamic_pointer = ccore;
172
        
173
        
174
    def __del__(self):
175
        """!
176
        @brief Destructor of the dynamic of the legion network.
177
        
178
        """
179
        if (self.__ccore_legion_dynamic_pointer is not None):
180
            wrapper.legion_dynamic_destroy(self.__ccore_legion_dynamic_pointer);
181
182
183
    def __len__(self):
184
        """!
185
        @brief Returns length of output dynamic.
186
        
187
        """
188
        if (self.__ccore_legion_dynamic_pointer is not None):
189
            return wrapper.legion_dynamic_get_size(self.__ccore_legion_dynamic_pointer);
190
        
191
        return len(self._time);
192
193
194
    def allocate_sync_ensembles(self, tolerance = 0.1):
195
        """!
196
        @brief Allocate clusters in line with ensembles of synchronous oscillators where each synchronous ensemble corresponds to only one cluster.
197
        
198
        @param[in] tolerance (double): Maximum error for allocation of synchronous ensemble oscillators.
199
        
200
        @return (list) Grours of indexes of synchronous oscillators, for example, [ [index_osc1, index_osc3], [index_osc2], [index_osc4, index_osc5] ].
201
        
202
        """
203
204
        if (self.__ccore_legion_dynamic_pointer is not None):
205
            self.__output = wrapper.legion_dynamic_get_output(self.__ccore_legion_dynamic_pointer);
206
            
207
        return allocate_sync_ensembles(self.__output, tolerance);
208
209
210
class legion_network(network):
211
    """!
212
    @brief Local excitatory global inhibitory oscillatory network (LEGION) that uses relaxation oscillator
213
           based on Van der Pol model. 
214
           
215
    @details The model uses global inhibitor to de-synchronize synchronous ensembles of oscillators.
216
    
217
    Example:
218
    @code
219
        # Create parameters of the network
220
        parameters = legion_parameters();
221
        parameters.Wt = 4.0;
222
        
223
        # Create stimulus
224
        stimulus = [1, 1, 0, 0, 0, 1, 1, 1];
225
        
226
        # Create the network (use CCORE for fast solving)
227
        net = legion_network(len(stimulus), parameters, conn_type.GRID_FOUR, ccore = True);
228
        
229
        # Simulate network - result of simulation is output dynamic of the network
230
        output_dynamic = net.simulate(1000, 750, stimulus);
231
        
232
        # Draw output dynamic
233
        draw_dynamics(output_dynamic.time, output_dynamic.output, x_title = "Time", y_title = "x(t)");
234
    @endcode
235
    
236
    """
237
238
    def __init__(self, num_osc, parameters = None, type_conn = conn_type.ALL_TO_ALL, type_conn_represent = conn_represent.MATRIX, ccore = False):
239
        """!
240
        @brief Constructor of oscillatory network LEGION (local excitatory global inhibitory oscillatory network).
241
        
242
        @param[in] num_osc (uint): Number of oscillators in the network.
243
        @param[in] parameters (legion_parameters): Parameters of the network that are defined by structure 'legion_parameters'.
244
        @param[in] type_conn (conn_type): Type of connection between oscillators in the network.
245
        @param[in] type_conn_represent (conn_represent): Internal representation of connection in the network: matrix or list.
246
        @param[in] ccore (bool): If True then all interaction with object will be performed via CCORE library (C++ implementation of pyclustering).
247
        
248
        """
249
        
250
        self._params = None;                 # parameters of the network
251
    
252
        self.__ccore_legion_pointer = None;
253
        self._params = parameters;
254
        
255
        # set parameters of the network
256
        if (self._params is None):
257
            self._params = legion_parameters();
258
        
259
        if (ccore is True):
260
            self.__ccore_legion_pointer = wrapper.legion_create(num_osc, type_conn, self._params);
261
            
262
        else: 
263
            super().__init__(num_osc, type_conn, type_conn_represent);
264
                
265
            # initial states
266
            self._excitatory = [ random.random() for _ in range(self._num_osc) ];
267
            self._inhibitory = [0.0] * self._num_osc;
268
            self._potential = [0.0] * self._num_osc;
269
            
270
            self._coupling_term = None;      # coupling term of each oscillator
271
            self._global_inhibitor = 0;      # value of global inhibitory
272
            self._stimulus = None;           # stimulus of each oscillator
273
            
274
            self._dynamic_coupling = None;   # dynamic connection between oscillators
275
            self._coupling_term = [0.0] * self._num_osc;
276
            self._buffer_coupling_term = [0.0] * self._num_osc;
277
                
278
            # generate first noises
279
            self._noise = [random.random() * self._params.ro for i in range(self._num_osc)];
280
281
282
    def __del__(self):
283
        """!
284
        @brief Default destructor of LEGION.
285
        
286
        """
287
        if (self.__ccore_legion_pointer is not None):
288
            wrapper.legion_destroy(self.__ccore_legion_pointer);
289
            self.__ccore_legion_pointer = None;
290
291
292
    def __len__(self):
293
        """!
294
        @brief (uint) Returns size of LEGION.
295
        
296
        """
297
        
298
        if (self.__ccore_legion_pointer is not None):
299
            return wrapper.legion_get_size(self.__ccore_legion_pointer);
300
        
301
        return self._num_osc;
302
303
304
    def __create_stimulus(self, stimulus):
305
        """!
306
        @brief Create stimulus for oscillators in line with stimulus map and parameters.
307
        
308
        @param[in] stimulus (list): Stimulus for oscillators that is represented by list, number of stimulus should be equal number of oscillators.
309
        
310
        """
311
        
312
        if (len(stimulus) != self._num_osc):
313
            raise NameError("Number of stimulus should be equal number of oscillators in the network.");
314
        else:
315
            self._stimulus = [];
316
             
317
            for val in stimulus:
318
                if (val > 0): self._stimulus.append(self._params.I);
319
                else: self._stimulus.append(0);
320
    
321
    
322
    def __create_dynamic_connections(self):
323
        """!
324
        @brief Create dynamic connection in line with input stimulus.
325
        
326
        """
327
        
328
        if (self._stimulus is None):
329
            raise NameError("Stimulus should initialed before creation of the dynamic connections in the network.");
330
        
331
        self._dynamic_coupling = [ [0] * self._num_osc for i in range(self._num_osc)];
332
        
333
        for i in range(self._num_osc):
334
            neighbors = self.get_neighbors(i);
335
            
336
            if ( (len(neighbors) > 0) and (self._stimulus[i] > 0) ):
337
                number_stimulated_neighbors = 0.0;
338
                for j in neighbors:
339
                    if (self._stimulus[j] > 0):
340
                        number_stimulated_neighbors += 1.0;
341
                
342
                if (number_stimulated_neighbors > 0):
343
                    dynamic_weight = self._params.Wt / number_stimulated_neighbors;
344
                    
345
                    for j in neighbors:
346
                        self._dynamic_coupling[i][j] = dynamic_weight;
347
    
348
    
349
    def simulate(self, steps, time, stimulus, solution = solve_type.RK4, collect_dynamic = True):
350
        """!
351
        @brief Performs static simulation of LEGION oscillatory network.
352
        
353
        @param[in] steps (uint): Number steps of simulations during simulation.
354
        @param[in] time (double): Time of simulation.
355
        @param[in] stimulus (list): Stimulus for oscillators, number of stimulus should be equal to number of oscillators,
356
                   example of stimulus for 5 oscillators [0, 0, 1, 1, 0], value of stimulus is defined by parameter 'I'.
357
        @param[in] solution (solve_type): Method that is used for differential equation.
358
        @param[in] collect_dynamic (bool): If True - returns whole dynamic of oscillatory network, otherwise returns only last values of dynamics.
359
        
360
        @return (list) Dynamic of oscillatory network. If argument 'collect_dynamic' = True, than return dynamic for the whole simulation time,
361
                otherwise returns only last values (last step of simulation) of dynamic.
362
        
363
        """
364
        
365
        if (self.__ccore_legion_pointer is not None):
366
            pointer_dynamic = wrapper.legion_simulate(self.__ccore_legion_pointer, steps, time, solution, collect_dynamic, stimulus);
367
            return legion_dynamic(None, None, None, pointer_dynamic);
368
        
369
        # Check solver before simulation
370
        if (solution == solve_type.FAST):
371
            raise NameError("Solver FAST is not support due to low accuracy that leads to huge error.");
372
        
373
        elif (solution == solve_type.RKF45):
374
            raise NameError("Solver RKF45 is not support in python version. RKF45 is supported in CCORE implementation.");
375
        
376
        # set stimulus
377
        self.__create_stimulus(stimulus);
378
            
379
        # calculate dynamic weights
380
        self.__create_dynamic_connections();
381
        
382
        dyn_exc = None;
383
        dyn_time = None;
384
        dyn_ginh = None;
385
        
386
        # Store only excitatory of the oscillator
387
        if (collect_dynamic == True):
388
            dyn_exc = [];
389
            dyn_time = [];
390
            dyn_ginh = [];
391
            
392
        step = time / steps;
393
        int_step = step / 10.0;
394
        
395
        for t in numpy.arange(step, time + step, step):
396
            # update states of oscillators
397
            self._calculate_states(solution, t, step, int_step);
398
            
399
            # update states of oscillators
400
            if (collect_dynamic == True):
401
                dyn_exc.append(self._excitatory);
402
                dyn_time.append(t);
403
                dyn_ginh.append(self._global_inhibitor);
404
            else:
405
                dyn_exc = self._excitatory;
406
                dyn_time = t;
407
                dyn_ginh = self._global_inhibitor;
408
        
409
        return legion_dynamic(dyn_exc, dyn_ginh, dyn_time); 
410
    
411
    
412
    def _calculate_states(self, solution, t, step, int_step):
0 ignored issues
show
Unused Code introduced by
The argument solution seems to be unused.
Loading history...
413
        """!
414
        @brief Calculates new state of each oscillator in the network.
415
        
416
        @param[in] solution (solve_type): Type solver of the differential equation.
417
        @param[in] t (double): Current time of simulation.
418
        @param[in] step (double): Step of solution at the end of which states of oscillators should be calculated.
419
        @param[in] int_step (double): Step differentiation that is used for solving differential equation.
420
        
421
        """
422
        
423
        next_excitatory = [0.0] * self._num_osc;
424
        next_inhibitory = [0.0] * self._num_osc;
425
        
426
        next_potential = [];
427
        if (self._params.ENABLE_POTENTIONAL is True):
428
            next_potential = [0.0] * self._num_osc;
429
        
430
        # Update states of oscillators
431
        for index in range (0, self._num_osc, 1):
432
            if (self._params.ENABLE_POTENTIONAL is True):
433
                result = odeint(self._legion_state, [self._excitatory[index], self._inhibitory[index], self._potential[index]], numpy.arange(t - step, t, int_step), (index , ));
434
                [ next_excitatory[index], next_inhibitory[index], next_potential[index] ] = result[len(result) - 1][0:3];
435
                
436
            else:
437
                result = odeint(self._legion_state_simplify, [self._excitatory[index], self._inhibitory[index] ], numpy.arange(t - step, t, int_step), (index , ));
438
                [ next_excitatory[index], next_inhibitory[index] ] = result[len(result) - 1][0:2];
439
            
440
            # Update coupling term
441
            neighbors = self.get_neighbors(index);
442
            
443
            coupling = 0
444
            for index_neighbor in neighbors:
445
                coupling += self._dynamic_coupling[index][index_neighbor] * heaviside(self._excitatory[index_neighbor] - self._params.teta_x);
446
            
447
            self._buffer_coupling_term[index] = coupling - self._params.Wz * heaviside(self._global_inhibitor - self._params.teta_xz);
448
        
449
        # Update state of global inhibitory
450
        result = odeint(self._global_inhibitor_state, self._global_inhibitor, numpy.arange(t - step, t, int_step), (None, ));
451
        self._global_inhibitor = result[len(result) - 1][0];
452
        
453
        self._noise = [random.random() * self._params.ro for i in range(self._num_osc)];
454
        self._coupling_term = self._buffer_coupling_term[:];
455
        self._inhibitory = next_inhibitory[:];
456
        self._excitatory = next_excitatory[:];
457
        
458
        if (self._params.ENABLE_POTENTIONAL is True):
459
            self._potential = next_potential[:];
460
            
461
    
462
    
463
    def _global_inhibitor_state(self, z, t, argv):
0 ignored issues
show
Unused Code introduced by
The argument t seems to be unused.
Loading history...
Unused Code introduced by
The argument argv seems to be unused.
Loading history...
464
        """!
465
        @brief Returns new value of global inhibitory
466
        
467
        @param[in] z (dobule): Current value of inhibitory.
468
        @param[in] t (double): Current time of simulation.
469
        @param[in] argv (tuple): It's not used, can be ignored.
470
        
471
        @return (double) New value if global inhibitory (not assign).
472
        
473
        """
474
        
475
        sigma = 0.0;
476
        
477
        for x in self._excitatory:
478
            if (x > self._params.teta_zx):
479
                sigma = 1.0;
480
                break;
481
        
482
        return self._params.fi * (sigma - z);
483
    
484
    
485
    def _legion_state_simplify(self, inputs, t, argv):
0 ignored issues
show
Unused Code introduced by
The argument t seems to be unused.
Loading history...
486
        """!
487
        @brief Returns new values of excitatory and inhibitory parts of oscillator of oscillator.
488
        @details Simplify model doesn't consider oscillator potential.
489
        
490
        @param[in] inputs (list): Initial values (current) of oscillator [excitatory, inhibitory].
491
        @param[in] t (double): Current time of simulation.
492
        @param[in] argv (uint): Extra arguments that are not used for integration - index of oscillator.
493
        
494
        @return (list) New values of excitatoty and inhibitory part of oscillator (not assign).
495
        
496
        """
497
        
498
        index = argv;
499
        
500
        x = inputs[0];  # excitatory
501
        y = inputs[1];  # inhibitory
502
        
503
        dx = 3.0 * x - x ** 3.0 + 2.0 - y + self._stimulus[index] + self._coupling_term[index] + self._noise[index];
504
        dy = self._params.eps * (self._params.gamma * (1.0 + math.tanh(x / self._params.betta)) - y);
505
        
506
        neighbors = self.get_neighbors(index);
507
        potential = 0.0;
0 ignored issues
show
Unused Code introduced by
The variable potential seems to be unused.
Loading history...
508
        
509
        for index_neighbor in neighbors:
510
            potential += self._params.T * heaviside(self._excitatory[index_neighbor] - self._params.teta_x);
511
        
512
        return [dx, dy];    
513
    
514
    
515
    def _legion_state(self, inputs, t, argv):
516
        """!
517
        @brief Returns new values of excitatory and inhibitory parts of oscillator and potential of oscillator.
518
        
519
        @param[in] inputs (list): Initial values (current) of oscillator [excitatory, inhibitory, potential].
520
        @param[in] t (double): Current time of simulation.
521
        @param[in] argv (uint): Extra arguments that are not used for integration - index of oscillator.
522
        
523
        @return (list) New values of excitatoty and inhibitory part of oscillator and new value of potential (not assign).
524
        
525
        """
526
        
527
        index = argv;
528
        
529
        x = inputs[0];  # excitatory
530
        y = inputs[1];  # inhibitory
531
        p = inputs[2];  # potential
532
        
533
        potential_influence = heaviside(p + math.exp(-self._params.alpha * t) - self._params.teta);
534
        
535
        dx = 3.0 * x - x ** 3.0 + 2.0 - y + self._stimulus[index] * potential_influence + self._coupling_term[index] + self._noise[index];
536
        dy = self._params.eps * (self._params.gamma * (1.0 + math.tanh(x / self._params.betta)) - y);
537
        
538
        neighbors = self.get_neighbors(index);
539
        potential = 0.0;
540
        
541
        for index_neighbor in neighbors:
542
            potential += self._params.T * heaviside(self._excitatory[index_neighbor] - self._params.teta_x);
543
        
544
        dp = self._params.lamda * (1.0 - p) * heaviside(potential - self._params.teta_p) - self._params.mu * p;
545
        
546
        return [dx, dy, dp];