forked from manpen/netsci
-
Notifications
You must be signed in to change notification settings - Fork 0
/
erdos_renyi.tex
1355 lines (1126 loc) · 80.6 KB
/
erdos_renyi.tex
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
% !TeX root = skript.tex
Ein Zufallsgraph ist ein Modell einer Graphfamilie -- etwa wie eine Zufallsvariable ein Modell einer Zufallsgröße ist.
Der Zufallsgraph beschreibt also ein Ensemble von Graphen, das so entworfen sein kann, dass bestimmte Eigenschaften in den betrachteten Graphen verstärkt auftreten.
Solche maßgeschneiderte Zufallsgraphen werden oft als \emph{Netzwerkmodell} bezeichnet.
Das Forschungsfeld der Zufallsgraphen nahm mit den Arbeiten von Gilbert~\cite{gilbert_1959}, sowie Erd\H{o}s und R\'enyi~\cite{erdos_renyi_1960} Anfang der 1960er Jahre an Fahrt auf.
Die beiden Arbeiten definieren Zufallsgraphen, die zunächst für abstrakte Untersuchungen (z.\,B. die probabilistische Methode) verwendet wurden.
Der Fokus auf die Modellierung von beobachteten Netzwerken kam erst später hinzu.
Dennoch spielen die Modelle in der Netzwerkforschung bis heute eine wichtige Rolle.
Wir werden sie daher bald genauer betrachten, fangen jedoch mit einem noch einfacheren Modell an.
\section{Der Zufallsgraph $\Gn$}
Ein Zufallsgraph $(\mathbb G, f)$ ist eine Wahrscheinlichkeitsverteilung $f\colon \mathbb G \to [0, 1]$ über einer Menge von Graphen $\mathbb G$.
Oftmals wird der Grundraum $\mathbb G$ durch eine Parametrisierung eingeschränkt; auch $f$ kann parametrisiert sein.
Das \aside{$\Gn$-Graphen} einfachste Modell ist $\Gn$, welches die Gleichverteilung über alle Graphen mit $n$ Knoten beschreibt.
Hierbei ist zu beachten, dass wir mit \emph{alle Graphen} in der Regel entweder alle \emph{gerichteten} oder alle \emph{ungerichteten} Graphen meinen.
Die Details der Analysen in den folgenden Kapiteln hängen von dieser Entscheidung ab -- allerdings ergeben sich meist keine qualitativen Unterschiede.
Daher werden wir oft den für uns einfacheren Fall wählen.
Je nach Wahl ist der Grundraum von $\Gn$
\begin{align}
\mathbb G_\text{ger}(n) & =
\twoset{G = (V,E)}{|V| = n \land E \subseteq V \times V} \\\label{eq:gerichtet_gn}
\mathbb G_\text{unge}(n) & =
\twoset{G = (V,E)}{|V| = n \land E \subseteq \twoset{\set{u,v}}{u,v \in V \text{ mit } u \ne v}}.
\end{align}
\noindent Die Wahrscheinlichkeitsverteilung $f$ folgt dann als $f_{\gGn}(G) = 1 / | \gGn |$ oder konkret:
\begin{align}
f_{\mathbb G_\text{ger}(n)}(G) & = \frac{1}{| \mathbb G_\text{ger}(n) |} = 2^{-n^2}\label{eq:gleichverteilt_gerichtet_gn} \\
f_{\mathbb G_\text{unge}(n)}(G) & = \frac{1}{| \mathbb G_\text{unge}(n) |} = 2^{-\binom n 2} = 2^{-\frac{n(n-1)}{2}}.
\end{align}
Die geschlossene Form $| \mathbb G_\text{ger}(n) | = 2^{n^2}$ in \cref{eq:gleichverteilt_gerichtet_gn} ergibt sich daraus, dass ein gerichteter Graph $n^2$ potentielle Kanten hat (die Adjazenzmatrix hat $n \times n$ Einträge).
Für jede dieser Kanten haben wir unabhängig genau zwei Möglichkeiten: Sie existiert oder eben nicht.
Analog sieht es für die $\binom n 2$ möglichen Kanten in einem ungerichteten Graphen aus.
Beachte aber auch, dass formal \cref{eq:gerichtet_gn} und \cref{eq:gleichverteilt_gerichtet_gn} widersprüchlich scheinen.
In \cref{eq:gerichtet_gn} fordern wir nur, dass die Knotenmenge~$V$ die Kardinalität $|V| = n$ hat.
Offensichtlich gibt es unbeschränkt viele dieser Wahlen; für $n=3$ z.\,B. $\set{1,2,3}$, $\set{2,3,4}$, \ldots, $\set{k, k+1, k+2}$ für alle $k \in \mathbb N$.
Die Wahl der Knotenbezeichner hat jedoch keinen Einfluss auf die Netzwerkeigenschaften.
Daher werden wir in dieser Veranstaltung grundsätzlich annehmen, dass die Knotenmengen in irgendeiner Form fixiert sind, z.\,B. als $V = \set{1, \ldots, n} = [n]$ oder $V = \set{v_1, \ldots, v_n}$.
\bigskip
Ist das $\Gn$-Modell aber nun realistisch?
Das hängt davon ab, was wir mit ihm bezwecken und was wir mit \emph{realistisch} meinen.
Es ist aber sicherlich nicht geeignet, um gängige Netzwerke zu beschreiben.
Dies liegt unter anderem daran, dass $\Gn$ oft Graphen mit vielen Kanten erzeugt;
intuitiv hat \glqq jeder zweite Graph\grqq{} mindestens die Hälfte aller möglichen Kanten.
\begin{figure}[t]
\begin{center}
\includegraphics[width=0.5\textwidth]{data/network-rep-edges.pdf}%
\includegraphics[width=0.5\textwidth]{data/network-rep-edges-thin.pdf}%
\end{center}
\caption{
Die Verteilung der Kantenanzahl in verschiedenen Netzwerken auf~\cite{networkrepository}.
Die Farben der einzelnen Punkte geben an, aus welchem Bereich die Netzwerke stammen.
}
\label{fig:kantenanzahl}
\end{figure}
\begin{observation}
Sei $G = (V,E)$ ein Graph, der zufällig aus $\Gn$ mit $n > 1$ gezogen wurde.
Dann gilt für gerichtete Graphen $\prob{|E| \ge n^2 / 2} \ge 1/2$ und für ungerichtete Graphen $\prob{|E| \ge \binom{n}{2} / 2} \ge 1/2$.
\end{observation}
\begin{proof}
Im Folgenden betrachten wir nur gerichtete Graphen; der Beweis läuft analog für ungerichtete Graphen.
Stellen wir uns einen beliebigen Graphen~$G = (V, E)$ vor.
Dann sei $\bar G = (V, \bar E)$ sein Komplement, d.\,h. für alle möglichen Kanten gilt:
\begin{equation}
\forall e \in V\times V\colon \quad\quad e \in \bar E \Leftrightarrow e \notin E
\end{equation}
Beobachte, dass es eine Bijektion zwischen allen Graphen in $\mathbb G$ und ihren Komplementen gibt: Jeder Graph hat ein eineindeutiges Komplement.
Per Konstruktion gilt außerdem:
\begin{eqnarray}
E \cup \bar E &=& V \times V\\
\Rightarrow |E| + |\bar E| &\ge& n^2\\
\Rightarrow \max(|E|, |\bar E|) &\ge& n^2 / 2.
\end{eqnarray}
Für jeden Graph~$G$ gilt also, dass entweder $G$ selbst oder sein Komplement~$\bar G$ mindestens $n^2 / 2$ Kanten hat.
Da wir $G$ und $\bar G$ je mit gleicher Wahrscheinlichkeit ziehen, gilt $\prob{|E| \ge n^2 / 2} \ge 1/2$.
\end{proof}
\section{Kantenanzahl in beobachteten Netzwerken}\label{sec:kanten-in-beobachteten-netzen}
Wie viele Kanten haben echte Netzwerke? Hierzu führen wir ein Experiment durch:
Wir nutzen die Datenbank \url{https://networkrepository.com/}, die über 5000 Netzwerke aus unterschiedlichen Bereichen enthält~\cite{networkrepository}.
In \cref{fig:kantenanzahl} zeichnen wir die Kantenanzahl als Funktion der Knotenanzahl.
Zwei Eigenschaften fallen direkt auf:
\begin{enumerate}
\item Die Kantendichte hängt vom Netzwerktyp ab.
\item Es sind fast immer deutlich weniger als die Hälfte der Kanten vorhanden.
\end{enumerate}
\subsection{Netzwerktypen haben unterschiedliche Kantendichten}
Betrachten wir die erste Beobachtung genauer, indem wir Straßennetze und Freundschaftsnetze vergleichen.
Wir modellieren ein Straßennetz dadurch, dass Adressen (Häuser, Kreuzung usw.) als Knoten und Straßen als Kanten dargestellt werden.
Diese Netze \aside{Straßennetze} sind im wesentlichen ein zweidimensionales Konstrukt.
Wenn wir Tunnel, Brücken und dergleichen ignorieren, verlaufen Straßen nicht über einander.
Daher erwarten wir, dass die Graphen von Straßennetzen fast planar sind.
Nach dem Euler'schen Polyedersatz erfüllen einfache, \aside{planare Graphen} planare und zusammenhängende Graphen:
\begin{equation}
|E| \le 3 |V| - 6
\end{equation}
Knoten in einem Straßennetz sollten also im Schnitt höchstens sechs Nachbarn haben.
In \aside{Freundschaftsnetze} sozialen Netzwerken ist die Situation anders --
stellen wir uns etwa einen Freundschaftsgraphen vor, in dem Knoten die Nutzer eines sozialen Netzwerks sind und Kanten eine Freundschaft anzeigen.
Im Jahre 2014 hatten Facebook-Nutzer im Schnitt mehr als 300 Freunde (mehr dazu später).
Dies ist offensichtlich deutlich mehr als in planaren Graphen möglich wäre.
Ganz ähnlich sieht es mit anderen sozialen Netzen aus: Ein durchschnittlicher Erwachsener kennt deutlich mehr als sechs andere Menschen persönlich (oft werden Zahlen zwischen 100 und 300 genannt).
\subsection{Die meisten Netzwerke sind dünn}
In \cref{fig:kantenanzahl} hat nur ein verschwindend geringer Anteil der Netzwerke mindestens die Hälfte aller Kanten (d.\,h. ist oberhalb der roten Linie).
Wie erklärt sich das?
In der Regel verursacht eine Kante Kosten:
Eine Straße muss gebaut werden, eine Freundschaft muss aufrecht erhalten werden (Zeitinvestment), eine Nachricht muss geschrieben werden etc.
Daher gibt es in den meisten Netzwerken einen gewissen Selektionsdruck, der dazu führt, dass jeder Knoten nur ausgewählte Nachbarn besitzt.
Wir klassifizieren Netzwerktypen, die auffallend viele oder wenig Kanten haben:
\begin{definition}
Graphen \aside{dünne und dichte Graphen} mit $n$ Knoten und $m$ Kanten heißen \emph{dünn} (engl. sparse), wenn $m = \Oh{n \log n}$ gilt, und \emph{dicht} wenn $m = \Omega(n^2)$.
\end{definition}
\begin{remark}
Aufgrund der asymptotischen Definition verwendet man \emph{dünn} und \emph{dicht} eher für Graphfamilien (so auch Zufallsgraphen) als für einzelne Instanzen.
Manche abweichende Definitionen fordern für \emph{dünne} Graphen $m = \Oh{n}$.
Wir nutzen hier $\Oh{n\log n}$, um soziale Netzwerke mit einzubeziehen.
\end{remark}
\section{Die Zufallsgraphen $\Gnm$ und $\Gnp$}
Das $\Gn$-Modell verfügt über keinen Mechanismus um Kanten auszudünnen.
Wir benötigen also Prozesse, die weniger dichte Graphen erzeugen können --
am besten parametrisiert, damit wir unterschiedlichen Netzwerktypen Rechnung tragen können.
Im Folgenden betrachten wir zwei solcher Modelle.
Das \Gnm-Modell \aside{Erd\H{o}s-R\'enyi Graphen \Gnm} von P.~Erd\H{o}s und A.~R\'enyi beschreibt die Gleichverteilung über allgemeinen Graphen mit $n$ Knoten und $m$ Kanten, d.\,h. wir betrachten die Grundmenge
\begin{equation}
\gGnm = \twoset{G}{G=(V,E) \in \mathbb G(n) \textnormal{ und } |E| = m}.
\end{equation}
Da gleichverteilt gewählt wird, gilt für die Wahrscheinlichkeitsverteilung~{$f\colon \mathbb G \to [0,1]$}
\begin{equation}
f(G) = \frac{1}{| \gGnm |} \quad \text{ für alle } G \in \gGnm.
\end{equation}
\begin{exercise}
Berechne $|\gGnm|$ für gerichtete und ungerichtete Graphen.
\end{exercise}
E.~Gilbert \aside{Gilbert-Graphen \Gnp} beschreibt ein ähnliches Modell, das \Gnp-Modell -- das oft fälschlicherweise \glqq Erd\H{o}s-R\'enyi-Modell \grqq{} genannt wird.
Um das Modell zu beschreiben, weichen wir von der bisherigen expliziten Definition der Grundmenge und Verteilung ab.
Stattdessen spezifizieren wir eine randomisierte Konstruktionsvorschrift:
\begin{enumerate}
\item Erzeuge $n$ Knoten $V = \{v_1, \ldots, v_n\}$.
\item Setze $E = \emptyset$.
\item Für jedes Paar von Knoten $(v_i, v_j)$ führe ein unabhängiges Bernoulli Experiment durch.
Mit Wahrscheinlichkeit $p$ füge $(v_i, v_j)$ als Kante zu $E$ hinzu.
\item Gebe den Graphen $G=(V, E)$ zurück.
\end{enumerate}
\noindent
Graphisch kann man sich also \Gnp als Adjazenzmatrix vorstellen, in der die Einträge unabhängig voneinander mit Wahrscheinlichkeit $p$ auf $1$ gesetzt werden:
\begin{center}
\begin{tikzpicture}
\node (mat) at (0,0) {
$\begin{pmatrix}
1 & 0 & 1 & 0 & 0 \\
0 & 1 & 1 & 0 & 1 \\
0 & 1 & 0 & 1 & 0 \\
1 & 0 & 1 & 0 & 0 \\
1 & 0 & 0 & 0 & 1 \\
\end{pmatrix}$
};
\node[anchor=west, align=left, xshift=4em] (label) at (mat.east) {
$\begin{cases*}
1 & mit Wahrscheinlichkeit $p$ \\
0 & mit Wahrscheinlichkeit $1-p$
\end{cases*}$
};
\path[draw, thick, bend right, ->] (label.west) to ($(mat.center)+(0.45em,0)$);
\end{tikzpicture}
\end{center}
\begin{exercise}
Die genannte Konstruktion erzeugt gerichtete Graphen.
Zeige, wie die Konstruktion so angepasst werden kann, dass gerichtete Graphen ohne Eigenschleifen erzeugt werden.
Wie verändert sich dann die Adjazenzmatrix?
Wie verhält es sich mit ungerichteten Graphen?
\end{exercise}
Durch ihre einfache Konstruktion sind beide Zufallsgraphen bis heute sehr verbreitete Modelle in der Netzwerkforschung.
Wir werden jedoch sehen, dass viele Eigenschaften von echten Netzwerken auch von \Gnp- oder \Gnm-Graphen nicht beschrieben werden können.
\subsection{Anzahl von Kanten in \Gnp}\label{subsec:anzahl_kanten_in_gnp}
Während bei \Gnm-Graphen die Anzahl der Kanten durch den Parameter~$m$ fixiert ist, ist $|E|$ bei \Gnp Graphen eine Zufallsvariable.
\begin{lemma}\label{lemma:erwartete_kanten_in_gnp}
Die \aside{erwartete Kantenanzahl $\expv{|E|}$ in \Gnp} erwartete Kantenanzahl~$m$ in einem gerichteten \Gnp Graphen ist \begin{equation*} \expv{m} = p n^2. \qedhere \end{equation*}
\end{lemma}
\begin{proof}
Fixiere einen Graphen~$G=(V,E)$ aus \Gnp.
Für jede Kante $(u,v)$ definiere die Indikatorvariable $I_{u,v}$, die anzeigt, ob die Kante $(u,v) \in E$ enthalten ist:
\begin{equation}
I_{u,v} = \begin{cases*}
1 & falls $(u,v) \in E$, \\
0 & sonst.
\end{cases*}
\end{equation}
\noindent Somit folgt die Kantenanzahl~$m$ in $G$ als Summe über die Indikatorvariablen
\begin{equation}
|E| = \smashoperator{\sum_{u,v \in V}} I_{u,v} = \smashoperator{\sum_{(u,v) \not\in E}} 0 \quad + \quad \smashoperator{\sum_{(u,v) \in E}} 1.
\end{equation}
\noindent Per Definition von \Gnp gilt $\prob{I_{u,v} {=} 1} = p$ und $\prob{I_{u,v} {=} 0} = 1-p$.
Somit folgt für den Erwartungswert jeder Indikatorvariable
\begin{equation}
\expv{I_{u,v}} = 1 \cdot p + 0 \cdot (1-p) = p \quad \text{unabhängig für alle } u, v \in V.
\end{equation}
\noindent Durch die Linearität des Erwartungswertes folgt schließlich
\begin{equation}
\expv{|E|} = \smashoperator{\sum_{u,v \in V}} \expv{I_{u,v}} = \smashoperator{\sum_{u,v \in V}} p = p n^2. \qedhere
\end{equation}
\end{proof}
\begin{exercise}
Zeige, dass die erwartete Anzahl an Kanten in einem ungerichteten \Gnp-Graphen $\expv{|E|} = \binom{n}{2} p = p n(n-1)/2$ beträgt.
\end{exercise}
\begin{exercise}
Zeige, dass $\Gn = \mathcal G(n, 1/2)$.
\end{exercise}
\bigskip
Wie wir am Beweis von \cref{lemma:erwartete_kanten_in_gnp} sehen, ergibt sich die Kantenanzahl~$m$ als Summe von unabhängigen Bernoulli-Zufallsvariablen;
sie ist also selbst eine Zufallsvariable und binomial verteilt.
Da uns die Binomialverteilung regelmäßig begegnen wird, wollen wir uns diese kurz in Erinnerung rufen.
\begin{definition}
Die \aside{Binomialverteilung.} Binomialverteilung $B_{N, p}(k)$ beschreibt die Wahrscheinlichkeit, dass bei $N$~unabhängigen Bernoulli Experimenten mit Wahrscheinlichkeit $p$ genau $k$ Experimente erfolgreich sind.
Es gilt \aside{\ \\ \ \\ \ \\ Binomialkoeffizient $\binom{n}{k} = \frac{n!}{(n-k)!k!}$}
\begin{eqnarray*}
\prob{B_{N,p} {=} k} &=& \binom{N}{k} p^k (1-p)^{N-k} \\
\expv{B_{N,p}} &=& Np \\
\varv{B_{N,p}} &=& Np(1-p),
\end{eqnarray*}
wobei $\expv{B_{N,p}}$ und $\varv{B_{N,p}}$ den Erwartungswert und die Varianz der Binomialverteilung beschreiben.
\end{definition}
Die Standardabweichung $\sigma = \sqrt{\varv{B_{N,p}}} \le \sqrt{\expv{B_{N,p}}}$ ist also relativ klein.
Wir können daher davon ausgehen, dass für hinreichend großes $N$ die Binomialverteilung recht stark um ihren Erwartungswert $Np$ konzentriert ist.
Daher sagen wir, dass \Gnm und \Gnp mit $p=m/n^2$ asymptotisch (d.\,h. für $n \to \infty$) äquivalent sind.
Häufig ist es jedoch einfacher \Gnp-Graphen zu analysieren, da -- im Gegensatz zu \Gnm-Graphen -- alle Kanten unabhängig gezogen werden.
\section{Effizientes Ziehen von \Gnp-Graphen}
Zufallsgraphen sind nicht nur in der theoretischen Analysen von Prozessen und Algorithmen nützlich, sondern auch für empirische Untersuchungen.
Nehmen wir etwa an, dass wir einen Algorithmus implementiert haben und dessen Geschwindigkeit oder Qualität vermessen möchten.
Dann können Zufallsgraphen unter anderem aus folgenden Gründen nützlich sein:
\begin{itemize}
\item Wenn die Generatorsoftware vorhanden ist, können wir synthetische Instanzen in quasi unbeschränkter Menge generieren.
\item Parametrisierte Generatoren erlauben, den Einfluss gewisser Eigenschaften detailliert zu analysieren, z.\,B. Skalierungsexperimente.
\item Beobachtete Graphen haben oft \glqq Rauschen\grqq, d.\,h. nicht verstandene oder insignifikante Strukturen, die Messungen verfälschen können.
Zufallsgraphen hingegen haben i.\,d.\,R. eine gut verstandene Struktur, die es uns oft ermöglicht, vor dem Testen Hypothesen aufzustellen.
\item Einige Generatoren können anhand eines fixierten Startzustands (Random-Seed) dieselben Graphen wiederholt erzeugen.
Wir müssen also die Eingaben nicht speichern/transferieren und können dennoch reproduzierbare Experimente durchführen.
\end{itemize}
\noindent Da in verschiedenen wissenschaftlichen Communities der Begriff \glqq Generator\grqq{} unterschiedlich benutzt wird, hier noch eine Definition:
\begin{definition}
Ein \aside{Graphgenerator} Graphgenerator ist ein Algorithmus, der einen Graphen aus einem Netzwerkmodell erzeugt.
Dieser respektiert die Verteilung des zugrundeliegenden Modells.
\end{definition}
\subsection{Ein naiver Ansatz}
Betrachten wir folgenden Graphgenerator für \Gnp-Graphen, der im wesentlichen eine Eins-zu-eins-Implementierung der Definition von \Gnp ist:
\begin{algorithm}[H]
\KwIn{Anzahl der Knoten~$n$ und Verbindungswahrscheinlichkeit~$p$}
\KwOut{Adjazenzmatrix eines zufälligen Graphen $G \follows \Gnp$}
Allokiere eine Matrix $A[1..n, 1..n]$\;
\For{$1 \le i \le n$}{
\For{$1 \le j \le n$}{
Setze $A[i,j] \gets \begin{cases*}
0 & mit Wahrscheinlichkeit $1-p$ \\
1 & mit Wahrscheinlichkeit $p$
\end{cases*}$\;
}
}
Gebe $A$ zurück
\caption{Naiver Graphgenerator für \Gnp-Graphen}
\label{alg:naive-gnp}
\end{algorithm}
Um \aside{Annahme: Aus einfachen Verteilungen kann in konstanter Zeit gezogen werden.} die Laufzeit dieses Algorithmus bestimmen zu können, treffen wir in dieser Veranstaltung folgende Annahme.
Wir können unter anderem aus folgenden Zufallsverteilungen in konstanter Zeit ziehen:
\begin{itemize}
\item Ganzzahl uniform aus $[a, b]$ für $a, b \in \mathbb{Z}$
\item Gleitkommazahl uniform aus $[a, b]$ für $a, b \in \mathbb{R}$
\item Bernoulli (folgt aus vorherigem Punkt)
\item Binomial, Geometrisch (siehe \cref{{sec:inversionsmethode}}), Hypergeometrisch, Poisson
\end{itemize}
Tatsächlich können die meisten dieser Verteilungen nur in \emph{erwartet} konstanter Zeit gezogen werden, praktisch sind die Fluktuationen jedoch so klein, dass wir sie vernachlässigen können.
Das erlaubt es uns, uns auf die Eigenschaften der Graphgeneratoren zu konzentrieren.
\begin{lemma}
\label{lem:naive-gnp}
Der Algorithmus \ref{alg:naive-gnp} erzeugt einen \Gnp-Graphen in $\Theta(n^2)$ Zeit.
\end{lemma}
\begin{exercise}
Beweise das Lemma.
\end{exercise}
Nun könnte man argumentieren, dass dieser naive Algorithmus bereits optimal ist:
Die Ausgabe hat Größe $\Theta(n^2)$ und somit benötigen wir $\Omega(n^2)$ Zeit, um die Ausgabe zu erzeugen.
Außerdem hat ein Graph mit $p=\Omega(1)$ in Erwartung $\Theta(n^2)$ Kanten, was wiederum die Laufzeit von unten beschränkt.
Diese Worst-Case-Schranken passen jedoch nicht zu unserer Beobachtung in \cref{sec:kanten-in-beobachteten-netzen}, dass Graphen in der freien Wildbahn meist dünn sind.
Es wäre schön, wenn wir diese Graphen schneller erzeugen könnten, jedoch ist \cref{alg:naive-gnp} für alle $p$ gleich \glqq schnell\grqq.
Hier \aside{ausgabesensitive Algorithmen} helfen \emph{ausgabesensitive Algorithmen} (output-sensitive algorithm), deren Laufzeit maßgeblich durch die Größe der Ausgabe bestimmt wird;
die meisten Generatoren, die wir betrachten werden, fallen genau in diese Kategorie.
Ausgabesensitive \aside{Kantenlisten} Algorithmen müssen jedoch eine Ausgabedatenstruktur nutzen, die sich an die Größe anpasst.
Eine Matrix hat immer dieselbe Größe, unabhängig davon, wie viele Einsen (Kanten) vorhanden sind.
Wir werden im Folgenden fast immer eine Kantenliste annehmen, die nur $\Oh{1}$ Speicherworte pro Kante benötigt.
\subsection{Ziehen mit zufälligen Sprüngen}\label{subsec:gnp_zufaellige_spruenge}
Wie können wir nun \cref{alg:naive-gnp} verändern, so dass er -- analog zu Kantenlisten -- keine Arbeit für Nichtkanten verrichtet?
Idee: Wir müssen diese einfach überspringen!
Aber woher wissen wir, wo \glqq Nullen\grqq{} gezogen werden, bzw. wie viele Nullen wir überspringen müssen?
Zählen können wir sie nicht, ohne wieder direkt quadratische Arbeit zu verrichten.
Ganz \aside{effizientes Ziehen durchs Überspringen von Nullen} einfach: Wir fragen uns, wie viele Nullen~$\ell$ \cref{alg:naive-gnp} gesetzt hätte und überspringen diese.
Zwischen diese zufälligen Sprünge setzen wir dann die Einsen ein.
Wenn wir nun $\ell$ aus einer geeigneten Verteilung ziehen, haben wir zwei unterschiedliche Zufallsprozesse, deren Ausgabe aber nicht unterscheidbar ist.
Generalisierungen dieser Methode werden uns noch häufiger begegnen.
Zur \aside{$n^2$-dimensionaler Vektor statt $n \times n$-Matrix} Vereinfachung stellen wir uns die Adjazenzmatrix~$A = (a_{ij})_{ij}$ als einen Vektor~$B = (b_k)_k$ vor, der dadurch entsteht, dass wir die Zeilen hinter einander anfügen (row-major layout).
Eintrag $a_{ij}$ wird also mit $b_{(i-1)n + j}$ identifiziert.
Dann besteht unsere Aufgabe darin, genau alle $b_k$ zu finden, die den Wert $1$ haben.
\begin{equation}
\Big(\quad \cdots \overbrace{1}^{b_k} \ \underbrace{0 \ \ 0\ \ 0 \ \cdots\ \ 0\ \ 0\ \ 0}_{\text{$\ell$ zusammenhängende Nullen}}\ \overbrace{1}^{b_{k + 1+ \ell}} \cdots \quad \Big)
\end{equation}
\noindent
Sei $X$ die \aside{Anzahl Nullen $X$ geometrisch verteilt} Zufallsvariable, welche die Anzahl der Nullen nach einer Eins modelliert.
\begin{eqnarray}
\prob{X = 0} &=& \prob{\textnormal{Nächster Versuch: 1}} = p \\
\prob{X = 1} &=& \prob{\textnormal{Nächster Versuch: 0, dann 1}} = (1 - p)p \\
\prob{X = 2} &=& \prob{\textnormal{Nächster Versuch: 0, dann 0, dann 1}} = (1 - p)^2 p \\
&\vdots& \\
\prob{X = \ell} &=& (1 - p)^\ell p
\end{eqnarray}
\noindent Die Zufallsvariable ist also \emph{geometrisch verteilt} mit Parameter $p$.
Es ergibt sich also der folgende einfache Algorithmus:
\begin{algorithm}[H]
\KwIn{Anzahl der Knoten~$n$ und Verbindungswahrscheinlichkeit~$p$}
\KwOut{Kantenliste $E$ eines zufälligen Graphen $G \follows \Gnp$}
Initialisiere eine leere Kantenliste $E$\;
$k \gets -1$\;
\While{true}{
$\ell \gets \text{ziehe die Sprungweite aus einer geometrischen Verteilung mit Parameter $p$}$\;
$k \gets k + \ell + 1$\;
\If{$k < n^2$}{
$(i, j) \gets (\lfloor k / n \rfloor, k \bmod n)$\;
$E \gets E \cup \{(i+1, j+1)\}$\;
} \Else {
Gebe $E$ zurück und beende\;
}
}
\caption{Generator für \Gnp-Graphen mit zufälligen Sprüngen}
\label{alg:linear-gnp}
\end{algorithm}
\begin{lemma}
\label{lem:linear-gnp}
Algorithmus \ref{alg:linear-gnp} erzeugt Graph $G(V, E) \follows \Gnp$ in $\Theta(|E|)$ Zeit.
\end{lemma}
\begin{exercise}
Beweise das Lemma.
\end{exercise}
\begin{exercise}
Passe \cref{alg:linear-gnp} für ungerichtete Graphen an.
\end{exercise}
Jetzt müssen wir nur noch herausfinden, wie wir geometrische Sprünge effizient ziehen können.
Hierzu kommt meist die Inversionsmethode zum Einsatz.
\subsection{Inversionsmethode: Ziehen aus der geometrischen Verteilung}
\label{sec:inversionsmethode}
\begin{figure}[t]
\begin{center}
\includegraphics[width=0.5\textwidth]{data/geometric-distr.pdf}%
\includegraphics[width=0.5\textwidth]{data/geometric-distr-inv.pdf}
\end{center}
\caption{
Wahrscheinlichkeitsverteilung und kumulative Verteilungsfunktion einer geometrischen Verteilung mit Parameter $p = 0.2$ (links).
Die rechte Abbildung zeigt die \glqq Umkehrfunktion\grqq{} $F^{-1}(p)$ der kumulativen Verteilungsfunktion $F(x)$.
}
\label{fig:geometric-distr}
\end{figure}
Die \aside{Inversionsmethode: effizientes Ziehen, falls Umkehrfunktion bekannt} Inversionsmethode (inversion transform method) ist ein allgemeines Verfahren, um eine uniforme Zufallsvariable in eine andere Verteilung zu übersetzen.
Wir betrachten die Methode hier am Beispiel der geometrischen Verteilung.
In \cref{fig:geometric-distr} ist die Wahrscheinlichkeitsverteilung $f(\ell)$ und die kumulative Verteilungsfunktion $F(\ell)$ einer geometrischen Verteilung mit Parameter $p = 0.2$ gezeigt.
%
%Da $f(\ell)$ ist eine Wahrscheinlichkeitsverteilung ist, muss $\sum_{\ell=0}^{\infty} f(\ell) = 1$ gelten.
Betrachten wir konkret die ersten drei Werte der Funktionen:
\begin{center}
\begin{tabular}{l|p{0.2\textwidth}p{0.25\textwidth}p{0.35\textwidth}}
& $\ell = 0$ & $\ell = 1$ & $\ell = 2$ \\\hline\hline
$f(\ell)$ & $p =0.2$ & $(1 {-} p)p = 0.16$ & $(1 {-} p)^2p = 0.128$ \\
$F(\ell)$ & $p =0.2$ & $p{+}(1 {-} p)p = 0.36$ & $p{+}(1 {-} p)p{+}(1 {-} p)^2p = 0.488$
\end{tabular}
\end{center}
\vspace{1em}
Ein Generator sollte also z.\,B. den Wert $0$ mit Wahrscheinlichkeit 20\,\% ausgeben.
Da wir einen deterministischen Algorithmus konstruieren möchten, brauchen wir eine Quelle für den Zufall:
Wir erwarten eine uniforme Zufallsvariable $U \in [0, 1)$ als Eingabe.
Wir können also z.B. prüfen:
\begin{itemize}
\item Wenn $U < 0.2$ ist, dann gebe den Wert $0$ zurück.
Da $U$ uniform verteilt ist, ist die Wahrscheinlichkeit für $U < 0.2$ genau $0.2$.
\item Falls nicht, prüfen wir, ob $U < 0.36 = f(0) + f(1) = F(1)$ ist. Falls ja, geben wir den Wert $1$ zurück.
Da $U$ uniform verteilt ist und wir wissen, dass $U$ nicht kleiner als $0.2$ ist, ist die Wahrscheinlichkeit $\prob{U < 0.36 \mid U > 0.2} = 0.16 = f(1)$.
\item Diesen Prozess setzen wir fort, bis wir das passende Intervall gefunden haben.
\end{itemize}
Die graphische Interpretation ist also die folgende:
In \cref{fig:geometric-distr} (rechts) ist die Umkehrfunktion~$F^{-1}(p)$ der kumulativen Verteilungsfunktion $F(x)$ gezeigt.
Wir werfen nun mittels der zufälligen Eingabe $U$ einen Punkt auf der $x$-Achse und schauen, welchen Wert $F^{-1}(U)$ hat -- dies ist unsere Ausgabe.
Im folgenden Theorem wird die Inversionsmethode formalisiert, wobei wir vereinfachend einige Annahmen treffen.
Die Methode lässt sich aber auch auf andere $\Omega$ (inklusive kontinuierliche) und nicht streng monotone $F(x)$ verallgemeinern.
\begin{theorem}
Sei $X \in \Omega$ eine diskrete Zufallsvariable und $f(x)$ und $F(x)$ ihre Wahrscheinlichkeitsverteilung sowie die kumulative Verteilungsfunktion.
Vereinfachend sei $\Omega$ total geordnet und $F(x)$ streng monoton steigend, sodass es eine Umkehrfunktion $F^{-1}(x)$ mit $F^{-1}(F(x)) = x$ gibt.
Sei $U$ eine uniforme Zufallsvariable aus $[0, 1)$.
Dann ist $X' \coloneq F^{-1}(U)$ eine Zufallsvariable mit der Wahrscheinlichkeitsverteilung $f(x)$.
\end{theorem}
\begin{proof}
\noindent Da $F$ streng monoton ist, existiert die Umkehrfunktion $F^{-1}$ und es gilt:
\begin{equation}\label{eq:inverse-apply-inverse}
F(x) = p
\quad\Leftrightarrow\quad F^{-1}(F(x)) = F^{-1}(p)
\quad\Leftrightarrow\quad x = F^{-1}(p)
\end{equation}
\noindent Daraus folgt dann direkt:
\begin{eqnarray}
\prob{X' \le x} &\stackrel{\text{Inv.-Methode}} = & \prob{F^{-1}(U) \le x} \\
&\stackrel{(\ref{eq:inverse-apply-inverse})}{=}& \prob{U \le F(x)} \\
&\stackrel{\prob{U < z} = z \forall z \in [0, 1]}=& F(x) \\
&\stackrel{\text{Def.} F(X)}=& \prob{X \le x} \hspace{5em}\hfill \qedhere
\end{eqnarray}
\end{proof}
\bigskip
\bigskip
\begin{figure}
\begin{center}
\includegraphics[width=0.5\textwidth]{data/geometric-distr-inv_hat.pdf}
\end{center}
\caption{Die in \cref{eq:inverse-hat} berechnet Umkehrfunktion $\hat F^{-1}(U)$}
\end{figure}
Zurück \aside{Inversionsmethode für geometrische Verteilungen} zur geometrischen Verteilung mit $f(\ell) = (1-p)^\ell p$.
O.\,b.\,d.\,A. sei $0 < p < 1$, da $X$ sonst eine Konstante ist.
Die kumulative Verteilungsfunktion lautet
\begin{equation}
F(\ell)
\ =\ p \sum_{i=0}^\ell (1-p)^\ell
\ \stackrel{(*)}{=} \ p \frac{1 - (1-p)^{\ell+1}}{\underbrace{1 - (1-p)}_{=p}}
\ = \ 1 - (1-p)^{\ell+1},
\end{equation}
wobei wir in $(*)$ die geschlossene Form der $\ell$-ten Partialsumme der geometrischen Reihe nutzen: $\sum_{i=0}^\ell x^\ell = (1 - x^{\ell+1})/(1 - x)$.
Nun berechnen die Umkehrfunktion $\hat F^{-1}(U) = \ell$ indem wir $F(\ell) = U$ setzen und nach $U$ umformen.
Hierbei nehmen wir zunächst an, dass die Funktion stetig wäre (daher schreiben wir $\hat F ^{-1}$ statt $F^{-1}$):
\begin{eqnarray}
U &=& F(\ell) = 1 - (1-p)^{\ell+1} \\
\Leftrightarrow (1-p)^{\ell+1} &=& 1 - U \\
\Leftrightarrow (\ell+1) \log(1-p) &=& \log(1 - U) \\
\Leftrightarrow \hat F^{-1}(U) := \log_{1-p}(1 - U) - 1 &=& \ell \label{eq:inverse-hat}
\end{eqnarray}
Schließlich erhalten wir $F^{-1}(U) = \lceil \hat F^{-1} (U) \rceil$ durch Aufrunden der gerade berechneten Funktion.
Zusammenfassend erzeugt \cref{alg:sample-geometric} eine geometrische Zufallsvariable $X$ mit Parameter $p$ in Zeit $\Oh{1}$:
\begin{algorithm}[H]
\If{$p=0$}{Gebe $\infty$ zurück}
\ElseIf{$p=1$}{Gebe $0$ zurück}
\Else{$U \gets \text{ziehe uniform aus $[0, 1)$}$\;
Gebe $\lceil \log_{1-p}(1 - U) - 1\rceil$ zurück}
\caption{Ziehen einer geometrischen Zufallsvariable}
\label{alg:sample-geometric}
\end{algorithm}
\bigskip
\bigskip
\section{Effizientes Ziehen von \Gnm Graphen}\label{sec:sampling_gnm}
Trotz der strukturellen Ähnlichkeiten zwischen \Gnp und \Gnm Graphen unterscheiden sich ihre Generatoren signifikant.
Grund hierfür ist, dass während \Gnp \emph{in Erwartung} $n^2p$ Kanten erzeugen, müssen wir für \Gnm Generatoren sicherstellen, dass exakt $m$ zufällige Kanten produziert werden.
Ähnlich wie zuvor, reduzieren wir das Ziehen der richtigen Position von \emph{Einsen} in der Adjazenzmatrix auf den eindimensionalen Fall.
Im Folgenden werden wir also folgendes äquivalentes Problem lösen:
Ziehe uniform ohne Zurücklegen $k$ Elemente aus einer Menge $S = \set{s_1, \ldots, s_N}$ mit $|S| = N > k$.
\subsection{Reservoir Sampling}
Reservoir Sampling ist eine allgemeine Technik, die es uns erlaubt $k$ Elemente aus einem Datenstrom zu ziehen.
Hierbei müssen wir anfangs \emph{nicht} wissen, wie viele Elemente insgesamt im Strom sind.\footnote{
Ein praktischer Anwendungsfall sind z.B. Iteratoren oder Generatoren, die in vielen Programmiersprachen genutzt werden.
In Rust ist z.B. \texttt{rand::IteratorRandom::choose\_multiple} mittels Reservoir Sampling implementiert.
}
Im konkreten Fall besteht der Datenstrom also aus allen $s_i$ in beliebiger Reihenfolge.
Vereinfachend können wir also annehmen, dass der Strom mindestens $k$ Elemente enthält.
\begin{algorithm}[H]
\KwIn{Eingabe: Datenstrom~$S$, Stichprobengröße~$k$}
\KwOut{Array $R[1..k]$ mit $k$ zufällig ausgewählten Elementen}
Initialisiere ein anfangs leeres Array $R$ mit Kapazität $k$\;
\While{$|R| < k$}{
$x \gets \text{nächstes Element aus $S$}$\;
Füge $x$ an $R$ an\;
}
Setze $i = k$\;
\While{$S$ nicht leer}{
$x \gets \text{nächstes Element aus $S$}$\;
$i \gets i + 1$\;
$j \gets \text{ziehe uniform aus $[1, i]$}$\;
\If{$j \le k$}{
Ersetze $R[j] \gets x$\;
}
}
Gebe $R$ zurück
\caption{Reservoir Sampling}
\label{algo:reservoir-sampling}
\end{algorithm}
\bigskip
\bigskip
\noindent
In \cref{thm:reservoir-sampling} zeigen wir die Korrektheit von \cref{algo:reservoir-sampling} für nur $k = 1$:
\begin{theorem}\label{thm:reservoir-sampling}
Für $k=1$ liefert \cref{algo:reservoir-sampling} liefert ein uniform gewähltes Element $x\in S$.
\end{theorem}
\begin{proof}
Sei $N$ die Länge des Stroms (ist dem Algorithmus nicht bekannt).
Seien zudem $s_1, \ldots, s_N$ die Elemente in der Reihenfolge, in der sie aus dem Strom gelesen werden;
oBdA gelte $s_i = s_j \ \Leftrightarrow\ i = j$ (d.h. alle Elemente sind unterschiedlich; falls nicht, können wir die Tupel $(s_i, i)$ betrachten).
Sei $R$ die Zufallsvariable, welche die Ausgabe des Algorithmus modelliert.
Dann ist zu zeigen, dass $\prob{R = s_i} = 1 / N$ für alle $1 \le i \le N$.
\begin{eqnarray}
\prob{R = s_i} &=& \prob{\text{$s_i$ wird in Runde~$i$ gewählt}} \cdot \prob{\text{$s_i$ wird nicht ersetzt}} \\
&=& \frac{1}{i} \cdot \underbrace{\prod_{j=i+1}^N \frac{j-1}{j}}_{\text{Teleskop-Produkt}}
= \frac{1}{i} \cdot \frac{i}{N} = \frac{1}{N}
\end{eqnarray}
\noindent Insbesondere gilt also, dass $\prob{R = s_i} = 1/N$ unabhängig von $i$ oder $s_i$ ist.
\end{proof}
\begin{exercise}
Zeige die Korrektheit von \cref{algo:reservoir-sampling} indem du den Beweis von \cref{thm:reservoir-sampling} auf $k \ge 1$ verallgemeinerst.
Dies lässt sich etwa über Induktion über die Datenstromlänge $|S|$ erreichen.
\end{exercise}
Wir können also \Gnm Graphen mittels Reservoir Sampling erzeugen.
Allerdings entspricht jede potentielle Kante einem Element in $S = V \times V$, für das wir je $\Theta(1)$ Arbeit verrichten müssen;
somit folgt eine Gesamtlaufzeit von $\Theta(n^2)$.
\subsection{Ziehen mit Lookups}
Im Folgenden skizzieren wir einen alternativen Ansatz um $k$~Elemente ohne Zurücklegen aus dem Universum~$S$ zu ziehen; dieser basiert auf \cite{batagelj2005efficient}.
Dieser nutzt das Grundgerüst in~\cref{algo:basis-ziehen-ohne-zuruecklegen}, das wir um einige Details erweitern müssen, um eine praktisch schnelle Implementierung zu erreichen:
\begin{algorithm}[H]
\KwIn{Eingabe: Universum~$S$, Stichprobengröße~$k$}
\KwOut{Menge $R$ mit $k$ zufällig ausgewählten Elementen}
Initialisiere $R = \emptyset$\;
\While{$|R| < k$}{
$x \gets $ zufälliges Element aus $S$\;
\If{$x \not\in R$}{
$R \gets R \cup \set{x}$
}
}
Gebe $R$ zurück
\caption{Basisalgorithmus für Ziehen ohne Zurücklegen}
\label{algo:basis-ziehen-ohne-zuruecklegen}
\end{algorithm}
\bigskip
Wir messen die Geschwindigkeit des Algorithmus zunächst nur indirekt, indem wir die Anzahl der Iterationen der \texttt{while}-Schleife zählen:
\begin{lemma}\label{lemma:basis-ziehen-ohne-zuruecklegen-versuche}
\Cref{algo:basis-ziehen-ohne-zuruecklegen} verwendet in Erwartung $N \log(N/(N-k+1)) + \Oh{1}$ Versuche (Iterationen der \texttt{while}-Schleife) um aus $S$ mit $|S| = N$ exakt $1 \le k < N$ Elemente zu ziehen.
\end{lemma}
\begin{proof}
Der Algorithmus baut $R$ inkrementell auf: wir beginnen mit $|R| = 0$ und fügen dann schrittweise je ein Element hinzu.
Sei $X$ eine Zufallsvariable, die beschriebt, wie viele Iterationen der Algorithmus benötigt um $R$ zu erzeugen.
Ferner seien $X_1, \ldots, X_k$ Zufallsvariablen, wobei $X_i$ beschriebt, wie viele Versuche benötigt werden um das $i$-te Elemente zu finden (d.h. um $R$ von $|R| = i - 1$ auf $|R| = i$ zu heben).
Per Konstruktion gilt $X = \sum_{i=1}^k X_i$ und aufgrund der Linearität des Erwartungswertes folgt
\begin{equation}
\expv{X} = \expv{\sum_{i=1}^k X_i} = \sum_{i=1}^k \expv{X_i}. \label{eq:basis-ziehen-ohne-zuruecklegen_expvX}
\end{equation}
Wir können also zunächst die Erwartungswerte $\expv{X_i}$ isoliert betrachten.
Unmittelbar bevor wir das $i$-te Element ziehen, haben wir $i-1$ Elemente bereits bestimmt, d.h. $|R| = i - 1$.
Mit Wahrscheinlichkeit $q_i = (i-1) / N$ erwischen wir also ein bereits gezogenes Element und versuchen es erneut.
Mit Wahrscheinlichkeit $p_i = 1 - q_i = (N - (i - 1)) / N$ ziehen wir aber ein Element, dass noch nicht Teil von $R$ ist.
Es gilt also
\begin{align*}
\prob{X_i = 1} & = p_i \\
\prob{X_i = 2} & = (1 -p_i) p_i \\
\vdots \quad\quad & = \quad\quad \vdots \\
\prob{X_i = 1 + \ell} & = (1- p_i)^{\ell} p_i.
\end{align*}
Die Zufallsvariable $X_i - 1$ ist also geometrisch verteilt mit Erfolgswahrscheinlichkeit~$p_i$.
Es gilt daher $\expv{X_i} = 1 / p_i = N / (N - i + 1)$. Durch Einsetzen in \cref{eq:basis-ziehen-ohne-zuruecklegen_expvX} erhalten wir:
\begin{align}
\expv{X} & = \sum_{i=1}^k \frac{N}{N-(i-1)} = \sum_{i=0}^{k-1} \frac{N}{N - i} \\
& = N \sum_{i=N-k+1}^N \frac{1}{i} \\
& = N \left[ \left( \sum_{i=1}^N \frac{1}{i} \right) - \left( \sum_{i=1}^{N-k} \frac{1}{i} \right) \right] \\
& = N \left[ H_N - H_{N-k} \right]
\end{align}
Im letzten Schritt nutzen wir die Definition $H_t = \sum_{i=1}^t 1/i$, also die $t$-te Partialsumme der Harmonischen Reihe.
Es gilt, dass $H_t = \ln t + \gamma + 1/(2t) + \varepsilon_t$, wobei $\gamma \approx 0.5772$ die Euler–Mascheroni Konstante ist und $0 \le \varepsilon_t \le N^{-2}/8$ ein verschwindender Fehlerterm ist.
Es folgt also:
\begin{align}
\expv{X} & = N \left[ H_N - H_{N-k} \right] \\
& \le N \left[ \ln N - \ln (N - k + 1) \right] + \Oh{1} \\
& = N \ln \frac{N}{N - k + 1} + \Oh{1}
\end{align}
\end{proof}
\noindent
\begin{figure}
\begin{center}
\includegraphics[width=0.5\textwidth]{data/gnm_iterationen.pdf}%
\includegraphics[width=0.5\textwidth]{data/gnm_iterationen_hist.pdf}
\end{center}
\caption{Simulation von \cref{algo:basis-ziehen-ohne-zuruecklegen} (ohne Komplementbildung).
\textbf{Links}: $N=10\,000$ und $2 \le k \le 3N/4$; für jeden Parameter wurden 1\,000 unabhängige Simulationen ausgeführt.
\textbf{Rechts}: $N=1\,000$, $k = N/2$. Es wurden 10\,000 unabhängige Simulationen ausgeführt.
}
\label{fig:ziehen-ohne-zuruecklegen-iterationen}
\end{figure}
Beobachte, dass die Analyse im Beweis von \cref{lemma:basis-ziehen-ohne-zuruecklegen-versuche} sehr eng ist.
Wie wir in \cref{fig:ziehen-ohne-zuruecklegen-iterationen} sehen, gilt für relativ kleine $N$ schon $\expv{X} \approx N \ln \frac{N}{N - k + 1}$.
\begin{exercise}
Zeige, dass $N \ln \frac{N}{N - k + 1} = \Oh{k}$ für $k \le N / 2$.
\end{exercise}
\bigskip
\Cref{lemma:basis-ziehen-ohne-zuruecklegen-versuche} besagt, dass \cref{algo:basis-ziehen-ohne-zuruecklegen} für $k \le N / 2$ erwartet $\Oh{k}$ Iterationen ausführt und somit (in Erwartung) optimal ist.
Nur für $k \to N$ wird der $\log$-Faktor asymptotisch relevant, und wir bekommen eine suboptimale Laufzeit von $\Omega(N \log N)$.
Auf \aside{Komplementbildung stellt sicher: $k \ge N /2$} der einen Seite sollte das für \emph{übliche} Graphen kein Problem sein (sie sind dünn und daher $k \ll N$); dennoch lässt sich das Problem einfach vermeiden.
Für $k > N / 2$ können wir zunächst das Komplement~$\bar R$ mit $| \bar R | = N - k < N/2$ ziehen und dann alle Elemente außer solche in $\bar R$ ausgeben.
Zwar benötigt die Komplementbildung $R = \bar{\bar R}$ während der Ausgabe $\Theta(N)$ Zeit, allerdings hat die Ausgabe ebenfalls Größe $\Omega(N)$, wodurch das asymptotisch optimal bleibt.
\subsection{Konzentration um den Erwartungswert}
Bisher haben wir uns darauf konzentriert Performance \emph{in Erwartung} zu analysieren.
In der Praxis kann das aber manchmal zu wenig sein:
wir stellen uns eine Implementierung vor, die in 99~\% der Fälle in 1\,s durchläuft aber in 1~\% der Ausführungen 101\,s benötigt.
Die durchschnittliche (\qq{erwartete}) Laufzeit wären also 2\,s.
Wenn wir unsere Pipelines darauf auslegen, kann das zu (im schlimmsten Fall kaskadierenden) Timeouts führen.
\Cref{fig:ziehen-ohne-zuruecklegen-iterationen} (rechts) scheint zu zeigen, dass \cref{algo:basis-ziehen-ohne-zuruecklegen} nicht in diese Kategorie fällt.
Unter 10\,000 Ausführung ist die größte Abweichungen vom Mittelwert nur rund 12\,\%.
Um Beobachtungen wie diese zu formalisieren, nutzt man häufig folgende Definition:
\begin{definition}
Eine \aside{mit hoher Wahrscheinlichkeit} Eigenschaft~$X$ gilt \emph{mit hoher Wahrscheinlichkeit} (with high probability, whp),
wenn für einen hinreichend großen Parameter $n$ die Gegenwahrscheinlichkeit $\prob{\lnot X} \le 1/n$ nicht übersteigt.
\end{definition}
In anderen Worten: je \qq{größer} unser Problem wird, desto unwahrscheinlicher sind extreme Ausreißer.
In der Regel ist der Parameter aus dem Kontext klar und wird daher in der Literatur oft nicht erwähnt (im folgenden Lemma wäre es $k$):
\begin{lemma}
\Cref{algo:basis-ziehen-ohne-zuruecklegen} verwendet mit hoher Wahrscheinlichkeit \Oh{k} Iterationen für $k \le N/2$.
\end{lemma}
\begin{proof}
Wir verwenden dieselbe Modellierung wie im Beweis von \cref{lemma:basis-ziehen-ohne-zuruecklegen-versuche}:
Sei $X = \sum_{i=1}^k X_i$ die Gesamtanzahl an Iterationen und $X_i - 1$ geometrisch verteilt mit Erfolgswahrscheinlichkeit~$p_i$.
Da $k \le N/2$, ziehen wir zu jedem Zeitpunkt mindestens mit Wahrscheinlichkeit $1/2$ ein neues Element, d.h. $p_i \ge 1/2$ für alle $i \ge 1$.
Um die harmonische Reihe zu vermeiden, nutzen wir eine schlechtere Abschätzung, bei der wir alle $p_i = 1/2$ auf den minimalen Erfolgswert setzen.
Sei $X' = \sum_{i=1}^k X'_i$ diese Verteilung und $X'_i - 1$ geometrisch mit Erfolgswahrscheinlichkeit $1/2$.
\noindent
Für einen beliebigen Grenzwert~$t$ gilt also
\begin{equation}\label{eq:xprimet_groesser_xt}
\prob{X > t} \le \prob{X' > t}.
\end{equation}
Wir können nun unsere Worst-Case Abschätzung $X'$ auch alternativ interpretieren:
Wir werfen solange Münzen bis wir $k$ mal \qq{Kopf} gesehen haben.
Die Anzahl der Versuche entspricht dann $X'$.
Wir formalisieren diese alternative Idee:
Sei $Y^{(t)} = \sum_{i=1}^t Y_i$ wobei $Y_i$ unabhängige Zufallsvariablen mit $\prob{Y_i = 1} = 1/2$ sind.
Per Definition gilt
\begin{equation}
\prob{X > t}
\quad \stackrel{\text{Eq. \ref{eq:xprimet_groesser_xt}}}{\le} \quad
\prob{X' > t}
\quad = \quad
\prob{Y^{(t)} < k}
\quad \le \quad
\prob{Y^{(t)} \le k}
\end{equation}
\noindent
Außerdem gilt $\expv{Y^{(t)}} = t/2$.
Da $Y_i$ unabhängige Bernoullizufallsvariablen sind, hält die Chernoff-Ungleichung:
\begin{align}
& & \prob{Y^{(t)} \le (1-\delta) \expv{Y_t}} & \le \exp\left(-\delta^2 \expv{Y^{(t)}} / 2\right) \\
& \stackrel{\delta:=1/2}{\Rightarrow} & \prob{Y^{(t)} \le t / 4} & \le \exp\left(-(1/2)^2 (t/2) / 2\right) = \exp\left( -t/16 \right)
\end{align}
\noindent Im letzten Schritt wählten wir $\delta = 1/2$; nun setzen wir $t = 4k$:
\begin{align}
& & \prob{Y^{(4k)} \le k} & \le \exp(-k/4) \le 1/k \quad \text{für } k \ge 9 \\
& \Rightarrow & \prob{X > 4k} & \le 1/k
\end{align}
\noindent
Es folgt $\prob{X \le 4k} \ge 1 - 1/k$. Somit gilt $X < 4k$ mit hoher Wahrscheinlichkeit.
\end{proof}
Welche Datenstruktur benötigen wir für $R$? Sie sollte Existenzanfragen ($x \in R$?) und Einfügen ($R \gets R \cup \set{x}$) effizient durchführen können.
Ein Hashset unterstützt beide Operationen in erwartet konstanter Zeit, wodurch sich eine erwartete Laufzeit von $\Oh{k}$ ergibt; dies ist optimal in Erwartung.
Für praktische Implementierungen kann es sich jedoch auch lohnen, das Hashset durch ein Bitset (d.h. ein Array in dem jeder Eintrag genau ein Bit groß ist) zu ersetzen.
Auf dem Papier benötigt ein Bitset $\Theta(N)$ viele Bits und ist somit in der Initialisierung teuer.
In der Praxis benötigt aber jeder Eintrag in einem Hashset mindestens 64 Bit, sodass ein Bitset bereits für $N / k > 64$ weniger Platz benötigt (das trifft etwa auf 30\,\% der Netzwerke in \cref{fig:kantenanzahl} zu).
Zusätzlich sind Existenzanfragen und Einfügeoperationen extrem schnell, wodurch ein Bitset oft einen lohnenswerten Kompromiss zwischen mehr Speicher für schnellere Ausführung darstellen kann.
\begin{widefigure}
\includegraphics[width=0.32\textwidth]{data/gnm_scale0.pdf}\hfill
\includegraphics[width=0.32\textwidth]{data/gnm_scale1.pdf}\hfill
\includegraphics[width=0.32\textwidth]{data/gnm_scale2.pdf}
\caption{
Laufzeit~$T$ pro Sample~$k$ für das Ziehen von $k$ Elementen aus $S = \set{1, \ldots, N}$ als Funktion von $|S|$.\\
\textbf{Links:} $k=10$, \textbf{Mitte: } $k = \sqrt{N}$, \textbf{Rechts: } $k = N / 4$.
}
\label{fig:benchmark_gnm_scale}
\end{widefigure}
In \cref{fig:benchmark_gnm_scale} stellen wir Laufzeitmessungen für folgendes Experiment dar.
Wir variieren die Größe unseres Universums $|S| = N$ von $2^8$ bis $2^{32}$ und betrachten drei Szenarien:
\begin{itemize}
\item $k=10$: Trotz steigender Universumsgröße ziehen wir immer nur 10 Samples.
\item $k=\sqrt{N}$: Die Stichprobengröße steigt langsam in der Universumsgröße (dieser Fall entspricht einem $\Gnm$ Graphen mit $m = \Theta(n)$).
\item $k=N/4$: Jedes vierte Element des Universums wird ausgewählt.
\end{itemize}
\begin{exercise}
Welche Effekte kannst du in \cref{fig:benchmark_gnm_scale} beobachten? Wie erklären sich diese?
\end{exercise}
Wir möchten uns aber nur auf einen zentralen Effekt in den Abbildungen konzentrieren, da dieser bei großen zufälligen Daten oft auftritt.\footnote{foreshadowing\ldots}
In \cref{fig:benchmark_gnm_scale} (rechts) beobachteten wir einen sprunghaften Anstieg der Zeit pro Element --- teils um fast einen Faktor 10.
Unsere Analysen sagen jedoch eine konstante Laufzeit pro gezogenem Element für $k = N/4$ voraus; einen sprunghaften Anstieg können sie nicht erklären.
Wir sind also im Algorithm-Engineering-Zyklus auf eine Diskrepanz zwischen Theorie und Praxis gestoßen!
Ist unsere Analyse falsch? Jaein.
In der Bestimmung der Laufzeit nutzen wir eine Unit-Cost Random-Access-Machine (d.h. wir gingen davon aus, dass jede Operation gleich viel Zeit benötigt).
Moderne Computer verfügen jedoch über ausgeklügelte Techniken, die gewisse Aspekte besonders beschleunigen.
Hierzu gehören die sog. Speicherhierarchien: Ein Prozessor kann nur auf Daten operieren, die sich auf dem physischen Chip befinden, genauer in den Registern.
Aus physikalischen und ökonomischen Gründen, gibt es aber nur sehr wenig \qq{Register-Speicher}.
Die meisten aktiven Daten werden daher im Arbeitsspeicher vorgehalten; dieser ist aber recht langsam.
Eine typische CPU kann in der Zeit, die es dauert ein Datum aus dem Arbeitsspeicher zu beziehen, etwa 100 bis 1000 Operationen ausführen.
Um dies aufzufangen, werden mehrere Caches eingefügt: die CPU \qq{rät} welche Daten in naher Zukunft benötigt werden und hält diese im schnelleren Cache vor.
Dies klappt aber mit Zufallsdaten oft nicht (sie sollen ja gerade \emph{nicht} vorhersehbar sein).
In \cref{fig:benchmark_gnm_scale} (rechts) sehen wir also genau den Punkt ab dem die Datenstrukturen nicht mehr in den Cache passen.
Für Hashset und Reservoir-Sampling nutzen diese Datenstrukturen etwa $64 k = 16 N$ Bits (das Hashset aufgrund des Load-Factors < 1 auch ein bisschen mehr).
Das Bitset benötigt $N$ Bits --- daher kommt hier der Performanceeinbruch auch später.
Um über diese qualitative Erklärung hinaus, das Verhalten analytisch zu erklären, benötigen wir also ein Maschinenmodell, das auch Speicherzugriffe berücksichtigt.
Dem widmen wir uns aber erst später.
Im Folgenden betrachten wir stattdessen eine Parallelisierung des Algorithmus --- dieselbe Strategie führt jedoch auch zu einem Generator, der zufällige Speicherzugriffe auf den Cache beschränken kann.
Wir brauchen zunächst ein Maschinenmodell für parallele Berechnungen.
\subsection{Fork-Join Parallelismus}
Das Fork-Join-Modell ist eine Erweiterung der Random-Access-Machine, welche die RAM um die zwei namesgebenenden Instruktionen erweitert.
Eine Berechnung startet als gewöhnliches sequentielles Programm, einem sog. Task.
Angenommen wir führen gerade $t_0$ aus, dann können wir die Berechnung in zwei \glqq Kindertasks\grqq{} $t_1$ und $t_2$ teilen \aside{fork} (\emph{forken}).
Beide Tasks werden unabhängig von einander ausgeführt, während $t_0$ ruht.
Sobald $t_1$ und $t_2$ fertig sind, kommt es zur Wiedervereinigung \aside{join} (\emph{join}):
$t_1$ und $t_2$ werden \glqq gelöscht \grqq, $t_0$ wieder geweckt und setzt seine Berechnung fort.
Betrachten wir das parallele Summieren $\sum_{i=1}^n x_n$ von Zahlen $X = (x_1, \ldots, x_n)$ als einfaches Beispiel:
\begin{algorithm}
\SetKwProg{Fn}{Function}{}{end}
\SetKwFunction{ParSum}{ParSum}
\SetKwFunction{Fork}{Fork}
\SetKwFunction{Join}{Join}
\Fn{\ParSum{$X=(x_1, \ldots, x_n)$}}{
\If{$n = 1$}{
gebe $x_1$ zurück und beende Task\;
}
Berechne Mitte~$m \gets \lfloor n/2\rfloor$\;
\Fork{$
\underbrace{\ParSum{$X_L = (x_1, \ldots, x_m)$}}_{\text{Task $t_1$}}, \ \
\underbrace{\ParSum{$X_R = (x_{m + 1}, \ldots, x_n)$}}_{\text{Task $t_2$}}
$}\;
$(s_1, s_2) \gets \Join{$t_1, t_2$}$\;
gebe $s_1 + s_2$ zurück und beende Task
}
\caption{Parallele Summe im Fork-Join Modell}
\label{algo:parallel_sum}
\end{algorithm}
In \cref{algo:parallel_sum} teilen wir die Eingabe rekursiv in zwei gleich große Teile (modulo Rundung) und berechnen die Teilsummen für beide rekursiv.
Sobald die Teilergebnisse vorliegen, summieren wir diese in $\Oh{1}$ Zeit und geben es als Gesamtergebnis zurück.
Dieses Schema funktioniert analog für alle assoziativen Operationen (z.B. Produkt, Min, Max).
In der Praxis findet man einige Software-Bibliotheken, die Fork-Join-Parallelismus zur Verfügung stellen;
für C++ etwa \texttt{Cilk} oder Intels \texttt{oneTBB}, für Rust \texttt{rayon}.
Diese praktischen Lösungen nutzen einen sog. Work-Stealing-Scheduler.
Grob \aside{Work Stealing} kann man annehmen, dass der Prozessor~$P_0$, der $t_0$ ausführte den Task $t_2$ als verfügbar bekannt gibt und dann selbst $t_1$ ausführt.
Wenn während der Ausführung von $t_1$ ein anderer Prozessor verfügbar ist, übernimmt dieser $t_2$.
In diesem Fall sind $t_1$ und $t_2$ also nebenläufig.
Wenn sich kein anderer Prozessor findet, wird $P_0$ den Task $t_2$ übernehmen, sobald $t_1$ beendet wurde.
Wir erhalten also ein dynamisches System, das (i) nicht wissen muss, wie viele Prozessoren es gibt, und (ii) auch relativ gut damit umgehen kann, dass $t_1$ und $t_2$ evtl. sehr unterschiedliche Laufzeiten haben.
Praktisch ist es kein Problem viel mehr Tasks zu erzeugen als es Prozessoren gibt.
Wir messen die Güte eines Fork-Join-Algorithmus mit zwei Eigenschaften:
\begin{enumerate}
\item
Die \aside{Arbeit / work} \emph{Arbeit} (engl. work) ist die Summe aller Instruktionen, die das Programm ausführt.
Dies entspricht also der Laufzeit, wenn nur ein Prozessor zur Verfügung steht.
Für \cref{algo:parallel_sum} ergibt sich die Arbeit als
\begin{equation}
W(n) = \begin{cases}
\Oh{1} & \text{falls } n = 1 \\
2 \cdot W(n/2) + \Oh{1} & \text{falls } n > 1
\end{cases}
= \Oh{n}.
\end{equation}
Der \aside{arbeits-optimal} Algorithmus leistet also asymptotisch nicht mehr Arbeit als die Laufzeit der besten sequentiellen Lösung; wir sagen er ist \emph{arbeits-optimal}.
\item
Der \aside{Span / Parallel Depth} Span (engl. span oder parallel depth) misst die Anzahl der Instruktionen auf einem kritischen Pfad (d.h. eine längste Folge abhängiger Schritte).
Dies entspricht der Laufzeit, wenn wir unbeschränkt viele Prozessoren zur Verfügung haben.
In \cref{algo:parallel_sum} müssen wir also analysieren, wie viele Instruktionen vor und nach dem Fork/Join Paar erforderlich sind und wie tief die Rekursion läuft:
\begin{equation}
S(n) = \underbrace{\Oh{1}}_\text{Arbeit vor/nach Fork/Join} + \underbrace{S(n/2)}_\text{Rekursion} = \Oh{\log n}.
\end{equation}
\end{enumerate}
Wenn wir ausreichend Prozessoren zur Verfügung haben, können wir also $n$ Zahlen in Zeit $\Oh{\log n}$ summieren.
Zu ähnlichen Ergebnissen kommt man auch in anderen parallelen Maschinenmodellen (z.B. P-Rams).
\subsection{Rekursives Ziehen}
In \cref{algo:parallel_sum} haben wir bereits ein gutes Grundgerüst für parallele Algorithmen gesehen:
wenn wir es schaffen das Problem \emph{schnell} in unabhängige Hälften zu teilen und dann rekursiv zu bearbeiten, sind wir auf einem guten Weg.
Bezogen auf unser Problem, $k$~Stichproben aus $N$~Elementen zu ziehen, haben wir zwei Dimensionen:
\begin{enumerate}
\item Teile $S$ deterministisch in $S_1$ und $S_2$: Wir können das Universum~$S$ mit $|S| = N$ in zwei Hälften $S_1$ und $S_2$ partitionieren mit $|S_1| = \lfloor N / 2 \rfloor$.
Dann stellt sich die Frage: was ist die Anzahl~$k_1$ der Samples aus $S_1$, bzw. $k_2$ aus $S_2$. Klar ist nur $k_1 + k_2 = k$.
\item Teile $k$ deterministisch in $k_1$ und $k_2$: Wir können die Stichproben teilen, d.h. $k_1 = \lfloor k/n \rfloor$ wählen; dann müssen wir jedoch eine unbekannte Partitionierung von $S$ finden.
\end{enumerate}
Auf den ersten Blick wirkt die zweite Variante zielführender.
Es ist vorteilhaft, wenn beide Teilprobleme etwa die gleiche Stichprobenanzahl haben, da die Arbeit in dieser Anzahl skalieren soll.
Tatsächlich liefert aber die erste Variante einfachere Algorithmen, die auch keine signifikanten Probleme mit unbalancierten Teilproblemen haben.
\begin{figure}[t]
\begin{center}
\begin{tikzpicture}
\foreach \x in {-8, -4, -2, 3, 5, 7, 9} {
\node[inner sep=0, minimum size=1em, fill=black!20, anchor=west] at (\x em, 0) {};
}
\node[draw, minimum height=1em, minimum width=20em] at (0,0) {};
\node[anchor=north] at (-9.5em, -0.5em) {1};
\node[anchor=north] at ( 9.5em, -0.5em) {$N$};
\foreach \x in {-9,...,9} {
\path[draw, black!50] (\x em, 0.5em) to ++(0, -1em);
}
\path[draw, thick, dashed] (0, 2em) to (0, -3.5em);
\draw[decorate, decoration = {brace}] (-9.75em,3em) to node[above] {$k$} ++(19.5em,0);
\draw[red, decorate, decoration = {brace}] (-9.5em,1em) to node[above] {$k_1$} ++(9em,0);
\draw[blue, decorate, decoration = {brace}] ( 0.5em,1em) to node[above] {$k_2$} ++(9em,0);
\node[anchor=east] at (-12em, 1.8em) {zufällig:};
\draw[decorate, decoration = {brace, mirror}] (-9.75em,-4em) to node[below] {$N$} ++(19.5em,0);
\draw[red, decorate, decoration = {brace, mirror}] (-9.5em,-2em) to node[below] {$N / 2$} ++(9em,0);
\draw[blue, decorate, decoration = {brace, mirror}] ( 0.5em,-2em) to node[below] {$N / 2$} ++(9em,0);
\node[anchor=east] at (-12em, -2.8em) {deterministisch:};
\end{tikzpicture}
\end{center}
\caption{
Ziehen von $k=7$ Stichproben aus $N=20$ Elementen ohne Zurücklegen; die grauen Elemente seien ausgewählt.
Die linke Hälfte des Universums liefert $k_1$ Stichproben, die rechte $k_2$.
Im gezeigten Beispiel gilt $k_1 = 3$ und $k_2 = 4$.
}
\label{fig:rekursive_k_aus_N}
\end{figure}
Wie in \cref{fig:rekursive_k_aus_N}, teilen wir also das Universum~$S$ in zwei gleich große Hälften $S_1$ und $S_2$.
Angenommen, wir täten dies nicht, wie viele Stichproben $k_1$ würde ein anderer Algorithmus aus $S_1$ ziehen?
Wir wissen es nicht genau, da $k_1$ eine Zufallsvariable ist --- sie ist Hypergeometrisch verteilt:
aus einer Gesamtpopulation von $N$ mit $k$~Treffern (d.h. die Elemente, die ein anderer Algorithmus gezogen hätte), ziehen wir $N_1$ Stichproben und erhalten $k_1$~Treffer.
In Erwartung sind es also $\expv{k_1} = k |S_1| / |S|$.
Analog zu den zufälligen Sprüngen in \cref{subsec:gnp_zufaellige_spruenge}, ziehen wir $k_1$ aus einer Hypergeometrisch Verteilung.
Dann setzen wir den Prozess bedingt auf $k_1$ fort:
\begin{algorithm}
\SetKwProg{Fn}{Function}{}{end}
\SetKwFunction{ParSelect}{ParSample}
\SetKwFunction{Fork}{Fork}
\SetKwFunction{Join}{Join}
\Fn{\ParSelect{$S=\set{s_1, \ldots, s_N}, k$}}{
\uIf{$k = 0$}{
gebe $\emptyset$ zurück und beende Task\;
} \uElseIf{$k = 1$}{
$x \gets$ zufällig uniform aus $S$\;
gebe $\set{x}$ zurück und beende Task\;
} \ElseIf{$k = N$}{
gebe $S$ zurück und beende Task\;
}
\BlankLine
$N_L \gets \lfloor N / 2 \rfloor$\;
$k_L \gets$ zufällig hypergeometrisch: ziehe $N_L$ Element aus $N$ mit $k$ Treffern\;
\BlankLine
\Fork{$
\underbrace{\ParSelect{$S_L = \set{s_1, \ldots, s_{N_L}}, k_L$}}_{\text{Task $t_1$}}, \ \
\underbrace{\ParSelect{$S_R = \set{s_{N_L + 1}, \ldots, s_N}, k - k_L$}}_{\text{Task $t_2$}}
$}\;
$(R_1, R_2) \gets \Join{$t_1, t_2$}$\;
gebe $R_1 \cup R_2$ zurück und beende Task
}
\caption{Paralleles Ziehen von $k$ Stichproben aus $S$ ohne Zurücklegen}
\label{algo:parallel_k_aus_N}
\end{algorithm}