AgentBasedModelling/GlobalSA.py at main · Nickvs99/AgentBasedModelling · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
from mesa.batchrunner import BatchRunner
from mesa import Model
import matplotlib.pyplot as plt
from gridmodel import GridModel as model
from SALib.sample import saltelli
from mesa.batchrunner import BatchRunner
from SALib.analyze import sobol
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from itertools import combinations

# define parameters and bounds
problem = {
    'num_vars': 8,
    'names': ['density', 'init_neutral', 'similar_wanted', 'radius','size','network_p', 'randomize_part', 'decrease_intolerance'],
    'bounds': [[0.1, 0.9], [0, 0.9], [0, 1], [1, 5], [5, 100], [0, 0.3], [0, 0.3], [0.9, 0.99]]
}

# set the repetitions, the amount of steps, and the amount of distinct values per variable
replicates = 10
max_steps = 100
distinct_samples = 8

model_reporters = {'Entropy': lambda m: m.calc_entropy(),
                   'Happy agents': lambda m: m.happiness()}

# get samples
param_values = saltelli.sample(problem, distinct_samples, calc_second_order = False)

batch = BatchRunner(model,
                    max_steps=max_steps,
                    variable_parameters={name:[] for name in problem['names']},
                    model_reporters=model_reporters)

count = 0
data = pd.DataFrame(index=range(replicates*len(param_values)),
                                columns=['density', 'init_neutral', 'similar_wanted', 'radius', 'size', 'network_p',
                                 'randomize_part', 'decrease_intolerance'])
data['Run'], data['Entropy'], data['Happy agents'] = None, None, None

for i in range(replicates):
    for vals in param_values:
        variable_parameters = {}

        for name, val in zip(problem['names'], vals):
            if name == 'size' or name == 'radius':
                val = int(val)
            variable_parameters[name] = val

        batch.run_iteration(variable_parameters, tuple(vals), count)
        iteration_data = batch.get_model_vars_dataframe().iloc[count]
        iteration_data['Run'] = count
        data.iloc[count, 0:8] = vals
        data.iloc[count, 8:11] = iteration_data
        count += 1

        for i in range(10, 101, 10):
            if count / (len(param_values) * (replicates)) * 100 == i:
                print(f'{i}% done!')

print(data)

# unhash to retrieve data as csv file
# data.to_csv('GlobalSA_data.csv')

Si_happy_agents = sobol.analyze(problem, data['Happy agents'].values, calc_second_order=False, print_to_console=True)
Si_entropy = sobol.analyze(problem, data['Entropy'].values, calc_second_order=False, print_to_console=True)

def plot_index(s, params, i, title=''):
    """
    Creates a plot for Sobol sensitivity analysis that shows the contributions
    of each parameter to the global sensitivity.

    Args:
        s (dict): dictionary {'S#': dict, 'S#_conf': dict} of dicts that hold
            the values for a set of parameters
        params (list): the parameters taken from s
        i (str): string that indicates what order the sensitivity is.
        title (str): title for the plot
    """

    if i == '2':
        p = len(params)
        params = list(combinations(params, 2))
        indices = s['S' + i].reshape((p ** 2))
        indices = indices[~np.isnan(indices)]
        errors = s['S' + i + '_conf'].reshape((p ** 2))
        errors = errors[~np.isnan(errors)]
    else:
        indices = s['S' + i]
        errors = s['S' + i + '_conf']
        plt.figure()

    l = len(indices)

    plt.title(title)
    plt.ylim([-0.2, len(indices) - 1 + 0.2])
    plt.yticks(range(l), params)
    plt.errorbar(indices, range(l), xerr=errors, linestyle='None', marker='o')
    plt.axvline(0, c='k')


order_labels = ['1', 'T']
title_labels = ['First order sensitivity', 'Total order sensitivity']
Si_labels = ['Happy agents', 'Entropy']

# very simple code to save the sensivity analysis plots to the desired working directory
# save_results_to = 'C:/Users/ysijp/OneDrive/Bureaublad/Agent-Based Modelling/GroupProject/Figures/'

for i in range(len(order_labels)):
    plot_index(Si_happy_agents, problem['names'], order_labels[i], title_labels[i] + " " "(Happy agents)")
    # plt.savefig(save_results_to + title_labels[i] + "_Happy agents" + '.png', bbox_inches="tight",  dpi = 300)
    plt.show()

for i in range(len(order_labels)):
    plot_index(Si_entropy, problem['names'], order_labels[i], title_labels[i] + " " "(Entropy)")
    # plt.savefig(save_results_to + title_labels[i] + "_Entropy" + '.png', bbox_inches="tight", dpi = 300)
    plt.show()