.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "auto_reliability_sensitivity/sensitivity_analysis/plot_hsic_estimators_ishigami.py"
.. LINE NUMBERS ARE GIVEN BELOW.
.. only:: html
.. note::
:class: sphx-glr-download-link-note
:ref:`Go to the end <sphx_glr_download_auto_reliability_sensitivity_sensitivity_analysis_plot_hsic_estimators_ishigami.py>`
to download the full example code.
.. rst-class:: sphx-glr-example-title
.. _sphx_glr_auto_reliability_sensitivity_sensitivity_analysis_plot_hsic_estimators_ishigami.py:
The HSIC sensitivity indices: the Ishigami model
================================================
.. GENERATED FROM PYTHON SOURCE LINES 5-11
.. code-block:: Python
import openturns as ot
import openturns.viewer as otv
from openturns.usecases import ishigami_function
.. GENERATED FROM PYTHON SOURCE LINES 12-15
This example is a brief overview of the HSIC sensitivity indices classes and how to call them.
HSIC estimators rely on a reproducing kernel of a Hilbert space. We can use them to compute sensitivity
indices. We present the methods on the :ref:`Ishigami function<use-case-ishigami>`.
.. GENERATED FROM PYTHON SOURCE LINES 18-22
Definition of the model
-----------------------
We load the model from the usecases module.
.. GENERATED FROM PYTHON SOURCE LINES 22-24
.. code-block:: Python
im = ishigami_function.IshigamiModel()
.. GENERATED FROM PYTHON SOURCE LINES 25-26
We generate an input sample of size 100 (and dimension 3).
.. GENERATED FROM PYTHON SOURCE LINES 26-29
.. code-block:: Python
size = 100
X = im.inputDistribution.getSample(size)
.. GENERATED FROM PYTHON SOURCE LINES 30-31
We compute the output by applying the Ishigami model to the input sample.
.. GENERATED FROM PYTHON SOURCE LINES 31-34
.. code-block:: Python
Y = im.model(X)
.. GENERATED FROM PYTHON SOURCE LINES 35-45
Setting the covariance models
-----------------------------
The HSIC algorithms use reproducing kernels defined on Hilbert spaces to estimate independence.
For each input variable we choose a covariance kernel.
Here we choose a :class:`~openturns.SquaredExponential`
kernel for all input variables.
They are all stored in a list of :math:`d+1` covariance kernels where :math:`d` is the number of
input variables. The remaining one is for the output variable.
.. GENERATED FROM PYTHON SOURCE LINES 45-47
.. code-block:: Python
covarianceModelCollection = []
.. GENERATED FROM PYTHON SOURCE LINES 48-54
.. code-block:: Python
for i in range(3):
Xi = X.getMarginal(i)
inputCovariance = ot.SquaredExponential(1)
inputCovariance.setScale(Xi.computeStandardDeviation())
covarianceModelCollection.append(inputCovariance)
.. GENERATED FROM PYTHON SOURCE LINES 55-56
Likewise we define a covariance kernel associated to the output variable.
.. GENERATED FROM PYTHON SOURCE LINES 56-60
.. code-block:: Python
outputCovariance = ot.SquaredExponential(1)
outputCovariance.setScale(Y.computeStandardDeviation())
covarianceModelCollection.append(outputCovariance)
.. GENERATED FROM PYTHON SOURCE LINES 61-66
The Global HSIC estimator
-------------------------
In this paragraph, we perform the analysis on the raw data: that is
the global HSIC estimator.
.. GENERATED FROM PYTHON SOURCE LINES 68-80
Choosing an estimator
^^^^^^^^^^^^^^^^^^^^^
After having defined the covariance kernels one has to select an
appropriate estimator for the computations.
Two estimators are proposed:
- an unbiased estimator through the :class:`~openturns.HSICUStat` class
- a biased, but asymptotically unbiased, estimator through the :class:`~openturns.HSICVStat` class
Beware that the conditional analysis used later cannot be performed with the unbiased estimator.
.. GENERATED FROM PYTHON SOURCE LINES 80-83
.. code-block:: Python
estimatorType = ot.HSICUStat()
.. GENERATED FROM PYTHON SOURCE LINES 84-85
We now build the HSIC estimator:
.. GENERATED FROM PYTHON SOURCE LINES 85-89
.. code-block:: Python
globHSIC = ot.HSICEstimatorGlobalSensitivity(
covarianceModelCollection, X, Y, estimatorType
)
.. GENERATED FROM PYTHON SOURCE LINES 90-91
We get the R2-HSIC indices:
.. GENERATED FROM PYTHON SOURCE LINES 91-95
.. code-block:: Python
R2HSICIndices = globHSIC.getR2HSICIndices()
print("\n Global HSIC analysis")
print("R2-HSIC Indices: ", R2HSICIndices)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
Global HSIC analysis
R2-HSIC Indices: [0.0937949,0.00487779,0.062993]
.. GENERATED FROM PYTHON SOURCE LINES 96-97
and the HSIC indices:
.. GENERATED FROM PYTHON SOURCE LINES 97-100
.. code-block:: Python
HSICIndices = globHSIC.getHSICIndices()
print("HSIC Indices: ", HSICIndices)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
HSIC Indices: [0.00809327,0.000415869,0.00514516]
.. GENERATED FROM PYTHON SOURCE LINES 101-102
The p-value by permutation.
.. GENERATED FROM PYTHON SOURCE LINES 102-105
.. code-block:: Python
pvperm = globHSIC.getPValuesPermutation()
print("p-value (permutation): ", pvperm)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
p-value (permutation): [0,0.336634,0.00990099]
.. GENERATED FROM PYTHON SOURCE LINES 106-107
We have an asymptotic estimate of the value for this estimator.
.. GENERATED FROM PYTHON SOURCE LINES 107-110
.. code-block:: Python
pvas = globHSIC.getPValuesAsymptotic()
print("p-value (asymptotic): ", pvas)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
p-value (asymptotic): [0.000204667,0.292796,0.00231475]
.. GENERATED FROM PYTHON SOURCE LINES 111-112
We vizualise the results.
.. GENERATED FROM PYTHON SOURCE LINES 112-125
.. code-block:: Python
graph1 = globHSIC.drawHSICIndices()
view1 = otv.View(graph1)
graph2 = globHSIC.drawPValuesAsymptotic()
view2 = otv.View(graph2)
graph3 = globHSIC.drawR2HSICIndices()
view3 = otv.View(graph3)
graph4 = globHSIC.drawPValuesPermutation()
view4 = otv.View(graph4)
.. rst-class:: sphx-glr-horizontal
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_001.png
:alt: HSIC indices
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_001.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_002.png
:alt: Asymptotic p-values
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_002.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_003.png
:alt: R2-HSIC indices
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_003.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_004.png
:alt: p-values by permutation
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_004.png
:class: sphx-glr-multi-img
.. GENERATED FROM PYTHON SOURCE LINES 126-131
The Target HSIC estimator
-------------------------
We now perform the target analysis which consists in using a filter function over the
output.
.. GENERATED FROM PYTHON SOURCE LINES 134-140
Defining a filter function
^^^^^^^^^^^^^^^^^^^^^^^^^^
We define a filter function on the output variable for the target
analysis. In this example we use the function :math:`\exp{(-d/s)}` where :math:`d` is the distance
to a well-chosen interval.
.. GENERATED FROM PYTHON SOURCE LINES 142-143
We first define a critical domain: in our case that is the :math:`[5,+\infty[` interval.
.. GENERATED FROM PYTHON SOURCE LINES 143-145
.. code-block:: Python
criticalDomain = ot.Interval(5, float("inf"))
.. GENERATED FROM PYTHON SOURCE LINES 146-148
We have access to the distance to this domain thanks to the
:class:`~openturns.DistanceToDomainFunction` class.
.. GENERATED FROM PYTHON SOURCE LINES 148-150
.. code-block:: Python
dist2criticalDomain = ot.DistanceToDomainFunction(criticalDomain)
.. GENERATED FROM PYTHON SOURCE LINES 151-152
We define the empirical parameter values in our function from the output sample
.. GENERATED FROM PYTHON SOURCE LINES 152-154
.. code-block:: Python
s = 0.1 * Y.computeStandardDeviation()[0]
.. GENERATED FROM PYTHON SOURCE LINES 155-157
We now define our filter function by composition of the parametrized function and
the distance function.
.. GENERATED FROM PYTHON SOURCE LINES 157-161
.. code-block:: Python
f = ot.SymbolicFunction(["x", "s"], ["exp(-x/s)"])
phi = ot.ParametricFunction(f, [1], [s])
filterFunction = ot.ComposedFunction(phi, dist2criticalDomain)
.. GENERATED FROM PYTHON SOURCE LINES 162-163
We modify the output covariance kernel so as to adapt it to the filtered output
.. GENERATED FROM PYTHON SOURCE LINES 163-167
.. code-block:: Python
Y_filtered = filterFunction(Y)
outputCovariance.setScale(Y_filtered.computeStandardDeviation())
covarianceModelCollection[-1] = outputCovariance
.. GENERATED FROM PYTHON SOURCE LINES 168-169
We choose an unbiased estimator
.. GENERATED FROM PYTHON SOURCE LINES 169-171
.. code-block:: Python
estimatorType = ot.HSICUStat()
.. GENERATED FROM PYTHON SOURCE LINES 172-173
and build the HSIC estimator
.. GENERATED FROM PYTHON SOURCE LINES 173-177
.. code-block:: Python
targetHSIC = ot.HSICEstimatorTargetSensitivity(
covarianceModelCollection, X, Y, estimatorType, filterFunction
)
.. GENERATED FROM PYTHON SOURCE LINES 178-179
We get the R2-HSIC indices:
.. GENERATED FROM PYTHON SOURCE LINES 179-183
.. code-block:: Python
R2HSICIndices = targetHSIC.getR2HSICIndices()
print("\n Target HSIC analysis")
print("R2-HSIC Indices: ", R2HSICIndices)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
Target HSIC analysis
R2-HSIC Indices: [0.0650641,-0.000717287,0.00854149]
.. GENERATED FROM PYTHON SOURCE LINES 184-185
and the HSIC indices:
.. GENERATED FROM PYTHON SOURCE LINES 185-188
.. code-block:: Python
HSICIndices = targetHSIC.getHSICIndices()
print("HSIC Indices: ", HSICIndices)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
HSIC Indices: [0.00908026,-9.89096e-05,0.00112837]
.. GENERATED FROM PYTHON SOURCE LINES 189-190
The p-value by permutation.
.. GENERATED FROM PYTHON SOURCE LINES 190-193
.. code-block:: Python
pvperm = targetHSIC.getPValuesPermutation()
print("p-value (permutation): ", pvperm)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
p-value (permutation): [0.019802,0.326733,0.227723]
.. GENERATED FROM PYTHON SOURCE LINES 194-195
We have an asymptotic estimate of the value for this estimator.
.. GENERATED FROM PYTHON SOURCE LINES 195-198
.. code-block:: Python
pvas = targetHSIC.getPValuesAsymptotic()
print("p-value (asymptotic): ", pvas)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
p-value (asymptotic): [0.00296907,0.392551,0.201913]
.. GENERATED FROM PYTHON SOURCE LINES 199-200
We vizualise the results.
.. GENERATED FROM PYTHON SOURCE LINES 200-213
.. code-block:: Python
graph5 = targetHSIC.drawHSICIndices()
view5 = otv.View(graph5)
graph6 = targetHSIC.drawPValuesAsymptotic()
view6 = otv.View(graph6)
graph7 = targetHSIC.drawR2HSICIndices()
view7 = otv.View(graph7)
graph8 = targetHSIC.drawPValuesPermutation()
view8 = otv.View(graph8)
.. rst-class:: sphx-glr-horizontal
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_005.png
:alt: HSIC indices
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_005.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_006.png
:alt: Asymptotic p-values
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_006.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_007.png
:alt: R2-HSIC indices
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_007.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_008.png
:alt: p-values by permutation
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_008.png
:class: sphx-glr-multi-img
.. GENERATED FROM PYTHON SOURCE LINES 214-220
The Conditional HSIC estimator
------------------------------
In this last section we preprocess the input variables: that is the conditional
analysis. To do so, one has to work with a weight function.
Here the weight function is the filter function we used previously.
.. GENERATED FROM PYTHON SOURCE LINES 220-222
.. code-block:: Python
weightFunction = filterFunction
.. GENERATED FROM PYTHON SOURCE LINES 223-224
We revert to the covariance kernel associated to the unfiltered output
.. GENERATED FROM PYTHON SOURCE LINES 224-227
.. code-block:: Python
outputCovariance.setScale(Y.computeStandardDeviation())
covarianceModelCollection[-1] = outputCovariance
.. GENERATED FROM PYTHON SOURCE LINES 228-229
We have to select a biased -but asymptotically unbiased- estimator
.. GENERATED FROM PYTHON SOURCE LINES 229-232
.. code-block:: Python
estimatorType = ot.HSICVStat()
.. GENERATED FROM PYTHON SOURCE LINES 233-234
We build the conditional HSIC estimator
.. GENERATED FROM PYTHON SOURCE LINES 234-238
.. code-block:: Python
condHSIC = ot.HSICEstimatorConditionalSensitivity(
covarianceModelCollection, X, Y, weightFunction
)
.. GENERATED FROM PYTHON SOURCE LINES 239-240
We get the R2-HSIC indices:
.. GENERATED FROM PYTHON SOURCE LINES 240-244
.. code-block:: Python
R2HSICIndices = condHSIC.getR2HSICIndices()
print("\n Conditional HSIC analysis")
print("R2-HSIC Indices: ", R2HSICIndices)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
Conditional HSIC analysis
R2-HSIC Indices: [0.179912,0.0385186,0.0726262]
.. GENERATED FROM PYTHON SOURCE LINES 245-246
and the HSIC indices:
.. GENERATED FROM PYTHON SOURCE LINES 246-249
.. code-block:: Python
HSICIndices = condHSIC.getHSICIndices()
print("HSIC Indices: ", HSICIndices)
.. rst-class:: sphx-glr-script-out
.. code-block:: none
HSIC Indices: [0.00710788,0.00183856,0.00293997]
.. GENERATED FROM PYTHON SOURCE LINES 250-251
For the conditional estimator we only have access to the p-value by permutation:
.. GENERATED FROM PYTHON SOURCE LINES 251-255
.. code-block:: Python
pvperm = condHSIC.getPValuesPermutation()
print("p-value (permutation): ", pvperm)
print("")
.. rst-class:: sphx-glr-script-out
.. code-block:: none
p-value (permutation): [0,0.188119,0.0594059]
.. GENERATED FROM PYTHON SOURCE LINES 256-257
We can vizualise the results thanks to the various drawing methods.
.. GENERATED FROM PYTHON SOURCE LINES 257-267
.. code-block:: Python
graph9 = condHSIC.drawHSICIndices()
view9 = otv.View(graph9)
graph10 = condHSIC.drawR2HSICIndices()
view10 = otv.View(graph10)
graph11 = condHSIC.drawPValuesPermutation()
view11 = otv.View(graph11)
.. rst-class:: sphx-glr-horizontal
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_009.png
:alt: HSIC indices
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_009.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_010.png
:alt: R2-HSIC indices
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_010.png
:class: sphx-glr-multi-img
*
.. image-sg:: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_011.png
:alt: p-values by permutation
:srcset: /auto_reliability_sensitivity/sensitivity_analysis/images/sphx_glr_plot_hsic_estimators_ishigami_011.png
:class: sphx-glr-multi-img
.. GENERATED FROM PYTHON SOURCE LINES 268-269
.. code-block:: Python
otv.View.ShowAll()
.. _sphx_glr_download_auto_reliability_sensitivity_sensitivity_analysis_plot_hsic_estimators_ishigami.py:
.. only:: html
.. container:: sphx-glr-footer sphx-glr-footer-example
.. container:: sphx-glr-download sphx-glr-download-jupyter
:download:`Download Jupyter notebook: plot_hsic_estimators_ishigami.ipynb <plot_hsic_estimators_ishigami.ipynb>`
.. container:: sphx-glr-download sphx-glr-download-python
:download:`Download Python source code: plot_hsic_estimators_ishigami.py <plot_hsic_estimators_ishigami.py>`
.. container:: sphx-glr-download sphx-glr-download-zip
:download:`Download zipped: plot_hsic_estimators_ishigami.zip <plot_hsic_estimators_ishigami.zip>`