Bestimmung von U(0,1)-verteilten Zufallszahlen

• Lineare Kongruenzgeneratoren (LCG):
  • gegeben seien natürliche Zahlen m, a, c und Startwert X₀ (a, c, X₀ jeweils < m)
  • berechne Folge ganzer Zahlen nach
    • X_n+1 = (a X_n + c) mod m
  • zur Erinnerung: a mod m = Rest von a bei Teilen durch m
  • Zufallszahl zwischen 0 und 1 dann
    • U_n = X_n / m
  • Beispiel mit kleinen Zahlen (m = 64, a = 9, c = 1, X₀ = 42)
    • 
• Eigenschaften von LCGs:
  • Werte immer rational mit Nenner m
  • Folge wiederholt sich mit gewisser Periode ≤ m
    • maximale Periode m möglich (für gewisse Werte von a, c, m)
  • also: m möglichst groß und Zweierpotenz (→ mod-Operation sehr schnell)
  • einige Beispielwerte

    • m	a	c	Quelle
      231	1103515245	12345	Ansi-C
      224	1140671485	12820163	MS Visual Basic (bis Version 6)
      248	25214903917	11	Java (einfache Version)
      231	65539	0	RANDU (bis in die 70er weit verbreitet, aber ganz schlecht!)

• Tests der Zufallszahlen:
  • viele Tests denkbar und wichtig
  • sehr ausführliche und vielseitige Testsammlung [11]
  • genaue statistische Analyse braucht Testtheorie (s.u.)
  • Fehlerniveau α bedeutet
    • Wenn die zu testende Vermutung (z. B. "Zufallszahlen gleichverteilt") abgelehnt wird, ist mit Wahrscheinlichkeit α die Vermutung doch wahr.
  • i. f. zwei konkrete Beispiele
• Test auf Gleichverteilung:
  • Vorgehensweise
    • erzeuge n Zufallszahlen U_i (i. F. n=10⁴)
    • zerlege [0,1] in k=10 gleich große Teilintervalle
    • f_j (j = 1 .. k) ist Anzahl der Werte im j-ten Teilintervall
    • berechne
    • 
  • Idee: χ² sollte klein sein
  • konkret für α = 0.01: Annahme ablehnen, wenn
    • 
  • ist man mit 5% Fehlerwahrscheinlichkeit zufrieden, reicht
    • 
  • zur Erinnerung: χ²_f = Chi-Quadrat-Verteilung mit f Freiheitsgraden
  • konkret für Visual-Basic-RNG mit X₀ = 42
    • χ² = 8.92
    • Verteilung auf die Teilintervalle
    • 
• Test auf Gleichverteilung im Raum:
  • Idee
    • fasse je drei aufeinander folgende U_i zu einem 3d-Vektor V_i zusammen
    • prüfe die Gleichverteilung der V_i im Einheitswürfel
  • konkretes Beispiel
    • erzeuge n=10⁴ Zufallsvektoren mit dem RANDU-Generator
    • zerlege das Einheitsinterval in k = 10 gleich große Teile
    • → der Einheitswürfel wird in k³ = 1000 Teilwürfel zerlegt
    • bestimme die Zahl f_i,j,l der Vektoren im Teilwürfel (i,j,l)
    • berechne
    • 
  • Ergebnis
    • bei seed = 42 erhält man χ² = 1175.2
    • Tabellenwerte:

    • 1-α	χ2999
      0.99	1105.9
      0.9999	1173.9

    • Die Annahme der Gleichverteilung wird abgelehnt mit einer Fehlerwahrscheinlichkeit < 0.01%
  • Ursache
    • Vektoren V_i sind auf nur 15 Ebenen im Raum verteilt!
    • 
• LFSR-Generator (Linear Feedback Shift Register):
  • erzeuge Folge von Bits folgendermaßen
    • q Bits werden in Schieberegister gespeichert
    • wähle feste Position r, berechne b = b_r ⊕ b_q (⊕ = Summe mod 2 = XOR)
    • schiebe die Bits nach rechts und speichere b im ersten Register
    • gib jeweils Wert in b₁ aus
  • Beispiel q = 5, r = 3
    • 
    • Startwert 11111 liefert Werte der Bits

    • 11111
      01111
      00111
      00011
      10001
      11000
      01100
      10110

  • 1. Bit liefert Ausgangsfolge: 1 0 0 0 1 1 0 1
  • daraus Zahl in [0,1]
    • fasse jeweils l Bits zu einer Binärzahl W_i zusammen
    • Zufallszahl U_i = W_i/2^l
  • maximale Periode 2^q-1 (alle Zahlen außer 0) für geeignete q, r und l
• Erweiterungen:
  • GFSR (Generalized LFSR)
    • arbeitet auf ganzen l-Bit-Worten
    • X_i = X_i-r ⊕ X_i-q
    • ⊕ bitweise berechnen
  • Beispiel mit q = 5, r = 3, l = 4 und den ersten 5 Werten als Startwerten

    • 1	1011
      2	1010
      3	1010
      4	1100
      5	1110
      -
      6	0001
      7	0110
      8	0100
      9	1101
      10	1000

  • Problem: Güte hängt stark von den Anfangswerten ab
  • TGFSR (Twisted GFSR)
    • A sei binäre lxl-Matrix, X_i werde als l-Vektor aufgefasst, damit
    • X_i = X_i-r ⊕ (A*X_i-q)
    • A*X wird ebenfalls mod 2 gerechnet
  • Mersenne-Twister
    • spezieller TGFSR-Generator mit ein paar zusätzlichen Tricks
    • leicht zu initialisieren
    • Periode 2¹⁹⁹³⁷ - 1
    • besteht "Würfeltest" in bis zu 623 Dimensionen
  • weit verbreitet
    • Programmiersprachen (C++, PHP, Python, ...)
    • numerische Programme (Matlab, Maple, Scilab, Excel, ...)
    • Statistik-Software (SPSS, SAS, R, ...)
• Weitere Verfahren:
  • RNGs für mehrere Streams
    • erzeugen viele unabhängige Zufallsfolgen
  • kryptographisch sichere RNGs
    • Parameter des zugrundeliegenden RNG lassen sich aus den Zufallszahlen nicht erschließen (in endlicher Zeit)
  • Vorsicht: fehlerhafte RNGs haben schon manche Forschungsergebnisse ruiniert!
    • RNG in Excel vor 2003 fehlerhaft, sein Ersatz bis 2007 ebenfalls fehlerhaft [12]
    • seit 2010 Mersenne-Twister, aber in VBA immer noch defekter alter RNG [13]
• RNGs in Matlab:
  • Grundfunktion rand(N,M)
    • erzeugt NxM-Matrix aus U(0,1)-Zufallszahlen
    • verwendet standardmäßig Mersenne-Twister
  • Startwert festlegen
    • rng(seed) mit seed ≥ 0 ganzzahlig
    • rnd('shuffle') verwendet Uhrzeit zur Initialisierung
  • RNG-Typ ändern mit rnd(type)

    • 'twister'	Mersenne Twister
      'simdTwister'	SIMD-oriented Fast Mersenne Twister
      'combRecursive'	Combined Multiple Recursive (viele parallele Streams)
      'multFibonacci'	Multiplicative Lagged Fibonacci (viele parallele Streams)

• Aufgaben:
  • Aufgabe 19