Coverage for aixcalibuha/sensitivity_analysis/sobol.py: 98%
82 statements
« prev ^ index » next coverage.py v7.4.1, created at 2024-01-27 10:48 +0000
« prev ^ index » next coverage.py v7.4.1, created at 2024-01-27 10:48 +0000
1"""
2Adds the SobolAnalyzer to the available
3classes of sensitivity analysis.
4"""
5import pandas as pd
6from SALib.sample import sobol
7from SALib.analyze import sobol as analyze_sobol
8import numpy as np
9from aixcalibuha.sensitivity_analysis import SenAnalyzer
10from aixcalibuha import CalibrationClass
11from aixcalibuha.sensitivity_analysis.plotting import plot_single, heatmaps
14class SobolAnalyzer(SenAnalyzer):
15 """
16 Sobol method from SALib
17 https://salib.readthedocs.io/en/latest/api.html#sobol-sensitivity-analysis
18 A variance-based method which can compute the sensitivity measures
19 'S1', 'ST' and 'S2' with their confidence intervals.
21 Additional arguments:
23 :keyword bool calc_second_order:
24 Default True, used for the sobol-method
25 :keyword int seed:
26 Used for the sobol-method
27 """
28 __analysis_variables = ['S1', 'ST', 'S1_conf', 'ST_conf']
29 __analysis_variables_1 = ['S1', 'ST', 'S1_conf', 'ST_conf']
30 __analysis_variables_2 = ['S2', 'S2_conf']
32 def __init__(self, sim_api, **kwargs):
33 # Additional kwarg which changes the possible analysis_variables
34 self.calc_second_order = kwargs.get("calc_second_order", True)
35 if self.calc_second_order:
36 self.__analysis_variables = ['S1', 'ST', 'S1_conf', 'ST_conf', 'S2', 'S2_conf']
37 # separately store first and total order (1) and second order (2) analysis variables
38 self.av_1_selected = []
39 self.av_2_selected = []
40 # Set additional kwargs
41 self.seed = kwargs.pop("seed", None)
43 super().__init__(
44 sim_api=sim_api,
45 **kwargs)
47 @property
48 def analysis_variables(self):
49 """The analysis variables of the sobol method"""
50 return self.__analysis_variables
52 def analysis_function(self, x, y):
53 """
54 Use the SALib.analyze.sobol method to analyze the simulation results.
56 :param np.array x:
57 placeholder for the `X` parameter of the morris method not used for sobol
58 :param np.array y:
59 The NumPy array containing the model outputs
60 :return:
61 returns the result of the SALib.analyze.sobol method (from the documentation:
62 a dictionary with cols `S1`, `S1_conf`, `ST`, and `ST_conf`, where each entry
63 is a list of size D (the number of parameters) containing the indices in the same
64 order as the parameter file. If calc_second_order is True, the dictionary also
65 contains cols `S2` and `S2_conf`.)
66 """
67 return analyze_sobol.analyze(self.problem, y,
68 calc_second_order=self.calc_second_order)
70 def create_sampler_demand(self):
71 """
72 Function to create the sampler parameters for the morris method
73 """
74 return {'calc_second_order': self.calc_second_order}
76 def generate_samples(self):
77 """
78 Run the sampler for sobol and return the results.
80 :return:
81 The list of samples generated as a NumPy array with one row per sample
82 and each row containing one value for each variable name in `problem['names']`.
83 :rtype: np.ndarray
84 """
85 return sobol.sample(self.problem,
86 N=self.num_samples,
87 **self.create_sampler_demand())
89 def _save(self, result: tuple, time_dependent: bool = False):
90 """
91 Save the results of the run and run_time_dependent function of the SobolAnalyzer.
92 """
93 if not result[0].empty:
94 super()._save(result=result[0], time_dependent=time_dependent)
95 if time_dependent:
96 savepath_result_2 = self.cd.joinpath(
97 f'{self.__class__.__name__}_results_second_order_time.csv')
98 else:
99 savepath_result_2 = self.cd.joinpath(
100 f'{self.__class__.__name__}_results_second_order.csv')
101 if not result[1].empty:
102 result[1].to_csv(savepath_result_2)
103 self.reproduction_files.append(savepath_result_2)
105 def _conv_local_results(self, results: list, local_classes: list):
106 """
107 Convert the result dictionaries form SALib
108 of each class and goal into a tuple of two DataFrames.
109 First is the single order and second is the second order result.
110 If one of the results is not computed an empty list is returned.
111 """
112 _conv_results = []
113 _conv_results_2 = []
114 tuples = []
115 tuples_2 = []
116 for class_results, local_class in zip(results, local_classes):
117 for goal, goal_results in class_results.items():
118 for analysis_var in self.analysis_variables:
119 res_dict = self._get_res_dict(result=goal_results,
120 cal_class=local_class,
121 analysis_variable=analysis_var)
122 if analysis_var in self.__analysis_variables_1:
123 _conv_results.append(res_dict)
124 tuples.append((local_class.name, goal, analysis_var))
125 elif analysis_var in self.__analysis_variables_2:
126 for tuner_para, res_dict in res_dict.items():
127 _conv_results_2.append(res_dict)
128 tuples_2.append((local_class.name, goal, analysis_var, tuner_para))
129 index = pd.MultiIndex.from_tuples(tuples=tuples,
130 names=['Class', 'Goal', 'Analysis variable'])
131 index_2 = pd.MultiIndex.from_tuples(tuples=tuples_2,
132 names=['Class', 'Goal', 'Analysis variable',
133 'Interaction'])
134 df = pd.DataFrame(_conv_results, index=index)
135 df_2 = pd.DataFrame(_conv_results_2, index=index_2)
136 return df, df_2
138 def _get_res_dict(self, result: dict, cal_class: CalibrationClass, analysis_variable: str):
139 """
140 Convert the result object to a dict with the key
141 being the variable name and the value being the result
142 associated to analysis_variable.
143 For second oder analysis variables the result is converted to a
144 dict with the key being the variable name and the value being another dict
145 with the variable names as the keys and the result associated to analysis_variable
146 from the interaction between the two variables.
147 """
148 names = self.create_problem(cal_class.tuner_paras)['names']
149 if analysis_variable in self.__analysis_variables_1:
150 if result is None:
151 res_dict_1 = {var_name: np.abs(res_val)
152 for var_name, res_val in zip(names,
153 np.zeros(len(names)))}
154 else:
155 res_dict_1 = {var_name: np.abs(res_val)
156 for var_name, res_val in zip(names,
157 result[analysis_variable])}
158 return res_dict_1
159 if analysis_variable in self.__analysis_variables_2:
160 if result is None:
161 res_dict_2 = {var_name: dict(zip(names, np.abs(res_val)))
162 for var_name, res_val in zip(names,
163 np.zeros((len(names), len(names))))}
164 else:
165 result_av = result[analysis_variable]
166 for i, _ in enumerate(result_av):
167 for j, _ in enumerate(result_av):
168 if i > j:
169 result_av[i][j] = result_av[j][i]
170 res_dict_2 = {var_name: dict(zip(names, np.abs(res_val)))
171 for var_name, res_val in zip(names,
172 result_av)}
173 return res_dict_2
175 def plot(self, result):
176 """
177 Plot the results of the sensitivity analysis method from run().
179 :param pd.DataFrame result:
180 Dataframe of the results like from the run() function.
181 :return tuple of matplotlib objects (fig, ax):
182 """
183 plot_single(result=result[0])
184 heatmaps(result=result[1])
186 @staticmethod
187 def load_second_order_from_csv(path):
188 """
189 Load second order sensitivity results which were saved with the run() or
190 run_time_dependent() function.
191 """
192 result = pd.read_csv(path, index_col=[0, 1, 2, 3])
193 return result