.. 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 Click :ref:`here ` 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-10 .. code-block:: default import openturns as ot import openturns.viewer as otv from openturns.usecases import ishigami_function .. GENERATED FROM PYTHON SOURCE LINES 11-14 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`. .. 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:: default 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:: default size = 100 X = im.distributionX.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:: default 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:: default covarianceModelCollection = [] .. GENERATED FROM PYTHON SOURCE LINES 48-54 .. code-block:: default 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:: default 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:: default estimatorType = ot.HSICUStat() .. GENERATED FROM PYTHON SOURCE LINES 84-85 We now build the HSIC estimator: .. GENERATED FROM PYTHON SOURCE LINES 85-88 .. code-block:: default globHSIC = ot.HSICEstimatorGlobalSensitivity( covarianceModelCollection, X, Y, estimatorType) .. GENERATED FROM PYTHON SOURCE LINES 89-90 We get the R2-HSIC indices: .. GENERATED FROM PYTHON SOURCE LINES 90-94 .. code-block:: default R2HSICIndices = globHSIC.getR2HSICIndices() print("\n Global HSIC analysis") print("R2-HSIC Indices: ", R2HSICIndices) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Global HSIC analysis R2-HSIC Indices: [0.0639736,-0.000188039,0.0140331] .. GENERATED FROM PYTHON SOURCE LINES 95-96 and the HSIC indices: .. GENERATED FROM PYTHON SOURCE LINES 96-99 .. code-block:: default HSICIndices = globHSIC.getHSICIndices() print("HSIC Indices: ", HSICIndices) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none HSIC Indices: [0.0051169,-1.57797e-05,0.00108955] .. GENERATED FROM PYTHON SOURCE LINES 100-101 The p-value by permutation. .. GENERATED FROM PYTHON SOURCE LINES 101-104 .. code-block:: default pvperm = globHSIC.getPValuesPermutation() print("p-value (permutation): ", pvperm) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none p-value (permutation): [0,0.356436,0.148515] .. GENERATED FROM PYTHON SOURCE LINES 105-106 We have an asymptotic estimate of the value for this estimator. .. GENERATED FROM PYTHON SOURCE LINES 106-109 .. code-block:: default pvas = globHSIC.getPValuesAsymptotic() print("p-value (asymptotic): ", pvas) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none p-value (asymptotic): [0.00207356,0.420886,0.156227] .. GENERATED FROM PYTHON SOURCE LINES 110-111 We vizualise the results. .. GENERATED FROM PYTHON SOURCE LINES 111-124 .. code-block:: default 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 125-130 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 133-139 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 141-142 We first define a critical domain: in our case that is the :math:`[5,+\infty[` interval. .. GENERATED FROM PYTHON SOURCE LINES 142-144 .. code-block:: default criticalDomain = ot.Interval(5, float('inf')) .. GENERATED FROM PYTHON SOURCE LINES 145-147 We have access to the distance to this domain thanks to the :class:`~openturns.DistanceToDomainFunction` class. .. GENERATED FROM PYTHON SOURCE LINES 147-149 .. code-block:: default dist2criticalDomain = ot.DistanceToDomainFunction(criticalDomain) .. GENERATED FROM PYTHON SOURCE LINES 150-151 We define the parameters in our function from the output sample .. GENERATED FROM PYTHON SOURCE LINES 151-153 .. code-block:: default s = 0.1 * Y.computeStandardDeviation()[0] .. GENERATED FROM PYTHON SOURCE LINES 154-156 We now define our filter function by composition of the parametrized function and the distance function. .. GENERATED FROM PYTHON SOURCE LINES 156-161 .. code-block:: default 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 choose an unbiased estimator .. GENERATED FROM PYTHON SOURCE LINES 163-165 .. code-block:: default estimatorType = ot.HSICUStat() .. GENERATED FROM PYTHON SOURCE LINES 166-167 and build the HSIC estimator .. GENERATED FROM PYTHON SOURCE LINES 167-171 .. code-block:: default targetHSIC = ot.HSICEstimatorTargetSensitivity( covarianceModelCollection, X, Y, estimatorType, filterFunction ) .. GENERATED FROM PYTHON SOURCE LINES 172-173 We get the R2-HSIC indices: .. GENERATED FROM PYTHON SOURCE LINES 173-177 .. code-block:: default R2HSICIndices = targetHSIC.getR2HSICIndices() print("\n Target HSIC analysis") print("R2-HSIC Indices: ", R2HSICIndices) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Target HSIC analysis R2-HSIC Indices: [0.0321494,0.00293329,-0.00627235] .. GENERATED FROM PYTHON SOURCE LINES 178-179 and the HSIC indices: .. GENERATED FROM PYTHON SOURCE LINES 179-182 .. code-block:: default HSICIndices = targetHSIC.getHSICIndices() print("HSIC Indices: ", HSICIndices) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none HSIC Indices: [0.000192388,1.84164e-05,-3.64356e-05] .. GENERATED FROM PYTHON SOURCE LINES 183-184 The p-value by permutation. .. GENERATED FROM PYTHON SOURCE LINES 184-187 .. code-block:: default pvperm = targetHSIC.getPValuesPermutation() print("p-value (permutation): ", pvperm) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none p-value (permutation): [0,0.257426,0.584158] .. GENERATED FROM PYTHON SOURCE LINES 188-189 We have an asymptotic estimate of the value for this estimator. .. GENERATED FROM PYTHON SOURCE LINES 189-192 .. code-block:: default pvas = targetHSIC.getPValuesAsymptotic() print("p-value (asymptotic): ", pvas) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none p-value (asymptotic): [0.0339153,0.297131,0.589994] .. GENERATED FROM PYTHON SOURCE LINES 193-194 We vizualise the results. .. GENERATED FROM PYTHON SOURCE LINES 194-207 .. code-block:: default 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 208-214 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 214-216 .. code-block:: default weightFunction = filterFunction .. GENERATED FROM PYTHON SOURCE LINES 217-218 We have to select a biased -but asymptotically unbiased- estimator .. GENERATED FROM PYTHON SOURCE LINES 218-221 .. code-block:: default estimatorType = ot.HSICVStat() .. GENERATED FROM PYTHON SOURCE LINES 222-223 We build the conditional HSIC estimator .. GENERATED FROM PYTHON SOURCE LINES 223-227 .. code-block:: default condHSIC = ot.HSICEstimatorConditionalSensitivity( covarianceModelCollection, X, Y, estimatorType, weightFunction ) .. GENERATED FROM PYTHON SOURCE LINES 228-229 We get the R2-HSIC indices: .. GENERATED FROM PYTHON SOURCE LINES 229-233 .. code-block:: default R2HSICIndices = condHSIC.getR2HSICIndices() print("\n Conditional HSIC analysis") print("R2-HSIC Indices: ", R2HSICIndices) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Conditional HSIC analysis R2-HSIC Indices: [0.250897,0.00916689,0.229481] .. GENERATED FROM PYTHON SOURCE LINES 234-235 and the HSIC indices: .. GENERATED FROM PYTHON SOURCE LINES 235-238 .. code-block:: default HSICIndices = condHSIC.getHSICIndices() print("HSIC Indices: ", HSICIndices) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none HSIC Indices: [0.0121911,0.000544078,0.0117331] .. GENERATED FROM PYTHON SOURCE LINES 239-240 For the conditional estimator we only have access to the p-value by permutation: .. GENERATED FROM PYTHON SOURCE LINES 240-244 .. code-block:: default pvperm = condHSIC.getPValuesPermutation() print("p-value (permutation): ", pvperm) print("") .. rst-class:: sphx-glr-script-out Out: .. code-block:: none p-value (permutation): [0,0.930693,0] .. GENERATED FROM PYTHON SOURCE LINES 245-246 We can vizualise the results thanks to the various drawing methods. .. GENERATED FROM PYTHON SOURCE LINES 246-256 .. code-block:: default 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 257-258 .. code-block:: default otv.View.ShowAll() .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 2.805 seconds) .. _sphx_glr_download_auto_reliability_sensitivity_sensitivity_analysis_plot_hsic_estimators_ishigami.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_hsic_estimators_ishigami.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_hsic_estimators_ishigami.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_