Akzeptanz- und Verwerfungsmethode

• Wir diskutieren nun eine weitere Methode, mit der man Pseudozufallszahlen $y_1,y_2,\ldots$ erzeugen kann,
  • die als Realisierungen von unabhängigen und identisch verteilten Zufallsvariablen $Y_1,Y_2\ldots$ mit einer vorgegebenen Verteilungsfunktion $G$ aufgefasst werden können,
  • wenn dabei von (nicht notwendig $(0,1]$-gleichverteilten) Pseudozufallszahlen $x_1,x_2,\ldots$ ausgegangen wird,
  • die Realisierungen von unabhängigen und identisch verteilten Zufallsvariablen $X_1,X_2\ldots$ mit einer gewissen Verteilungsfunktion $F$ sind,
  • so dass $G$ absolutstetig bezüglich $F$ ist, wobei die Dichte $g(x)=dG(x)/dF(x)$ beschränkt ist,
  • d.h., es gibt eine Konstante $c>0$ mit

    $$g(x)\le c \qquad G(y)=\int\limits_{-\infty}^y g(x)\,dF(x)\,,\qquad \forall\, x,y\in\mathbb{R}\,.$$ (19)

    
    

• Wir betrachten zunächst den diskreten Fall.
  • Sei $a_j=j$ für jedes $j=0,1,\ldots$, und seien ${\mathbf{p}}=(p_0,p_1,\ldots)^\top$ bzw. ${\mathbf{q}}=(q_0,q_1,\ldots)^\top$ zwei beliebige Wahrscheinlichkeitsfunktionen, so dass für jedes $j=0,1,\ldots$ aus $p_j=0$ die Gültigkeit von $q_j=0$ folgt.
  • Sei $X:\Omega\to\{0,1,\ldots\}$ eine Zufallsvariable mit $P(X=j)=p_j$ für jedes $j=0,1,\ldots$,
  • und sei $c>0$ eine positive Zahl, so dass

    $$\Bigl(\;g(j)=\;\Bigr)\qquad\frac{q_j}{p_j}\le c \mbox{für jedes $j\ge 0$\ mit $p_j>0$.}$$ (20)

    
    

Theorem 3.5    

• Sei $(U_1,X_1),(U_2,X_2),\ldots$ eine Folge von unabhängigen und identisch verteilten Zufallsvektoren, deren Komponenten unabhängig sind, wobei $U_i$ eine $(0,1]$-gleichverteilte Zufallsvariable und $X_i$ gemäß ${\mathbf{p}}$ verteilt ist.
• Dann ist
  • die Zufallsvariable

    $$I=\min\Bigl\{k\ge 1: U_k<\frac{q_{X_k}}{c\,p_{X_k}}\Bigr\}$$ (21)

    
    
    geometrisch verteilt mit Erwartungswert $c$, d.h. $I\sim\,{\rm Geo}(c^{-1})$,
  • und die Zufallsvariable $Y=X_I$ ist gemäß ${\mathbf{q}}$ verteilt.

Beweis

 

• Aus der Definitionsgleichung (21) von $I$ ergibt sich, dass für jedes $j\ge 1$
  

  $$P(I=j) = P\Bigl(U_1\ge\frac{q_{X_1}}{c\,p_{X_1}}\,,\ldots,\, U_{j-1}\ge\frac{q_{X_{j-1}}}{c\,p_{X_{j-1}}}\,,\,U_j<\frac{q_{X_j}}{c\,p_{X_j}}\Bigr)$$

  $$= P\Bigl(U_1\ge\frac{q_{X_1}}{c\,p_{X_1}}\Bigr)\,\ldots\, P\Bigl(U_... ..._{X_{j-1}}}{c\,p_{X_{j-1}}}\Bigr)\, P\Bigl(U_j<\frac{q_{X_j}}{c\,p_{X_j}}\Bigr)$$

  $$= p\,q^{j-1}\,,$$

  
  
  • wobei $q=1-p$ und
    

    $$p = P\Bigl(U_1<\frac{q_{X_1}}{c\,p_{X_1}}\Bigr)$$

    $$= \sum\limits_{k:\, p_k>0} P\Bigl(U_1<\frac{q_{X_1}}{c\,p_{X_1}}\mid X_1=k\Bigr)\,P(X_1=k)$$

    $$= \sum\limits_{k:\, p_k>0} P\Bigl(U_1<\frac{q_k}{c\,p_k}\Bigr)\, p_k$$

    $$= \sum\limits_{k:\, p_k>0}\frac{q_k}{c\,p_k}\; p_k = \;\frac{1}{c}\;.$$

    
    

  • Hieraus folgt, dass $I\sim\,{\rm Geo}(c^{-1})$.
• Außerdem gilt für jedes $j\ge 1$ mit $p_j>0$
  

  $$P(Y=j) = P(X_I=j) \;=\; \sum\limits_{k=1}^\infty P(X_I=j,\, I=k)$$

  $$= \sum\limits_{k=1}^\infty P(X_k=j,\,I=k) \;=\; \sum\limits_{k=1}^\infty P(X_k=j)\, q^{k-1}\,P\Bigl(U_k\le\;\frac{q_j}{cp_j}\Bigr)$$

  $$= \sum\limits_{k=1}^\infty p_j\,q^{k-1}\;\frac{q_j}{cp_j} \;=\;\frac{q_j}{c}\; \sum\limits_{k=1}^\infty q^{k-1}$$

  $$= \frac{q_j}{c}\;\frac{1}{1-q}\; \;=\; q_j\,,$$

  
  
  und für jedes $j\ge 1$ mit $p_j=0$ gilt
  $$P(Y=j) \;=\; \sum\limits_{k=1}^\infty P(X_k=j,\, I=k)\;\le\; \sum\limits_{k=1}^\infty P(X_k=j)\;=\; 0\,.$$
  
   
    $\Box$
  
  
  

Beachte

 

• Aus Theorem 3.5 ergibt sich, dass im Mittel $c$ Pseudozufallszahlen erforderlich sind, die $F$-verteilt sind, um eine $G$-verteilte Pseudozufallszahl zu erzeugen.
• Wenn mehrere Alternativen für die Wahl der Verteilungsfunktion $F$ zur Verfügung stehen,
  • die ähnlich günstige Eigenschaften hinsichtlich der Erzeugung von $F$-verteilten Pseudozufallszahlen aufweisen,
  • dann sollte somit diejenige Verteilungsfunktion $F$ gewählt werden, für die $c$ am kleinsten ist.
• Außerdem ergibt sich aus Theorem 3.5,
  • dass die Werte $g(x)$ bzw. $g(j)$ der Dichte in (19) bzw. (20) nicht vollständig,
  • sondern nur bis auf einen konstanten Proportionalitätsfaktor bekannt sein müssen.

Im allgemeinen (d.h. nicht notwendig diskreten) Fall kann man ähnlich vorgehen, wobei das folgende Resultat als Grundlage zur Konstruktion von Akzeptanz- und Verwerfungsalgorithmen dient.

Theorem 3.6    

• Seien $F,G:\mathbb{R}\to[0,1]$ zwei beliebige Verteilungsfunktionen, so dass $(\ref{def.dic.geh})$ gilt.
• Sei $(U_1,X_1),(U_2,X_2),\ldots$ eine Folge von unabhängigen und identisch verteilten Zufallsvektoren, deren Komponenten unabhängig sind, wobei $U_i$ eine $(0,1]$-gleichverteilte Zufallsvariable und $X_i$ gemäß $F$ verteilt ist.
• Dann ist die Zufallsvariable

  $$I=\min\Bigl\{k\ge 1: U_k<\frac{g(X_k)}{c}\Bigr\}$$ (22)

  
  
  geometrisch verteilt mit Erwartungswert $c$, d.h. $I\sim\,{\rm Geo}(c^{-1})$, und die Zufallsvariable $Y=X_I$ ist gemäß $G$ verteilt.

Beweis

 

• Ähnlich wie im Beweis von Theorem 3.5 ergibt sich, dass $P(I=j)=p\,q^{j-1}$ für jedes $j\ge 1$ gilt, wobei
  

  $$p = P\Bigl(U_1<\frac{g(X_1)}{c}\Bigr)= \int\limits_\mathbb{R}P\Bigl(U_1<\frac{g(X_1)}{c}\mid X_1=x\Bigr)\,dF(x)$$

  $$= \int\limits_\mathbb{R}P\Bigl(U_1<\frac{g(x)}{c}\Bigr)\,dF(x)= \int\limits_\mathbb{R}\;\frac{g(x)}{c}\;dF(x)=\frac{1}{c}\;.$$

  
  

• Außerdem gilt für jedes $y\in\mathbb{R}$
  

  $$P(Y\le y) = P(X_I\le y)\;=\; \sum\limits_{k=1}^\infty P(X_I\le y,\, I=k)$$

  $$= \sum\limits_{k=1}^\infty P(X_k\le y,\, I=k)$$

  $$= \sum\limits_{k=1}^\infty \int\limits_{-\infty}^y P(I=k\mid X_k=v)\, dF(v)$$

  $$= \sum\limits_{k=1}^\infty q^{k-1}\int\limits_{-\infty}^y P\Bigl(U_k<\;\frac{g(v)}{c}\Bigr)\, dF(v)$$

  $$= \frac{1}{1-q}\;\int\limits_{-\infty}^y \;\frac{g(v)}{c}\; dF(v)$$

  $$= \int\limits_{-\infty}^y g(v) dF(v)\;=\; G(y)\,,$$

  
  
  wobei $1-q=p=c^{-1}$.
  
  $\Box$

Auf die gleiche Weise ergibt sich die folgende vektorielle Version von Theorem 3.6.

Theorem 3.7    

• Sei $m\ge 1$ eine beliebige, jedoch fest vorgegebene natürliche Zahl, und seien $F,G:\mathbb{R}^m\to[0,1]$ zwei beliebige Verteilungsfunktionen (von $m$-dimensionalen Zufallsvektoren), so dass es gibt eine Konstante $c>0$ mit

  $$g({\mathbf{x}})\le c \qquad G({\mathbf{y}})=\int\limits_{(-\infty,{\mathbf{y}}]} g({\m... ...dF({\mathbf{x}})\,,\qquad \forall\, {\mathbf{x}},{\mathbf{y}}\in\mathbb{R}^m\,.$$ (23)

  
  

• Sei $(U_1,{\mathbf{X}}_1),(U_2,{\mathbf{X}}_2),\ldots$ eine Folge von unabhängigen und identisch verteilten Zufallsvektoren, deren Komponenten unabhängig sind, wobei $U_i$ eine $(0,1]$-gleichverteilte Zufallsvariable und ${\mathbf{X}}_i$ gemäß $F$ verteilt ist.
• Dann ist die Zufallsvariable

  $$I=\min\Bigl\{k\ge 1: U_k<\frac{g({\mathbf{X}}_k)}{c}\Bigr\}$$ (24)

  
  
  geometrisch verteilt mit Erwartungswert $c$, d.h. $I\sim\,{\rm Geo}(c^{-1})$, und der Zufallsvektor ${\mathbf{Y}}={\mathbf{X}}_I$ ist gemäß $G$ verteilt.

Beispiele

 

1. Gleichverteilung in beschränkten Borel-Mengen 
   • Der Zufallsvektor ${\mathbf{X}}:\Omega\to\mathbb{R}^m$ (mit der Verteilungsfunktion $F$) sei gleichverteilt in dem Quader $(-1,1]^m$, und $B\in\mathcal{B}((-1,1]^m$ sei eine beliebige Borel-Menge in $(-1,1]^m$ mit positivem Lebesgue-Maß $\vert B\vert$.
   • Dann ist die Verteilungsfunktion $G:\mathbb{R}^m\to[0,1]$ mit
     $$G({\mathbf{y}})=\int\limits_{(-\infty,{\mathbf{y}}]}\;\frac{{1\hs... ...}{\vert B\vert}\;dF({\mathbf{x}})\,,\qquad\forall\,{\mathbf{y}}\in\mathbb{R}^m$$
     absolutstetig bezüglich $F$, und für die (Radon-Nikodym-) Dichte $g:\mathbb{R}^m\to[0,\infty)$ gilt
     $$g({\mathbf{x}})=\frac{{1\hspace{-1mm}{\rm I}}({\mathbf{x}}\in B)}{\vert B\vert}\le c=\vert B\vert^{-1}$$   bzw.$$\qquad \frac{g({\mathbf{x}})}{c}\;={1\hspace{-1mm}{\rm I}}({\mathbf{x}}\in B)\,, \qquad\forall\,{\mathbf{x}}\in\mathbb{R}^m\,.$$

   • Gemäß Theorem 3.7 können nun in $B$ gleichverteilte Pseudozufallsvektoren ${\mathbf{y}}_1,{\mathbf{y}}_2,\ldots$ auf die folgende Weise erzeugt werden.

     1.

     Generiere $m$ Pseudozufallszahlen $u_1,\ldots,u_m$, die im Intervall $(0,1]$ gleichverteilt sind.

     2.

     Falls $(2u_1-1,\ldots,2u_m-1)^\top\not\in B$ gilt, dann kehre zu Schritt 1 zurück.

     3.

     Setze ${\mathbf{y}}=(2u_1-1,\ldots,2u_m-1)^\top$.

   
   

2. Normalverteilung 
   • Als eine Alternative zu dem in Abschnitt 3.2.1 diskutierten Box-Muller-Algorithmus betrachten wir noch ein anderes Verfahren zur Erzeugung von normalverteilten Pseudozufallszahlen,
     • das in der Literatur die Polar-Methode genannt wird,
     • wobei nun die Berechnung der trigonometrischen Funktionen in (12) vermieden wird.
   • Der Zufallsvektor $(V_1,V_2)$ sei im Einheitskreis $B=\{(x_1,x_2)\in\mathbb{R}^2:\,x_1^2+x_2^2\le 1\}$ gleichverteilt.
   • Dann ist der Zufallsvektor $(Y_1,Y_2)$ mit
     $$Y_1=\sqrt{-2 \log(V_1^2+V_2^2)}\; \frac{V_1}{\sqrt{V_1^2+V_2^2}}\,,\qquad Y_2=\sqrt{-2 \log(V_1^2+V_2^2)}\; \frac{V_2}{\sqrt{V_1^2+V_2^2}}$$
     N $({\mathbf{o}},{\mathbf{I}})$-verteilt, d.h. $Y_1$, $Y_2$ sind unabhängige und N$(0,1)$-verteilte Zufallsvariablen,
     • denn mit der Substitution
       $$v_1=r\cos \theta\,,\qquad v_2=r\sin \theta\,,$$

     • d.h. durch Übergang zu Polarkoordinaten, ergibt sich für beliebige $y_1,y_2\in\mathbb{R}$
       

       $${P(Y_1\le y_1,\, Y_2\le y_2)}$$

       $$= \frac{1}{\pi}\;\int\limits_B {1\hspace{-1mm}{\rm I}}\Bigl(\frac{v... ...v_2 \sqrt{-2 \log(v_1^2+v_2^2)}}{\sqrt{v_1^2+v_2^2}}\,\le y_2\Bigr)\,d(v_1,v_2)$$

       $$= \frac{1}{\pi}\;\int\limits_0^{2\pi}\int\limits_0^1 r\,{1\hspace{-... ...\theta\,\le y_1\,,\, \sqrt{-2 \log(r^2)}\; \sin\theta\le y_2\Bigr)\,dr\,d\theta$$

       $$= \frac{1}{2}\;\frac{1}{2\pi}\;\int\limits_0^{2\pi}\int\limits_0^\i... ...}\cos\theta\le y_1,\, \sqrt{x}\sin\theta\le y_2\Bigr)\;e^{-x/2}\,dx\,d\theta\,,$$

       
       

     • wobei die letzte Gleichheit das Ergebnis der folgenden Substitution ist:
       $$x=-2 \log(r^2)$$   bzw.$$\qquad -\;\frac{1}{2}\; e^{-x/2}\,dx=2r\,dr\,.$$

     • Hieraus ergibt sich nun die Behauptung (genauso wie bei der Begründung von Formel (12) in Abschnitt 3.2.1).
   • Die Pseudozufallszahlen $y_1,\ldots,y_{2n}$ mit
     $$y_{2k-1}=\sqrt{-2 \log(v_{2k-1}^2+v_{2k}^2)}\; \frac{v_{2k-1}}{\s... ...sqrt{-2 \log(v_{2k-1}^2+v_{2k}^2)}\; \frac{v_{2k}}{\sqrt{v_{2k-1}^2+v_{2k}^2}}$$
     • können somit als Realisierungen von unabhängigen und N$(0,1)$-verteilten Zufallsvariablen aufgefasst werden,
     • falls $(v_1,v_2),\ldots,(v_{2n-1},v_{2n})$ die Realisierungen von unabhängigen und im Einheitskreis $B=\{(x_1,x_2)\in\mathbb{R}^2:\,x_1^2+x_2^2\le 1\}$ gleichverteilten Zufallsvektoren $(V_1,V_2),\ldots,(V_{2n-1},V_{2n})$ sind,
     • die (wie im vorhergehenden Beispiel beschrieben) mit der Akzeptanz- und Verwerfungsmethode erzeugt werden können.