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

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 

12 

13 

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. 

20 

21 Additional arguments: 

22 

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'] 

31 

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) 

42 

43 super().__init__( 

44 sim_api=sim_api, 

45 **kwargs) 

46 

47 @property 

48 def analysis_variables(self): 

49 """The analysis variables of the sobol method""" 

50 return self.__analysis_variables 

51 

52 def analysis_function(self, x, y): 

53 """ 

54 Use the SALib.analyze.sobol method to analyze the simulation results. 

55 

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) 

69 

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} 

75 

76 def generate_samples(self): 

77 """ 

78 Run the sampler for sobol and return the results. 

79 

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()) 

88 

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) 

104 

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 

137 

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 

174 

175 def plot(self, result): 

176 """ 

177 Plot the results of the sensitivity analysis method from run(). 

178 

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]) 

185 

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