The Higher Education and Research forge

Home My Page Projects Code Snippets Project Openings EMULSION public releases
Summary Activity Surveys SCM Listes Sympa

SCM Repository

77504782aebd99268da9fced61518da1cb19cbc2
1 """
2 .. module:: emulsion.agent.core.emulsion_agent
4 .. moduleauthor:: Sébastien Picault <sebastien.picault@inra.fr>
6 Part of this code was adapted from the PADAWAN framework (S. Picault,
7 Univ. Lille).
8 """
11 # EMULSION (Epidemiological Multi-Level Simulation framework)
12 # ===========================================================
13
14 # Contributors and contact:
15 # -------------------------
16
17 #     - Sébastien Picault (sebastien.picault@inra.fr)
18 #     - Yu-Lin Huang
19 #     - Vianney Sicard
20 #     - Sandie Arnoux
21 #     - Gaël Beaunée
22 #     - Pauline Ezanno (pauline.ezanno@inra.fr)
23
24 #     BIOEPAR, INRA, Oniris, Atlanpole La Chantrerie,
25 #     Nantes CS 44307 CEDEX, France
26
27
28 # How to cite:
29 # ------------
30
31 #     S. Picault, Y.-L. Huang, V. Sicard, P. Ezanno (2017). "Enhancing
32 #     Sustainability of Complex Epidemiological Models through a Generic
33 #     Multilevel Agent-based Approach", in: C. Sierra (ed.), 26th
34 #     International Joint Conference on Artificial Intelligence (IJCAI),
35 #     AAAI, p. 374-380. DOI: 10.24963/ijcai.2017/53
36
37
38 # License:
39 # --------
40
41 #    Copyright 2016 INRA and Univ. Lille
42
43 #    Inter Deposit Digital Number: IDDN.FR.001.280043.000.R.P.2018.000.10000
44
45 #    Agence pour la Protection des Programmes,
46 #    54 rue de Paradis, 75010 Paris, France
47
48 #    Licensed under the Apache License, Version 2.0 (the "License");
49 #    you may not use this file except in compliance with the License.
50 #    You may obtain a copy of the License at
51
52 #        http://www.apache.org/licenses/LICENSE-2.0
53
54 #    Unless required by applicable law or agreed to in writing, software
55 #    distributed under the License is distributed on an "AS IS" BASIS,
56 #    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
57 #    See the License for the specific language governing permissions and
58 #    limitations under the License.
60 import abc
62 import numpy                      as np
65 from   emulsion.tools.graph       import EdgeTypes
66 from   emulsion.model             import EDGE_KEYWORDS
67 from   emulsion.tools.misc        import rates_to_probabilities,\
68     aggregate_probabilities, probabilities_to_rates, count_population, AGENTS
70 from   emulsion.agent.core.abstract_agent  import AbstractAgent
73 #  ______                 _     _                                      _
74 # |  ____|               | |   (_)               /\                   | |
75 # | |__   _ __ ___  _   _| |___ _  ___  _ __    /  \   __ _  ___ _ __ | |_
76 # |  __| | '_ ` _ \| | | | / __| |/ _ \| '_ \  / /\ \ / _` |/ _ \ '_ \| __|
77 # | |____| | | | | | |_| | \__ \ | (_) | | | |/ ____ \ (_| |  __/ | | | |_
78 # |______|_| |_| |_|\__,_|_|___/_|\___/|_| |_/_/    \_\__, |\___|_| |_|\__|
79 #                                                      __/ |
80 #                                                     |___/
83 class EmulsionAgent(AbstractAgent):
84     """The EmulsionAgent is the base class for multi-level
85     epidemiologic agents. An EmulsionAgent can represent any entity
86     involved in the epidemiologic model (individual, compartment,
87     herd, [meta]population, etc.).
89     Each agent contains an exchange box (in & out) which is composed
90     of a list of messages. A message is generally a dictionary of
91     Source/Destination (In/Out) and content of exchange with
92     Source/Destination agent.
94     """
95     def __init__(self, model=None, level=None, name=None, host=None, **others):
96         """Initialize the unit with a health state and a name."""
97         super().__init__(**others)
98         # self.statevars.health_state = health_state
99         self.model = model
100         self.level = level
101         self._name = name
102         self._host = host
103         if 'step' not in self.statevars:
104             self.statevars.step = 0
105         # exchange inbox/outbox
106         self._inbox = []
107         self._outbox = []
110     def upper_level(self, init=True):
111         """Return the 'upper level' from this agent, i.e. the first host with
112         a not-None level attribute.
114         TAG: USER
115         """
116         if self.level is not None and not init:
117             return self
118         return self.get_host().upper_level(init=False)
120     @property
121     def name(self):
122         """Return the name of the unit. If no name was provided during
123         instantiation,  the class name is returned.
125         TAG: USER
126         """
127         return self._name if self._name is not None\
128             else super().__repr__() # self.__class__.__name__
130     def __repr__(self):
131         return self.name
134     def evolve(self, machine=None):
135         """This method is aimed at defining what has systematically to
136         be done in the unit at each time step (e.g. age change...). It
137         has to be overriden if needed in subclasses.
139         """
140         self.statevars.step += 1
142     def get_host(self, key=None):
143         """Return the host of the current unit.
145         TAG: USER
146         """
147         return self._host
149     @abc.abstractmethod
150     def get_content(self):
151         """Return either the population (number) of the current unit,
152         or the list of agents contained in the current unit. The
153         output is a tuple with either ('population', qty) or ('agents', list).
155         """
156         pass
158     def do_state_actions(self, event, state_machine, state_name,
159                          population=None, agents=None, **_):
160         """Perform the actions associated to the current state. If the
161         unit is a ViewAgent, actions are actually performed by each
162         unit of the specified agents list, in turn. Otherwise, the
163         actions are performed according to the population, which is
164         expected to be a number.
166         """
167         # if actions are actually associated to the current state of
168         # the state machine...
169            # ... and to the 'event' (enter/exit/stay)
170         if state_name in state_machine.state_actions\
171            and event in state_machine.state_actions[state_name]:
172             # print(f'Doing state actions {event} for {state_machine}\n\tfor {state_name}')
173             # retrieve the list of actions
174             l_actions = state_machine.state_actions[state_name][event]
175             # ask the current unit to perform the actions with the
176             # specified population
177             # print(f'\t{l_actions}')
178             for action in l_actions:
179                 action.execute_action(self, population=population, agents=agents)
181     def do_edge_actions(self, actions=None, population=None, agents=None):
182         """Perform the actions associated to a transition between
183         states. If the unit is a ViewCompartment, actions are actually
184         performed by each unit of the specified agents list, in
185         turn. Otherwise, the actions are performed according to the
186         population, which is expected to be a number.
188         """
189         # # if actions are actually associated to the current edge...
190         #    # ... and to the 'event' (cross)
191         # if from_ in state_machine.edge_actions\
192         #    and to_ in state_machine.edge_actions[from_]\
193         #    and event in state_machine.edge_actions[from_][to_]:
194         #     # retrieve the list of actions
195         #     l_actions = state_machine.edge_actions[from_][to_][event]
196         #     # ask the current unit to perform the actions with the
197         #     # specified population
198         for action in actions:
199 #            print(action)
200             action.execute_action(self, population=population, agents=agents)
203     def next_states_from(self, initial_state, state_machine):
204         """Return a list of tuples composed of:
206         - each of the possible states reachable from the specified
207         initial state (some depending on a condition)
209         - the transition rate, probability or amount to each state in
210         the specified state_machine (possible clauses: 'rate',
211         'proba', 'amount', 'amount-all-but' with the actual value)
213         - a tuple indicating who can cross the edge, depending on
214         conditions (either ('population', qty) or ('agents', list))
216         - the list of actions to perform when crossing the edge (if
217         any).
219         The conditions, rates and probabilities may depend on the state
220         variables or properties of the current unit.
222         """
223         result = []
224         # remove unfulfilled 'when' clauses if any
225         for (state, props) in state_machine.graph.edges_from(initial_state,
226                                                              type_id=EdgeTypes.TRANSITION):
227             if 'when' in props:
228                 when = state_machine.get_value(props['when'])
229                 fulfilled = when(self) if callable(when) else when
230                 if not fulfilled:
231                     continue
232             cond_result = self.get_content()
233             # if any condition, evaluate it
234             if 'cond' in props:
235                 cond = state_machine.get_value(props['cond'])
236                 # compute the content tuple (either ('population', qty)
237                 # or ('agents', lsit)) which fulfils the condition
238                 cond_result = self.evaluate_condition(cond)
239             # only edges with condition fulfilled are taken into account
240             if cond_result:
241                 # print(cond_result)
242                 flux = None
243                 for keyword in EDGE_KEYWORDS:
244                     if keyword in props:
245                         flux = keyword
246                         break
247                 value = state_machine.get_value(props[flux])
248                 actions = props['actions'] if 'actions' in props else []
249                 if callable(value):
250                     value = value(self)
251                 if value > 0 or flux == 'amount-all-but':
252                     result.append((state, flux, value, cond_result, actions))
253         return result
255     def production_from(self, initial_state, state_machine):
256         """Return a list of tuples composed of:
258         - each of the possible states that can be produced from the specified
259         initial state (some depending on a condition)
261         - the transition rate, probability or amount to each state in
262         the specified state_machine (possible clauses: 'rate',
263         'proba', 'amount', 'amount-all-but' with the actual value)
265         - a tuple indicating who can produce cross the edge, depending on
266         conditions (either ('population', qty) or ('agents', list))
268         - the prototypes of agents produced through the edge (if
269         any).
271         The conditions, rates and probabilities may depend on the state
272         variables or properties of the current unit.
274         """
275         result = []
276         # remove unfulfilled 'when' clauses if any
277         for (state, props) in state_machine.graph.edges_from(initial_state,
278                                                              type_id=EdgeTypes.PRODUCTION):
279             if 'when' in props:
280                 when = state_machine.get_value(props['when'])
281                 fulfilled = when(self) if callable(when) else when
282                 if not fulfilled:
283                     continue
284             cond_result = self.get_content()
285             # if any condition, evaluate it
286             if 'cond' in props:
287                 cond = state_machine.get_value(props['cond'])
288                 # compute the content tuple (either ('population', qty)
289                 # or ('agents', lsit)) which fulfils the condition
290                 cond_result = self.evaluate_condition(cond)
291             # only edges with condition fulfilled are taken into account
292             if cond_result:
293                 # print(cond_result)
294                 flux = None
295                 for keyword in EDGE_KEYWORDS:
296                     if keyword in props:
297                         flux = keyword
298                         break
299                 value = state_machine.get_value(props[flux])
300                 protos = props['prototype'] if 'prototype' in props else None
301                 if callable(value):
302                     value = value(self)
303                 if value > 0 or flux == 'amount-all-but':
304                     result.append((state, flux, value, cond_result, protos))
305         return result
308     def _compute_production(self, value, flux, pop_size, stochastic):
309         """Compute the values of production edges.
311         """
312         # ALL PRODUCTIONS COME FROM THE SAME POPULATION
313         # (reference_pop) The `values` can represent : 1) amounts - in
314         # that case no transformation occurs, 2) probabilities - in
315         # that case the values must be converted to rates if the
316         # target compartment is deterministic, otherwise the step
317         # duration must be taken into account; 3) rates - computed
318         # using Poisson distribution if the target compartment is
319         # stochastic
320         if flux == 'amount-all-but':
321             amount = max(0, pop_size - value)
322         elif flux == 'amount':
323             amount = value
324         elif flux == 'proba':
325             # if stochastic uses binomial
326             if stochastic:
327                 # aggregate probability wrt the time step duration
328                 value = aggregate_probabilities([value], self.model.delta_t)[0]
329                 amount = np.random.binomial(pop_size, value)
330             else:
331                 # otherwise treat proba as proportion
332                 amount = pop_size * value
333         elif flux == 'rate':
334             # if stochastic: use Poisson distribution
335             if stochastic:
336                 # amount = np.random.poisson(value * pop_size * self.model.delta_t)
337                 amount = np.random.binomial(pop_size, 1-np.exp(-value * self.model.delta_t))
338             else:
339                 # consider rate as a speed
340                 # S = S0 exp(rate * dt) => S += S * (exp(rate*dt) - 1)
341                 # amount = value * pop_size * self.model.delta_t
342                 amount = pop_size * (np.exp(value * self.model.delta_t) - 1)
343         else:                   # should not happen !
344             amount = None
345         return amount
348     def _compute_values_for_unique_population(self,
349                                               values,
350                                               flux,
351                                               reference_pop,
352                                               stochastic):
353         """Restructure the values according to the situation, for
354         edges concerning the same population.
356         """
357         # ALL TRANSITIONS AFFECT THE SAME POPULATION (reference_pop)
358         # The `values` can represent : 1) amounts - in that case no
359         # transformation occurs, 2) probabilities - in that case the
360         # values must be converted to rates if the target compartment
361         # is deterministic, otherwise the step duration must be taken
362         # into account; 3) rates - in that case the values must be
363         # converted to probabilities if the target compartment is
364         # stochastic
365         # print('IDENTICAL')
366         available_flux = set(flux)
367         # try:
368         assert len(available_flux) == 1 or available_flux == set(['amount', 'amount-all-but'])
369         # except:
370         #     print(available_flux)
372         method = None
373         if 'amount' in available_flux or 'amount-all-but' in available_flux:
374             # handle values as amounts
375             method = 'amount'
376             total_ref_pop = count_population(reference_pop)
377             # check that values are between 0 and the population size,
378             # if needed invert 'amount-all-but' values
379             values = tuple([max(0, min(total_ref_pop-v, total_ref_pop))\
380                               if f == 'amount-all-but'\
381                               else max(0, min(v, total_ref_pop))
382                             for (f, v) in zip(flux, values)])
383             # when finished the length of values is the number of
384             # outgoing edges
385             # print('AMOUNT', values)
386         elif 'proba' in available_flux:
387             if not stochastic:
388                 # then if the target compartment is deterministic,
389                 # probabilities must be converted into rates
390                 # print('PROBA -> RATES', values)
391                 values = probabilities_to_rates(values + (1 - sum(values),))
392                 # when finished the length of values is the number of
393                 # outgoing edges
394                 # print(values)
395             else:
396                 # aggregate probabilities wrt the time step duration
397                 values = aggregate_probabilities(values, self.model.delta_t)
398                 values = values + (1 - sum(values),)
399                 # when finished the length of values is the number of
400                 # outgoing edges + 1
401                 # print('PROBA', values)
402         elif not stochastic:
403             # print('RATES', values)
404             pass
405         else:
406             # otherwise values are transformed from rates to
407             # probabilities
408             values = rates_to_probabilities(sum(values),
409                                             values,
410                                             delta_t=self.model.delta_t)
411             # when finished the length of values is the number of
412             # outgoing edges + 1
413             # print("RATES -> PROBAS", values)
414         return values, method
416     def _compute_values_for_multiple_populations(self,
417                                                  values,
418                                                  flux,
419                                                  populations,
420                                                  stochastic):
421         """Restructure the values according to the situation, for
422         edges concerning distinct populations.
424         """
425         # IN THAT CASE, EACH POSSIBLE TRANSITION IS RELATED TO A
426         # SPECIFIC SUBGROUP OF THE COMPART.
427         ### IMPORTANT: we assume that all conditions are disjoint,
428         ### i.e. there is no intersection between the populations. IF
429         ### NOT THE CASE, this means that the model is not really
430         ### consistent, and the calculation of probabilities should be
431         ### done on each individual... thus why use a StructuredView
432         ### instead of a set of EvolvingAtom ???
433         # print('MULTIPLE')
434         pop_sets = [SortedSet(pop[AGENTS]) for pop in populations]
435         # pop_sets must be disjoint
436         assert(not any([pop_sets[i] & pop_sets[j]
437                         for i in range(len(pop_sets)-1)
438                         for j in range(i+1, len(pop_sets))]))
439         # binomial values must be computed for each transition to
440         # determine how many units among the candidates will cross the
441         # transition, thus we need probabilities. If all values are
442         # given as rates, they must be transformed
443         method = None
444         available_flux = set(flux)
445         # try:
446         assert len(available_flux) == 1 or available_flux == set(['amount', 'amount-all-but'])
447         # except:
448         #     print(available_flux)
450         if available_flux == {'rate'} and stochastic:
451             # values are transformed from rates to probabilities
452             values = rates_to_probabilities(sum(values),
453                                             values,
454                                             delta_t=self.model.delta_t)
455             # print('RATES -> PROBAS', values)
456         elif 'amount' in available_flux or 'amount-all-but' in available_flux:
457             method = 'amount'
458             pops = [len(p) for p in pop_sets]
459             values = tuple([max(0, min(pop-v, pop))\
460                               if f == 'amount-all-but'\
461                               else max(0, min(v, pop))
462                             for (f, v, pop) in zip(flux, values, pops)])
463             values = values #+ (1 - sum(values),)
464             # print('AMOUNTS:', values)
465         elif 'proba' in available_flux:
466             if stochastic:
467                 # values must be aggregated wrt the times step duration
468                 values = aggregate_probabilities(values,
469                                                  delta_t=self.model.delta_t)
470                 values = values + (1 - sum(values),)
471                 # print('PROBAS:', values)
472             else:
473                 # values must be transformed from probabilities to rates
474                 # print('PROBABILITIES -> RATES', values)
475                 values = probabilities_to_rates(values + (1 - sum(values),))
476                 # print(values)
477         else:
478             print('ERROR - DETERMINISTIC COMPARTMENT MODEL MUST NOT'
479                   ' TRY TO HANDLE SEPARATE SUBPOPULATIONS')
480         return values, method
483     def evaluate_condition(self, condition):
484         """Return the content (tuple either ('population', qty) or
485         ('agents', list)) if the condition is fulfilled, () otherwise.
487         """
488         if callable(condition):
489             condition = condition(self)
490         return self.get_content() if condition else ()
492     def evaluate_event(self, name):
493         """Evaluate if the specified event name is fulfilled, using a
494         calendar. The calendar is determined dynamically according to
495         the name of the event.
497         TAG: USER
498         """
499         # print(self, 'evaluating event', name, "at", self.statevars.step)
500         return self.model.get_calendar_for_event(name)[name](self.statevars.step)
502     def is_in_state(self, state_name):
503         """Evaluate if the agent is in the specified state. The state name is
504         expected to be a valid state of one of the state machines. If
505         so, the name of the state machine is expected to be one of the
506         statevars of the agent.
508         """
509         # retrieve name of the state machine associated with the state_name
510         varname = self.get_model_value(state_name).state_machine.machine_name
511         result = self.statevars[varname].name == state_name\
512                  if varname in self.statevars else False
513         # print(self, 'testing var:', varname, ' = ', result)
514         return result
516     # def time_before_event(self, event_name):
517     #     """Return the time to wait from the current time to the specified
518     #     *event_name*.
519     #     """
520     #     pass self.model.get_cal
522     def get_outbox(self):
523         """Return the content of the outbox.
525         TAG: USER
526         """
527         return self._outbox
529     def add_outbox(self, message):
530         """Appendd a message to the outbox.
532         TAG: USER
533         """
534         return self._outbox.append(message)
536     def reset_outbox(self):
537         """Reset the content of the outbox.
539         TAG: USER
540         """
541         self._outbox = []
543     def add_inbox(self, messages=[]):
544         """Add the specified list of messages in the inbox.
546         TAG: USER"""
547         self._inbox += messages
549     def clean_inbox(self):
550         """Remove non sticky messages from the inbox.
552         TAG: USER
553         """
554         self._inbox = [message for message in self._inbox if message.get('sticky')]
556     def checkout_inbox(self):
557         """Pick up agent's inbox content.
559         TAG: USER
560         """
561         for message in self._inbox:
562             for name, value in message.items():
563                 self.set_information(name, value)