Importance Sampling

Importance Sampling (im Deutschen manchmal auch Stichprobenentnahme nach Wichtigkeit, oder Stichprobenziehung nach Wichtigkeit^[1] genannt) ist ein Begriff aus der Statistik, der die Technik zur Erzeugung von Stichproben anhand einer Wahrscheinlichkeitsverteilung beschreibt. Importance Sampling ist eine von mehreren Möglichkeiten zur Varianzreduktion, also zur Steigerung der Effizienz von Monte-Carlo-Simulationen.

Monte-Carlo-Simulationen werden oft benutzt, um den Erwartungswert

${\displaystyle \operatorname {E} _{X}[f(X)]={\begin{cases}\displaystyle \sum _{x\in \Omega }f(x)p(x)&{\text{(diskreter Fall)}}\\\displaystyle \int _{\Omega }f(x)p(x)\mathrm {d} x&{\text{(stetiger Fall)}}\end{cases}}}$

einer reellen Zufallsvariablen ${\displaystyle f(X)}$ zu berechnen, wobei die Wahrscheinlichkeitsverteilung der Zufallsvariablen ${\displaystyle X}$ und die Funktion ${\displaystyle f}$ bekannt sind. Die Zufallsvariable ${\displaystyle X}$ nimmt Werte in der Ergebnismenge ${\displaystyle \Omega }$ an. Für eine Realisierung ${\displaystyle x\in \Omega }$ der Zufallsvariablen ${\displaystyle X}$ ist ${\displaystyle f(x)\in \mathbb {R} }$ eine Realisierung der Zufallsvariablen ${\displaystyle f(X)}$.

Im diskreten Fall ist ${\displaystyle \Omega \subset \mathbb {R} ^{d}}$ eine abzählbare Menge und die Wahrscheinlichkeitsverteilung von ${\displaystyle X}$ ist durch die Wahrscheinlichkeitsfunktion ${\displaystyle p}$ mit ${\displaystyle \sum _{x\in \Omega }p(x)=1}$ gegeben. Im stetigen Fall ist ${\displaystyle \Omega \subseteq \mathbb {R} ^{d}}$ typischerweise ein ${\displaystyle d}$-dimensionales Intervall und die Wahrscheinlichkeitsverteilung von ${\displaystyle X}$ ist durch eine Wahrscheinlichkeitsdichte ${\displaystyle p}$ mit der Eigenschaft ${\displaystyle \int _{\Omega }p(x)\mathrm {d} x=1}$ gegeben. Für eine Realisierung ${\displaystyle x\in \Omega }$ der Zufallsvariablen ${\displaystyle X}$ ist ${\displaystyle f(x)\in \mathbb {R} }$ eine Realisierung der Zufallsvariablen ${\displaystyle f(X)}$. Da die Ergebnismenge ${\displaystyle \Omega }$ im Allgemeinen hochdimensional sein kann, kann die Berechnung des Erwartungswertes sehr schwierig oder zeitaufwendig sein, so dass dann eine Approximation durch Monte-Carlo-Simulation sinnvoll ist.

Bei Anwendungen im Bereich der Physik ist z. B. die Ergebnismenge ${\displaystyle \Omega }$ der Phasenraum der Teilchen im System, die Zufallsvariablen ${\displaystyle X}$ und ${\displaystyle f(X)}$ sind interessierende Größen und ${\displaystyle p(x)}$ ist proportional zum Boltzmann-Faktor.

Statt den Erwartungswert analytisch zu berechnen oder numerisch zu approximieren, berechnet man im einfachsten Fall den Monte-Carlo-Schätzwert

${\displaystyle \operatorname {\hat {E}} _{X}[f(X)]={\frac {1}{n}}\sum _{i=1}^{n}f(x_{i}),}$

für den Erwartungswert ${\displaystyle \operatorname {E} _{X}[f(X)]}$, wobei die Zufallszahlen ${\displaystyle x_{1},\dots ,x_{n}}$ Realisierungen von stochastisch unabhängigen und identisch verteilten Zufallsvariablen ${\displaystyle X_{1},\dots ,X_{n}}$ mit der Verteilung von ${\displaystyle X\sim p}$ sind. In statistischer Terminologie ist ${\displaystyle (X_{1},\dots ,X_{n})}$ eine einfache Zufallsstichprobe mit Stichprobenumfang ${\displaystyle n}$. Der Schätzwert ${\displaystyle \operatorname {\hat {E}} _{X}[f(X)]}$ ist eine Realisierung des zugehörigen Schätzers

${\displaystyle \operatorname {\tilde {E}} _{X}[f(X)]={\frac {1}{n}}\sum _{i=1}^{n}f(X_{i})}$,

der eine Zufallsvariable ist, die nach dem starken Gesetz der großen Zahlen mit Wahrscheinlichkeit Eins gegen den zu schätzenden Erwartungswert konvergiert. Da die Berechnung durch Monte-Carlo-Simulation eine statistische Schätzung ist, gibt es einen Schätzfehler der z. B. durch die Varianz des Schätzer beschrieben werden kann. Beim Importance Sampling wird die zuvor beschriebene einfache Monte-Carlo-Schätzmethode mit dem Ziel modifiziert, die Varianz des Schätzer zu reduzieren bzw. die Genauigkeit bei gegebenem Stichprobenumfang zu erhöhen.

Der Standardansatz zur Approximation des Erwartungswertes ${\displaystyle \operatorname {E} [X]}$ einer Zufallsvariablen ${\displaystyle X}$ durch Monte-Carlo-Simulation besteht darin, Zufallszahlen ${\displaystyle x_{1},\dots ,x_{n}}$ als Realisierungen stochastisch unabhängiger und identisch verteilter Zufallsvariablen ${\displaystyle X_{1},\dots ,X_{n}}$ mit der Wahrscheinlichkeitsverteilung von ${\displaystyle X}$ zu erzeugen und dann den gesuchten Erwartungswert ${\displaystyle \operatorname {E} [X]}$ durch den Monte-Carlo-Schätzwert

${\displaystyle \operatorname {\hat {E}} _{\mathrm {MC} }[X]={\frac {1}{n}}\sum _{i=1}^{n}x_{i}}$

zu approximieren.

Beim Importance-Sampling werden stattdessen Zufallszahlen aus einer modifizierten Wahrscheinlichkeitsverteilung mit dem Ziel verwendet, durch eine Varianzreduktion zu rechentechnisch effizienteren Berechnung zu kommen.

Für eine diskrete Zufallsvariable ${\displaystyle X}$, die Werte ${\displaystyle x\in \Omega }$ mit den Wahrscheinlichkeiten ${\displaystyle p(x)=P(X=x)}$ annimmt, ist die Grundlage des Importance-Sampling eine Umdeutung des Erwartungswertes

${\displaystyle \operatorname {E} [X]=\sum _{x\in \Omega }xp(x)=\sum _{x\in \Omega }{\frac {xp(x)}{w(x)}}w(x)=\operatorname {E} \left[{\frac {Yp(Y)}{w(Y)}}\right]\;.}$

Dabei ist ${\displaystyle Y}$ eine diskrete Zufallsvariable, die Werte ${\displaystyle y\in \Omega }$ mit den positiven Wahrscheinlichkeiten ${\displaystyle w(y)=P(Y=y)}$ annimmt. Die Erzeugung von ${\displaystyle n}$ Zufallen ${\displaystyle y_{1},\dots ,y_{n}}$ aus der Verteilung von ${\displaystyle Y}$ führt dann zum Monte-Carlo-Schätzwert

${\displaystyle {\frac {1}{n}}\sum _{i=1}^{n}{\frac {y_{i}p(y_{i})}{w(y_{i})}}\approx \operatorname {E} \left[{\frac {Yp(Y)}{w(Y)}}\right]\;,}$

der als Importance-Sampling-Schätzwert

${\displaystyle \operatorname {\hat {E}} _{\mathrm {IS} }[X]={\frac {1}{n}}\sum _{i=1}^{n}{\frac {y_{i}p(y_{i})}{w(y_{i})}}}$

für ${\displaystyle \operatorname {E} [X]}$ interpretiert wird. Der Umweg über die Erzeugung aus Zufallszahlen mit einer anderen Verteilung kann lohnend sein, da durch eine geeignete Wahl der Verteilung von ${\displaystyle Y}$ ein Importance-Sampling-Schätzer eine erhebliche kleinere Varianz als der gewöhnliche Monte-Carlo-Schätzer haben kann.

Die Grundidee des Importance-Sampling sei an einem einfachen Beispiel mit ${\displaystyle \Omega =\{-10,0,10\}}$ veranschaulicht. Die Zufallsvariable ${\displaystyle X}$ mit der Wahrscheinlichkeitsfunktion

${\displaystyle p(x)=P(X=x)={\begin{cases}{\frac {1}{10}}&{\text{für }}x\in \{-10,10\}\\{\frac {8}{10}}&{\text{für }}x=0\end{cases}}}$

hat den Erwartungswert ${\displaystyle \operatorname {E} [X]=0}$ und die Varianz ${\displaystyle \mathrm {Var} [X]=20.}$ Die Wahrscheinlichkeitsfunktion ${\displaystyle p}$ sei als bekannt vorausgesetzt. Der Erwartungswert sei als unbekannt angenommen und soll durch Simulation bestimmt werden. Der gewöhnliche Monte-Carlo-Schätzer

${\displaystyle \operatorname {\tilde {E}} _{\mathrm {MC} }[X]={\frac {1}{n}}\sum _{i=1}^{n}X_{i}}$

basierend auf ${\displaystyle n}$ Zufallszahlen aus der Verteilung von ${\displaystyle X}$ hat dann den Erwartungswert Null und die Varianz

${\displaystyle \mathrm {Var} [{\tilde {E}}_{\mathrm {MC} }[X]]={\frac {\mathrm {Var} [X]}{n}}={\frac {20}{n}}\;.}$

Für das Importance-Sampling wird der Erwartungswert ${\displaystyle \operatorname {E} [X]}$ als

${\displaystyle \operatorname {E} [X]=\sum _{x\in \Omega }xp(x)=\sum _{x\in \Omega }{\frac {xp(x)}{w(x)}}{w(x)}=\operatorname {E} \left[{\frac {Yp(Y)}{w(Y)}}\right]}$

geschrieben und die positiven ${\displaystyle w(x)}$ werden so gewählt, dass sie sich zu Eins addieren und damit als Wahrscheinlichkeiten ${\displaystyle P(Y=y)=w(y)}$ für eine Zufallsvariable ${\displaystyle Y}$ interpretiert werden können. Der Importance-Sampling-Schätzer

${\displaystyle \operatorname {\tilde {E}} _{\mathrm {IS} }[X]={\frac {1}{n}}\sum _{i=1}^{n}{\frac {Y_{i}p(Y_{i})}{w(Y_{i})}},}$

der als Funktion der Zufallsvariablen ${\displaystyle Y_{1},\dots ,Y_{n}}$ selbst eine Zufallsvariable ist, hat in diesem Beispiel den Erwartungswert Null und die Varianz

${\displaystyle \mathrm {Var} [\operatorname {\hat {E}} _{\mathrm {IS} }[X]]={\frac {1}{n}}\left({\frac {1}{w(-10)}}+{\frac {1}{w(10)}}\right)}$

Die Varianz des Schätzers ${\displaystyle \operatorname {\tilde {E}} _{\mathrm {IS} }[X]}$ hängt also von den gewählten Wahrscheinlichkeiten ${\displaystyle w(y)}$ ab und kann bei geeigneter Wahl erheblich kleiner als die Varianz des gewöhnlichen Monte-Carlo-Schätzers sein. Z. B. ergibt sich für die Wahl ${\displaystyle w(y)=1/3}$ für ${\displaystyle y\in \Omega }$ die Varianz

${\displaystyle \mathrm {Var} [\operatorname {\tilde {E}} _{\mathrm {IS} }[X]]={\frac {6}{n}}<\mathrm {Var} [\operatorname {\tilde {E}} _{\mathrm {MC} }[X]]={\frac {20}{n}}\;.}$

Durch den Importance-Sampling-Schätzer mit der Verteilung von ${\displaystyle Y}$ lässt sich also der Erwartungswert ${\displaystyle \mathbb {E} [X]}$ mit kleinerer Varianz als durch den gewöhnlichen Monte-Carlo-Schätzer bestimmen. Damit erhält man bei gleicher Anzahl ${\displaystyle n}$ von erzeugten Zufallszahlen eine größerer Genauigkeit.

Im einfachsten Fall (einfache Stichprobenentnahme, englisch simple sampling) werden aus dem Ergebnisraum gleichverteilt zufällig Zustände ${\displaystyle Y_{i}\sim {\text{uniform}}(\Omega )}$ für die Stichprobe ${\displaystyle S}$ ausgewählt. Dann ergibt sich für den geschätzten Mittelwert:

${\displaystyle \operatorname {\hat {E}} _{X}[f(X)]={\frac {\sum _{y\in S}p(y)\,f(y)}{\sum _{y\in S}p(y)}},}$

wobei die Summation über die zufälligen Realisierungen ${\displaystyle y_{i}\sim Y_{i}}$in der Stichprobe läuft. ${\displaystyle p(y)}$ ist die ursprüngliche Wahrscheinlichkeit(sdichte) für die – durch Simple Sampling erzeugte – Realisierung ${\displaystyle y}$. ${\displaystyle f(y)}$ ist eine Funktionsauswertung

Die Methode des Simple Sampling ist meistens nicht sehr effizient, da oft nur wenige relevante Zustände in die Mittelwertbildung eingehen. Um dieses Problem zu umgehen und so die Standardabweichung des gemessenen Mittelwertes bei gleichem Stichprobenumfang zu reduzieren, versucht man Zustände mit einem größeren Gewicht häufiger in die Mittelwertbildung eingehen zu lassen als Zustände mit einem geringeren Gewicht: Der obigen Schätzer des Simple Sampling kann durch Erweitern mit ${\displaystyle 1=W(y)/W(y)}$ auch wie folgt ausgedrückt werden:

${\displaystyle \operatorname {E} _{X}[f(X)]={\frac {\sum _{y\in S}p(y)/W(y)\,f(y)W(y)}{\sum _{y\in S}p(y)/W(y)W(y)}}.}$

Werden Realisierungen ${\displaystyle y}$ mit der Wahrscheinlichkeit ${\displaystyle W(y)}$ erzeugt (Stichprobenentnahme nach Wichtigkeit englisch importance sampling, ${\displaystyle y\sim W}$), wird also eine andere Stichprobe S' erzeugt, so berechnet sich der geschätzte Mittelwert in der Folge einfach mithilfe von

${\displaystyle \operatorname {E} _{X}[f(X)]={\frac {\sum _{y\in S'}p(y)/W(y)f(y)}{\sum _{y\in S'}p(y)/W(y)}},}$ für ${\displaystyle y\sim W}$

Die Wahrscheinlichkeitsdichte ${\displaystyle W}$ wird auch als "biased distribution", "proposal distribution" oder "sample distribution" bezeichnet. Da die Stichprobe aus der Verteilung ${\displaystyle W}$ gezogen wurde, müssen die Erwartungswerte durch entsprechendes Reweighting mit ${\displaystyle W}$ errechnet werden (siehe Formel) und nicht einfach als arithmetischen Mittel. Dieses Reweighting wird beispielsweise beim "Temperature Reweighting" von Monte-Carlo Simulationen molekularer Systeme genutzt bei denen Aussagen über angrenzende Temperaturen gemacht werden sollen^[2].

Um eine Stichprobenentnahme nach Wichtigkeit in der Praxis zu erreichen, geht man von einer Startkonfiguration aus und erzeugt mithilfe des Metropolisalgorithmus eine Markow-Kette aus Systemzuständen.

Werden die Realisierungen ${\displaystyle y}$ mit einer Wahrscheinlichkeit ${\displaystyle W(y)}$ proportional zu ${\displaystyle p(y)}$ vorgeschlagen (das ist gerade die Metropoliswahl), so ergibt sich (wie beim Gesetz der großen Zahlen)

${\displaystyle \operatorname {\hat {E}} _{X}[f(X)]={\frac {1}{N}}\sum _{y\in S'}f(y).}$

Gerade, dass hier nur die Proportionalität ${\displaystyle W(y)\propto p(y)}$ erforderlich ist, ist ein Vorteil der Methode, da die Zustandssumme nicht ausgewertet werden muss.

Simple sampling erhält man für der Fall, dass die Vorschlagsdichte konstant ist: ${\displaystyle W(y)={\text{const}}}$.

Das multikanonische Ensemble kann mit dem Wang-Landau-Algorithmus simuliert werden, indem die Wahl ${\displaystyle W(x)=1/D(\operatorname {E} (x))}$ getroffen wird (wobei ${\displaystyle D(\operatorname {E} (x))}$ die Zustandsdichte der Energie ist, welche dem Zustand ${\displaystyle x}$ zugeordnet ist).

• W. K. Hastings: Monte Carlo Sampling Methods Using Markov Chains and Their Applications. In: Biometrika. Band 57, 1970, S. 97–109. 
• Thomas Müller-Gronbach, Erich Novak, Klaus Ritter: Monte Carlo-Algorithmen. Springer-Verlag, Berlin 2012, ISBN 978-3-540-89140-6, Abschnitt 5.4 Importance Sampling, S. 155–166, doi:10.1007/978-3-540-89141-3. 
• Surya T. Tokdar, Robert E. Kass: Importance sampling: a review. In: WIREs Computational Statistics. Band 2, Nr. 1, 2010, S. 54–60, doi:10.1002/wics.56. 
• Christian P. Robert, George Casella: Monte Carlo Statistical Methods (= Springer Texts in Statistics). 2. Auflage. Springer, 2004, ISBN 0-387-21239-6, Kap. 3.3 Importance Sampling, S. 90–107, doi:10.1007/978-1-4757-4145-2. 
• R. Srinivasan: Importance sampling – Applications in communications and detection. Springer-Verlag, Berlin 2002, ISBN 978-3-540-43420-7. 
• Suojin Wang: Importance Sampling. In: Samuel Kotz et al. (Hrsg.): Encyclopedia of Statistical Sciences. 2. Auflage. Band 5. Wiley, New York 2006, ISBN 978-0-471-15044-2, S. 3347–3353, doi:10.1002/0471667196. 

1. ↑ International Statistical Institute: Glossary of statistical terms.
2. ↑ Bachmann, M. (2014). Thermodynamics and Statistical Mechanics of Macromolecular Systems. Vereinigtes Königreich: Cambridge University Press. Seiten 104, 105 Google books