scipy.stats.qmc.

LatinHypercube#

class scipy.stats.qmc.LatinHypercube(d, *, scramble=True, strength=1, optimization=None, rng=None)[Quelle]#

Latin-Hypercube-Sampling (LHS).

Eine Latin-Hypercube-Stichprobe [1] erzeugt \(n\) Punkte in \([0,1)^{d}\). Jede univariat marginale Verteilung ist geschichtet, wobei genau ein Punkt in \([j/n, (j+1)/n)\) für \(j=0,1,...,n-1\) platziert wird. Sie sind immer noch anwendbar, wenn \(n << d\).

Parameter:
dint

Dimension des Parameterraums.

scramblebool, optional

Wenn False, werden die Stichproben innerhalb der Zellen eines mehrdimensionalen Gitters zentriert. Andernfalls werden die Stichproben zufällig innerhalb der Zellen des Gitters platziert.

Hinweis

Das Setzen von scramble=False gewährleistet keine deterministische Ausgabe. Verwenden Sie dazu den Parameter rng.

Standard ist True.

Hinzugefügt in Version 1.10.0.

optimization{None, “random-cd”, “lloyd”}, optional

Ob ein Optimierungsschema zur Verbesserung der Qualität nach dem Sampling verwendet werden soll. Beachten Sie, dass dies ein Nachbearbeitungsschritt ist, der nicht garantiert, dass alle Eigenschaften der Stichprobe erhalten bleiben. Standard ist None.

  • random-cd: zufällige Permutationen von Koordinaten zur Verringerung der zentrierten Diskrepanz. Die beste Stichprobe basierend auf der zentrierten Diskrepanz wird ständig aktualisiert. Das zentrierte dispunktbasierte Sampling zeigt eine bessere Raumfüllungsrobustheit gegenüber 2D- und 3D-Unterprojektionen im Vergleich zur Verwendung anderer Dispunktmaße.

  • lloyd: Stört Stichproben mithilfe eines modifizierten Lloyd-Max-Algorithmus. Der Prozess konvergiert zu gleichmäßig verteilten Stichproben.

Hinzugefügt in Version 1.8.0.

Geändert in Version 1.10.0: Fügt lloyd hinzu.

strength{1, 2}, optional

Stärke des LHS. strength=1 erzeugt einen einfachen LHS, während strength=2 einen orthogonalen Array-basierten LHS der Stärke 2 [7], [8] erzeugt. In diesem Fall können nur n=p**2 Punkte mit \(p\) einer Primzahl abgetastet werden. Es schränkt auch \(d <= p + 1\) ein. Standard ist 1.

Hinzugefügt in Version 1.8.0.

rngnumpy.random.Generator, optional

Zustand des Pseudozufallszahlengenerators. Wenn rng None ist, wird ein neuer numpy.random.Generator unter Verwendung von Entropie aus dem Betriebssystem erstellt. Typen außer numpy.random.Generator werden an numpy.random.default_rng übergeben, um einen Generator zu instanziieren.

Geändert in Version 1.15.0: Als Teil des SPEC-007-Übergangs von der Verwendung von numpy.random.RandomState zu numpy.random.Generator wurde dieses Schlüsselwort von seed zu rng geändert. Für eine Übergangszeit funktionieren beide Schlüsselwörter weiterhin, es kann jedoch jeweils nur eines angegeben werden. Nach der Übergangszeit werden Funktionsaufrufe, die das Schlüsselwort seed verwenden, Warnungen ausgeben. Nach einer Deprekationsfrist wird das Schlüsselwort seed entfernt.

Methoden

fast_forward(n)

Führt die Sequenz um n Positionen vorwärts.

integers(l_bounds, *[, u_bounds, n, ...])

Zieht n ganze Zahlen von l_bounds (inklusive) bis u_bounds (exklusive), oder wenn endpoint=True, von l_bounds (inklusive) bis u_bounds (inklusive).

random([n, workers])

Zieht n im halboffenen Intervall [0, 1).

reset()

Setzt die Engine auf den Grundzustand zurück.

Siehe auch

Quasi-Monte Carlo

Hinweise

Wenn LHS zur Integration einer Funktion \(f\) über \(n\) verwendet wird, ist LHS bei fast additiven Integranden sehr effektiv [2]. Mit einem LHS von \(n\) Punkten ist die Varianz des Integrals immer niedriger als bei reinem MC mit \(n-1\) Punkten [3]. Es gibt einen zentralen Grenzwertsatz für LHS bezüglich des Mittels und der Varianz des Integrals [4], aber nicht unbedingt für optimierte LHS aufgrund der Randomisierung.

\(A\) wird als orthogonale Matrix der Stärke \(t\) bezeichnet, wenn in jeder n-Zeilen-t-Spalten-Untermatrix von \(A\): alle \(p^t\) möglichen verschiedenen Zeilen gleich oft vorkommen. Die Elemente von \(A\) sind in der Menge \(\{0, 1, ..., p-1\}\), auch Symbole genannt. Die Einschränkung, dass \(p\) eine Primzahl sein muss, dient der Ermöglichung der modularen Arithmetik. Die Erhöhung der Stärke verleiht den Unterprojektionen einer Stichprobe eine gewisse Symmetrie. Mit Stärke 2 sind die Stichproben entlang der Diagonalen von 2D-Unterprojektionen symmetrisch. Dies kann unerwünscht sein, andererseits verbessert es die Stichprobenverteilung.

Stärke 1 (einfacher LHS) bringt einen Vorteil gegenüber Stärke 0 (MC), und Stärke 2 ist eine nützliche Steigerung gegenüber Stärke 1. Eine Steigerung auf Stärke 3 ist eine kleinere Steigerung, und durch Zufall (scrambled) generierte QMC-Verfahren wie Sobol' und Halton sind leistungsfähiger [7].

Um einen LHS der Stärke 2 zu erstellen, wird die orthogonale Matrix \(A\) randomisiert, indem eine zufällige, bijektive Abbildung der Symbolmenge auf sich selbst angewendet wird. Zum Beispiel könnten in Spalte 0 alle 0en zu 2 werden; in Spalte 1 könnten alle 0en zu 1 werden, usw. Dann wird für jede Spalte \(i\) und jedes Symbol \(j\) ein einfacher, eindimensionaler LHS der Größe \(p\) zum Unterarray hinzugefügt, wo \(A^i = j\). Die resultierende Matrix wird schließlich durch \(p\) geteilt.

Referenzen

[1]

Mckay et al., “A Comparison of Three Methods for Selecting Values of Input Variables in the Analysis of Output from a Computer Code.” Technometrics, 1979.

[2]

M. Stein, “Large sample properties of simulations using Latin hypercube sampling.” Technometrics 29, no. 2: 143-151, 1987.

[3]

A. B. Owen, “Monte Carlo variance of scrambled net quadrature.” SIAM Journal on Numerical Analysis 34, no. 5: 1884-1910, 1997

[4]

Loh, W.-L. “On Latin hypercube sampling.” The annals of statistics 24, no. 5: 2058-2080, 1996.

[5]

Fang et al. “Design and modeling for computer experiments”. Computer Science and Data Analysis Series, 2006.

[6]

Damblin et al., “Numerical studies of space filling designs: optimization of Latin Hypercube Samples and subprojection properties.” Journal of Simulation, 2013.

[7] (1,2)

A. B. Owen , “Orthogonal arrays for computer experiments, integration and visualization.” Statistica Sinica, 1992.

[8]

B. Tang, “Orthogonal Array-Based Latin Hypercubes.” Journal of the American Statistical Association, 1993.

[9]

Seaholm, Susan K. et al. (1988). Latin hypercube sampling and the sensitivity analysis of a Monte Carlo epidemic model. Int J Biomed Comput, 23(1-2), 97-112. DOI:10.1016/0020-7101(88)90067-0

Beispiele

Generiert Stichproben aus einem Latin-Hypercube-Generator.

>>> from scipy.stats import qmc
>>> sampler = qmc.LatinHypercube(d=2)
>>> sample = sampler.random(n=5)
>>> sample
array([[0.1545328 , 0.53664833], # random
        [0.84052691, 0.06474907],
        [0.52177809, 0.93343721],
        [0.68033825, 0.36265316],
        [0.26544879, 0.61163943]])

Berechnet die Qualität der Stichprobe anhand des Diskrepanzkriteriums.

>>> qmc.discrepancy(sample)
0.0196... # random

Stichproben können auf Grenzen skaliert werden.

>>> l_bounds = [0, 2]
>>> u_bounds = [10, 5]
>>> qmc.scale(sample, l_bounds, u_bounds)
array([[1.54532796, 3.609945 ], # random
        [8.40526909, 2.1942472 ],
        [5.2177809 , 4.80031164],
        [6.80338249, 3.08795949],
        [2.65448791, 3.83491828]])

Nachfolgend finden Sie weitere Beispiele, die alternative Möglichkeiten zur Konstruktion von LHS mit noch besserer Raumabdeckung zeigen.

Verwendung eines Basen-LHS als Referenz.

>>> sampler = qmc.LatinHypercube(d=2)
>>> sample = sampler.random(n=5)
>>> qmc.discrepancy(sample)
0.0196...  # random

Verwenden Sie das Schlüsselwortargument optimization, um einen LHS mit geringerer Diskrepanz bei höherem Rechenaufwand zu erzeugen.

>>> sampler = qmc.LatinHypercube(d=2, optimization="random-cd")
>>> sample = sampler.random(n=5)
>>> qmc.discrepancy(sample)
0.0176...  # random

Verwenden Sie das Schlüsselwortargument strength, um einen orthogonalen Array-basierten LHS der Stärke 2 zu erzeugen. In diesem Fall muss die Anzahl der Stichprobenpunkte das Quadrat einer Primzahl sein.

>>> sampler = qmc.LatinHypercube(d=2, strength=2)
>>> sample = sampler.random(n=9)
>>> qmc.discrepancy(sample)
0.00526...  # random

Optionen könnten kombiniert werden, um einen optimierten zentrierten, orthogonalen Array-basierten LHS zu erzeugen. Nach der Optimierung wäre nicht garantiert, dass das Ergebnis die Stärke 2 aufweist.

Realweltbeispiel

In [9] wurde eine Latin-Hypercube-Sampling-Strategie (LHS) verwendet, um einen Parameterraum zu untersuchen und die Bedeutung jedes Parameters eines epidemiologischen Modells zu analysieren. Eine solche Analyse wird auch als Sensitivitätsanalyse bezeichnet.

Da die Dimensionalität des Problems hoch ist (6), ist die Abdeckung des Raumes rechenintensiv. Wenn numerische Experimente kostspielig sind, ermöglicht QMC Analysen, die bei Verwendung eines Gitters nicht möglich wären.

Die sechs Parameter des Modells repräsentierten die Krankheitswahrscheinlichkeit, die Rückzugswahrscheinlichkeit und vier Kontaktwahrscheinlichkeiten. Die Autoren nahmen für alle Parameter gleichmäßige Verteilungen an und generierten 50 Stichproben.

Um das Protokoll mit scipy.stats.qmc.LatinHypercube zu replizieren, besteht der erste Schritt darin, eine Stichprobe im Einheits-Hyperwürfel zu erstellen

>>> from scipy.stats import qmc
>>> sampler = qmc.LatinHypercube(d=6)
>>> sample = sampler.random(n=50)

Dann kann die Stichprobe auf die entsprechenden Grenzen skaliert werden

>>> l_bounds = [0.000125, 0.01, 0.0025, 0.05, 0.47, 0.7]
>>> u_bounds = [0.000375, 0.03, 0.0075, 0.15, 0.87, 0.9]
>>> sample_scaled = qmc.scale(sample, l_bounds, u_bounds)

Eine solche Stichprobe wurde verwendet, um das Modell 50 Mal auszuführen, und eine polynomiale Response-Oberfläche wurde konstruiert. Dies ermöglichte es den Autoren, die relative Bedeutung jedes Parameters über den Bereich der Möglichkeiten jedes anderen Parameters zu untersuchen.

In diesem Computerexperiment zeigten sie eine 14-fache Reduzierung der benötigten Stichproben, um einen Fehler von unter 2 % auf ihrer Response-Oberfläche aufrechtzuerhalten, verglichen mit einer Gitterabtastung.

Auf dieser Seite