Sensitivity Analysis#

This is Example 2 of Bourinet (2017), which gives analytical solutions as:

\[\begin{split}\begin{align} \frac{\partial\beta}{\partial\mu_R} &= 0.5184 \\ \frac{\partial\beta}{\partial\sigma_R} &= −0.2548 \\ \frac{\partial\beta}{\partial\mu_S} &= −1.3629 \\ \frac{\partial\beta}{\partial\sigma_S} &= 0.0445 \\ \end{align}\end{split}\]

See: Bourinet (2017), FORM Sensitivities to Distribution Parameters with the Nataf Transformation, P. Gardoni (ed.), Risk and Reliability Analysis: Theory and Applications, Springer Series in Reliability Engineering, DOI 10.1007/978-3-319-52425-2_12

Import pystra as usual:

[1]:
import pystra as pr

Define the limit state function

[2]:
def lsf(R, S):
    return R - S

And create the relevant Pystra objects:

[3]:
limit_state = pr.LimitState(lsf)

model = pr.StochasticModel()
model.addVariable(pr.Lognormal("R", 5, 5))
model.addVariable(pr.Lognormal("S", 1, 1))
model.setCorrelation(pr.CorrelationMatrix([[1.0, 0.5], [0.5, 1.0]]))

Here we suppress the printed output from the analyses, as the sensitivity analysis runs FORM many times.

[4]:
options = pr.AnalysisOptions()
options.setPrintOutput(False)

Do a FORM analysis so we have a record of the base result.

[5]:
form = pr.Form(stochastic_model=model, limit_state=limit_state, analysis_options=options)
form.run()
form.showDetailedOutput()

======================================================
FORM
======================================================
Pf                       1.6927788976e-02
BetaHL                   2.1217876054
Model Evaluations        66
------------------------------------------------------
Variable            U_star             X_star        alpha
R                -0.966683           1.580985    -0.455543
S                 1.888784           1.580984    +0.890214
======================================================

Finally, we run the sensitivity analysis to find the sensitivity of the FORM result to the parameters of the stochastic model:

[6]:
sens = pr.SensitivityAnalysis(stochastic_model=model, limit_state=limit_state, analysis_options=options)
results = sens.run_form()
[7]:
print(results)
{'R': {'mean': 0.5155641831255604, 'std': -0.2545723002165454}, 'S': {'mean': -1.3623019484605414, 'std': 0.04167433447075417}}