Skip to content

Monte Carlo

machwave.montecarlo

MonteCarloSimulation

The MonteCarloSimulation class: - Stores data for a Monte Carlo simulation - Executes the simulation - Delegates plotting to plot_service.py

Source code in machwave/montecarlo/base.py 
 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
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
class MonteCarloSimulation:
    """
    The MonteCarloSimulation class:
    - Stores data for a Monte Carlo simulation
    - Executes the simulation
    - Delegates plotting to `plot_service.py`
    """

    def __init__(
        self,
        parameters: list[typing.Any],
        number_of_scenarios: int,
        simulation: type[Simulation],
    ) -> None:
        """
        Initializes a MonteCarloSimulation object.

        Args:
            parameters: list with the input parameters for a simulation class.
            number_of_scenarios: Number of scenarios to be simulated.
            simulation: Simulation class reference.
        """
        self.parameters = parameters
        self.number_of_scenarios = number_of_scenarios
        self.simulation = simulation

        self.scenarios: list = []
        self.results: list = []

        self._object_store = dict()  # for nested parameter processing

    def generate_scenario(self) -> list[typing.Any]:
        """
        Generates a Monte Carlo scenario in the form of a list of parameters.
        """
        new_scenario = []
        parameters_copy = deepcopy(self.parameters)

        for parameter in parameters_copy:
            if isinstance(parameter, MonteCarloParameter):
                parameter = parameter.get_random_value()
            else:
                self._process_nested_parameters(parameter)

            new_scenario.append(parameter)

        self.scenarios.append(new_scenario)
        return new_scenario

    def _process_nested_parameters(self, parameter: typing.Any) -> None:
        """
        Recursively processes an object's attributes to replace
        MonteCarloParameter instances with randomized values and store objects
        using UUIDs.

        Args:
            parameter: The object whose attributes will be processed.
        """
        parameter_uuid = uuid.uuid4()
        self._object_store[parameter_uuid] = parameter
        search_tree = {parameter_uuid: get_object_dict(parameter)}

        i = 0
        while search_tree and i < SEARCH_TREE_DEPTH_LIMIT:
            i += 1
            new_search_tree = {}
            for param_uuid, sub_params in search_tree.items():
                param = self._object_store[param_uuid]
                for name, attr in sub_params.items():
                    object_uuid = uuid.uuid4()

                    if isinstance(attr, MonteCarloParameter):
                        setattr(param, name, attr.get_random_value())
                    elif isinstance(attr, list):
                        for item in attr:
                            if isinstance(item, dict):
                                continue
                            self._object_store[object_uuid] = item
                            new_search_tree[object_uuid] = get_object_dict(item)
                    else:
                        object_uuid = uuid.uuid4()
                        self._object_store[object_uuid] = attr
                        new_search_tree[object_uuid] = get_object_dict(attr)

            search_tree = new_search_tree

    def run(self) -> None:
        """
        Executes `number_of_scenarios` runs of the underlying Simulation.
        """
        self.results = []
        for _ in range(self.number_of_scenarios):
            scenario = self.generate_scenario()
            self.results.append(self.simulation(*scenario).run())

    def retrieve_values_from_result(
        self,
        state_index: int,
        property_name: str,
    ) -> np.ndarray:
        """
        Retrieves a specific scalar property from all simulation results.
        Returns a NumPy array of length = number_of_scenarios.

        Args:
            state_index: Index of the state in the simulation result.
            property_name: Name of the property to retrieve from the results.
        Returns:
            NumPy array containing the values of the specified property
            across all simulation results.
        """
        return np.array(
            [
                getattr(sim_result[state_index], property_name)
                for sim_result in self.results
            ]
        )

    def get_property_stats(
        self, state_index: int, property_name: str
    ) -> dict[str, float]:
        """
        Compute descriptive statistics for a scalar property across all results.

        Metrics returned:
        mean, median, variance, std_dev, mode, skew
        (Fishers, unbiased), kurtosis (excess, unbiased),
        p5 (5th percentile), p95 (95th percentile).

        Args:
            state_index: Index of the state in the simulation result.
            property_name: Name of the property to retrieve from the
                results.
        Returns:
            Dictionary containing the mean, median, variance, and
            standard deviation of the specified property across all
            simulation results.
        """
        values = np.asarray(
            self.retrieve_values_from_result(state_index, property_name)
        )

        mean_val = np.mean(values)
        median_val = np.median(values)
        var_val = np.var(values)
        std_val = np.std(values)

        mode_val = float(scipy_stats.mode(values, nan_policy="omit").mode)
        skew_val = scipy_stats.skew(values, bias=False)  # unbiased Fisher skew
        kurt_val = scipy_stats.kurtosis(values, fisher=True, bias=False)
        p5, p95 = np.percentile(values, [5, 95])

        return {
            "mean": float(mean_val),
            "median": float(median_val),
            "variance": float(var_val),
            "std_dev": float(std_val),
            "mode": float(mode_val),
            "skew": float(skew_val),
            "kurtosis": float(kurt_val),
            "p5": float(p5),
            "p95": float(p95),
        }

    def plot_histogram(
        self,
        state_index: int,
        property_name: str,
        x_axes_title: str = "x",
        **plotly_kwargs,
    ) -> None:
        plot_service.plot_histogram(
            self.results, state_index, property_name, x_axes_title, **plotly_kwargs
        )

    def plot_histogram_with_kde(
        self,
        state_index: int,
        property_name: str,
        x_axes_title: str = "x",
        nbins: int = 30,
        kde_points: int = 200,
        **plotly_kwargs,
    ) -> None:
        plot_service.plot_histogram_with_kde(
            self.results,
            state_index,
            property_name,
            x_axes_title,
            nbins,
            kde_points,
            **plotly_kwargs,
        )

    def plot_cdf(
        self,
        state_index: int,
        property_name: str,
        x_axes_title: str = "x",
        **plotly_kwargs,
    ) -> None:
        plot_service.plot_cdf(
            self.results, state_index, property_name, x_axes_title, **plotly_kwargs
        )

    def plot_time_series_extremes(
        self,
        state_index: int,
        time_property: str,
        series_property: str,
        x_axes_title: str = "time",
        title: str | None = None,
        **plotly_kwargs,
    ) -> None:
        plot_service.plot_time_series_extremes(
            self.results,
            state_index,
            time_property,
            series_property,
            x_axes_title,
            title,
            **plotly_kwargs,
        )

__init__(parameters, number_of_scenarios, simulation)

Initializes a MonteCarloSimulation object.

Parameters:

Name Type Description Default
parameters list[Any]

list with the input parameters for a simulation class.

 required
number_of_scenarios int

Number of scenarios to be simulated.

 required
simulation type[Simulation]

Simulation class reference.

 required
Source code in machwave/montecarlo/base.py 
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
def __init__(
    self,
    parameters: list[typing.Any],
    number_of_scenarios: int,
    simulation: type[Simulation],
) -> None:
    """
    Initializes a MonteCarloSimulation object.

    Args:
        parameters: list with the input parameters for a simulation class.
        number_of_scenarios: Number of scenarios to be simulated.
        simulation: Simulation class reference.
    """
    self.parameters = parameters
    self.number_of_scenarios = number_of_scenarios
    self.simulation = simulation

    self.scenarios: list = []
    self.results: list = []

    self._object_store = dict()  # for nested parameter processing

generate_scenario()

Generates a Monte Carlo scenario in the form of a list of parameters.

Source code in machwave/montecarlo/base.py 
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
def generate_scenario(self) -> list[typing.Any]:
    """
    Generates a Monte Carlo scenario in the form of a list of parameters.
    """
    new_scenario = []
    parameters_copy = deepcopy(self.parameters)

    for parameter in parameters_copy:
        if isinstance(parameter, MonteCarloParameter):
            parameter = parameter.get_random_value()
        else:
            self._process_nested_parameters(parameter)

        new_scenario.append(parameter)

    self.scenarios.append(new_scenario)
    return new_scenario

get_property_stats(state_index, property_name)

Compute descriptive statistics for a scalar property across all results.

Metrics returned: mean, median, variance, std_dev, mode, skew (Fishers, unbiased), kurtosis (excess, unbiased), p5 (5th percentile), p95 (95th percentile).

Parameters:

Name Type Description Default
state_index int

Index of the state in the simulation result.

 required
property_name str

Name of the property to retrieve from the results.

 required

Returns: Dictionary containing the mean, median, variance, and standard deviation of the specified property across all simulation results.

Source code in machwave/montecarlo/base.py 
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
def get_property_stats(
    self, state_index: int, property_name: str
) -> dict[str, float]:
    """
    Compute descriptive statistics for a scalar property across all results.

    Metrics returned:
    mean, median, variance, std_dev, mode, skew
    (Fishers, unbiased), kurtosis (excess, unbiased),
    p5 (5th percentile), p95 (95th percentile).

    Args:
        state_index: Index of the state in the simulation result.
        property_name: Name of the property to retrieve from the
            results.
    Returns:
        Dictionary containing the mean, median, variance, and
        standard deviation of the specified property across all
        simulation results.
    """
    values = np.asarray(
        self.retrieve_values_from_result(state_index, property_name)
    )

    mean_val = np.mean(values)
    median_val = np.median(values)
    var_val = np.var(values)
    std_val = np.std(values)

    mode_val = float(scipy_stats.mode(values, nan_policy="omit").mode)
    skew_val = scipy_stats.skew(values, bias=False)  # unbiased Fisher skew
    kurt_val = scipy_stats.kurtosis(values, fisher=True, bias=False)
    p5, p95 = np.percentile(values, [5, 95])

    return {
        "mean": float(mean_val),
        "median": float(median_val),
        "variance": float(var_val),
        "std_dev": float(std_val),
        "mode": float(mode_val),
        "skew": float(skew_val),
        "kurtosis": float(kurt_val),
        "p5": float(p5),
        "p95": float(p95),
    }

retrieve_values_from_result(state_index, property_name)

Retrieves a specific scalar property from all simulation results. Returns a NumPy array of length = number_of_scenarios.

Parameters:

Name Type Description Default
state_index int

Index of the state in the simulation result.

 required
property_name str

Name of the property to retrieve from the results.

 required

Returns: NumPy array containing the values of the specified property across all simulation results.

Source code in machwave/montecarlo/base.py 
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
def retrieve_values_from_result(
    self,
    state_index: int,
    property_name: str,
) -> np.ndarray:
    """
    Retrieves a specific scalar property from all simulation results.
    Returns a NumPy array of length = number_of_scenarios.

    Args:
        state_index: Index of the state in the simulation result.
        property_name: Name of the property to retrieve from the results.
    Returns:
        NumPy array containing the values of the specified property
        across all simulation results.
    """
    return np.array(
        [
            getattr(sim_result[state_index], property_name)
            for sim_result in self.results
        ]
    )

run()

Executes number_of_scenarios runs of the underlying Simulation.

Source code in machwave/montecarlo/base.py 
130
131
132
133
134
135
136
137
def run(self) -> None:
    """
    Executes `number_of_scenarios` runs of the underlying Simulation.
    """
    self.results = []
    for _ in range(self.number_of_scenarios):
        scenario = self.generate_scenario()
        self.results.append(self.simulation(*scenario).run())