.. 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_probabilistic_modeling_distributions_plot_minimum_volume_level_set_1D.py: Draw minimum volume level set in 1D =================================== In this example, we compute the minimum volume level set of a univariate distribution. .. code-block:: default import openturns as ot import openturns.viewer as viewer from matplotlib import pylab as plt ot.Log.Show(ot.Log.NONE) With a Normal, minimum volume LevelSet -------------------------------------- .. code-block:: default n = ot.Normal() .. code-block:: default graph = n.drawPDF() view = viewer.View(graph) .. image:: /auto_probabilistic_modeling/distributions/images/sphx_glr_plot_minimum_volume_level_set_1D_001.png :alt: plot minimum volume level set 1D :class: sphx-glr-single-img We want to compute the minimum volume LevelSet which contains `alpha`=90% of the distribution. The `threshold` is the value of the PDF corresponding the `alpha`-probability: the points contained in the LevelSet have a PDF value lower or equal to this threshold. .. code-block:: default alpha = 0.9 levelSet, threshold = n.computeMinimumVolumeLevelSetWithThreshold(alpha) threshold .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 0.1031356403794346 The `LevelSet` has a `contains` method. Obviously, the point 0 is in the LevelSet. .. code-block:: default levelSet.contains([0.]) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none True .. code-block:: default def computeSampleInLevelSet(distribution, levelSet, sampleSize = 1000): ''' Generate a sample from given distribution. Extract the sub-sample which is contained in the levelSet. ''' sample = distribution.getSample(sampleSize) dim = distribution.getDimension() # Get the list of points in the LevelSet. inLevelSet = [] for x in sample: if levelSet.contains(x): inLevelSet.append(x) # Extract the sub-sample of the points in the LevelSet numberOfPointsInLevelSet = len(inLevelSet) inLevelSetSample = ot.Sample(numberOfPointsInLevelSet,dim) for i in range(numberOfPointsInLevelSet): inLevelSetSample[i] = inLevelSet[i] return inLevelSetSample .. code-block:: default def from1Dto2Dsample(oldSample): ''' Create a 2D sample from a 1D sample with zero ordinate (for the graph). ''' size = oldSample.getSize() newSample = ot.Sample(size,2) for i in range(size): newSample[i,0] = oldSample[i,0] return newSample .. code-block:: default def drawLevelSet1D(distribution, levelSet, alpha, threshold, sampleSize = 100): ''' Draw a 1D sample included in a given levelSet. The sample is generated from the distribution. ''' inLevelSample = computeSampleInLevelSet(distribution,levelSet,sampleSize) cloudSample = from1Dto2Dsample(inLevelSample) graph = distribution.drawPDF() mycloud = ot.Cloud(cloudSample) graph.add(mycloud) graph.setTitle("%.2f%% of the distribution, sample size = %d, " % (100*alpha, sampleSize)) return graph .. code-block:: default graph = drawLevelSet1D(n, levelSet, alpha, threshold) view = viewer.View(graph) .. image:: /auto_probabilistic_modeling/distributions/images/sphx_glr_plot_minimum_volume_level_set_1D_002.png :alt: 90.00% of the distribution, sample size = 100, :class: sphx-glr-single-img With a Normal, minimum volume Interval -------------------------------------- .. code-block:: default interval = n.computeMinimumVolumeInterval(alpha) interval .. raw:: html

[-1.64485, 1.64485]



.. code-block:: default def drawPDFAndInterval1D(distribution, interval, alpha): ''' Draw the PDF of the distribution and the lower and upper bounds of an interval. ''' xmin = interval.getLowerBound()[0] xmax = interval.getUpperBound()[0] graph = distribution.drawPDF() yvalue = distribution.computePDF(xmin) curve = ot.Curve([[xmin,0.],[xmin,yvalue],[xmax,yvalue],[xmax,0.]]) curve.setColor("black") graph.add(curve) graph.setTitle("%.2f%% of the distribution, lower bound = %.3f, upper bound = %.3f" % (100*alpha, xmin,xmax)) return graph The `computeMinimumVolumeInterval` returns an `Interval`. .. code-block:: default graph = drawPDFAndInterval1D(n, interval, alpha) view = viewer.View(graph) .. image:: /auto_probabilistic_modeling/distributions/images/sphx_glr_plot_minimum_volume_level_set_1D_003.png :alt: 90.00% of the distribution, lower bound = -1.645, upper bound = 1.645 :class: sphx-glr-single-img With a Mixture, minimum volume LevelSet --------------------------------------- .. code-block:: default m = ot.Mixture([ot.Normal(-5.,1.),ot.Normal(5.,1.)],[0.2,0.8]) .. code-block:: default graph = m.drawPDF() view = viewer.View(graph) .. image:: /auto_probabilistic_modeling/distributions/images/sphx_glr_plot_minimum_volume_level_set_1D_004.png :alt: plot minimum volume level set 1D :class: sphx-glr-single-img .. code-block:: default alpha = 0.9 levelSet, threshold = m.computeMinimumVolumeLevelSetWithThreshold(alpha) threshold .. rst-class:: sphx-glr-script-out Out: .. code-block:: none 0.04667473141178892 The interesting point is that a `LevelSet` may be non-contiguous. In the current mixture example, this is not an interval. .. code-block:: default graph = drawLevelSet1D(m, levelSet, alpha, threshold, 1000) view = viewer.View(graph) .. image:: /auto_probabilistic_modeling/distributions/images/sphx_glr_plot_minimum_volume_level_set_1D_005.png :alt: 90.00% of the distribution, sample size = 1000, :class: sphx-glr-single-img With a Mixture, minimum volume Interval --------------------------------------- .. code-block:: default interval = m.computeMinimumVolumeInterval(alpha) interval .. raw:: html

[-5.44003, 6.72227]



The `computeMinimumVolumeInterval` returns an `Interval`. The bounds of this interval are different from the previous `LevelSet`. .. code-block:: default graph = drawPDFAndInterval1D(m, interval, alpha) view = viewer.View(graph) plt.show() .. image:: /auto_probabilistic_modeling/distributions/images/sphx_glr_plot_minimum_volume_level_set_1D_006.png :alt: 90.00% of the distribution, lower bound = -5.440, upper bound = 6.722 :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 0.502 seconds) .. _sphx_glr_download_auto_probabilistic_modeling_distributions_plot_minimum_volume_level_set_1D.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_minimum_volume_level_set_1D.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_minimum_volume_level_set_1D.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_