Wie können sortierte, gleichmäßig verteilte Werte in einem Intervall effizient generiert werden?

12

Angenommen, ich möchte aus dem Intervall eine Reihe von Zufallszahlen generieren (a, b). Die generierte Sequenz sollte auch die Eigenschaft haben, dass sie sortiert ist. Ich kann mir zwei Möglichkeiten vorstellen, um dies zu erreichen.

Sei ndie Länge der zu erzeugenden Sequenz.

1. Algorithmus:

Let `offset = floor((b - a) / n)`
for i = 1 up to n:
   generate a random number r_i from (a, a+offset)
   a = a + offset
   add r_i to the sequence r

2. Algorithmus:

for i = 1 up to n:
    generate a random number s_i from (a, b)
    add s_i to the sequence s
sort(r)

Meine Frage ist, erzeugt Algorithmus 1 Sequenzen, die so gut sind wie die von Algorithmus 2 erzeugten?

Ultrajohn
quelle
Übrigens ist es bemerkenswert einfach, eine Liste sortierter Zufallszahlen in zu erstellen R. Um ein Array von Mengen von n Zufallszahlen über ein einheitliches Intervall [ a , b ] zu erzeugen, funktioniert der folgende Code : . kn[a,b]rand_array <- replicate(k, sort(runif(n, a, b))
RobertF

Antworten:

18

Der erste Algorithmus schlägt aus zwei Gründen schlecht fehl :

  1. Wenn Sie den Boden von kann dies drastisch reduziert werden. In der Tat, wenn b - a < n ist , ist es Null, was Ihnen eine Menge gibt, deren Werte alle gleich sind!(ein- -b)/.nb- -ein<n

  2. Wenn Sie nicht das Wort ergreifen, sind die resultierenden Werte zu gleichmäßig verteilt. Zum Beispiel gibt es in jeder einfachen Zufallsstichprobe von iid einheitlichen Variablen (etwa zwischen a = 0 und b = 1 ) eine ( 1 - 1 / n ) n1 / e 37 % Wahrscheinlichkeit, dass die größte nicht sein wird im oberen Intervall von 1 - 1 / n bis 1 . Mit Algorithmus 1 gibt es eine 100nein=0b=1(1- -1/.n)n1/.e37%.1- -1/.n1100%.Chance, dass das Maximum in diesem Intervall liegt. Für einige Zwecke ist diese Supergleichmäßigkeit gut, aber im Allgemeinen ist es ein schrecklicher Fehler, weil (a) viele Statistiken ruiniert werden, aber (b) es sehr schwierig sein kann, festzustellen, warum.

  3. Wenn Sie das Sortieren vermeiden möchten, generieren Sie stattdessen unabhängige exponentiell verteilte Variablen. Normalisieren Sie ihre kumulative Summe auf den Bereich ( 0 , 1 ), indem Sie durch die Summe dividieren. Löschen Sie den größten Wert (der immer 1 sein wird ). Skalieren Sie auf den Bereich ( a , b ) .n+1(0,1)1(ein,b)

Histogramme aller drei Algorithmen werden angezeigt. (Jedes zeigt die kumulativen Ergebnisse von unabhängigen Sätzen mit jeweils n = 100 Werten.) Das Fehlen einer sichtbaren Variation im Histogramm für Algorithmus 1 zeigt das Problem dort. Die Variation der beiden anderen Algorithmen ist genau das, was zu erwarten ist - und was Sie1000n=100 von einem Zufallszahlengenerator benötigen .

Weitere (amüsante) Möglichkeiten zum Simulieren unabhängiger gleichmäßiger Variablen finden Sie unter Simulieren von Zeichnungen aus einer gleichmäßigen Verteilung mithilfe von Zeichnungen aus einer Normalverteilung .

Abbildung: Histogramme

Hier ist der RCode, der die Figur erzeugt hat.

b <- 1
a <- 0
n <- 100
n.iter <- 1e3

offset <- (b-a)/n
as <- seq(a, by=offset, length.out=n)
sim.1 <- matrix(runif(n.iter*n, as, as+offset), nrow=n)
sim.2 <- apply(matrix(runif(n.iter*n, a, b), nrow=n), 2, sort)
sim.3 <- apply(matrix(rexp(n.iter*(n+1)), nrow=n+1), 2, function(x) {
  a + (b-a) * cumsum(x)[-(n+1)] / sum(x)
})

par(mfrow=c(1,3))
hist(sim.1, main="Algorithm 1")
hist(sim.2, main="Algorithm 2")
hist(sim.3, main="Exponential")
whuber
quelle
Was halten Sie von dem Algorithmus (basierend auf Rangordnungsstatistiken) in meiner Antwort? ;-)
Hat aufgehört - Anony-Mousse
@Anony Es ist eine weniger effiziente Version meines Algorithmus 3. (Ihre scheint eine Menge unnötiger Neuskalierungen zu beinhalten.) Sie generieren die exponentiellen Variablen, indem Sie Protokolle von Uniformen erstellen, was Standard ist.
whuber
6

Der erste Algorithmus erzeugt zu gleichmäßig verteilte Zahlen

Siehe auch Reihen mit geringer Diskrepanz .

[0;;1]]

(Wie bereits ausgeführt, ist dies eine gewünschte Eigenschaft zB für Schichtung sein kann. Low-Diskrepanz Serien wie Halton und Sobel haben ihre Fälle verwenden.)

Ein richtiger, aber teurer Ansatz (für echte Werte)

... soll Beta-verteilte Zufallszahlen verwenden. Die Rangordnungsstatistik der Gleichverteilung ist Beta-verteilt. Sie können dies verwenden, um zufällig die kleinste , dann die zweitkleinste, ... Wiederholung zu zeichnen .

[0;;1]]Beta[1,n]]n1- -X.Beta[n,1]]- -ln(1- -X.)Exponentiell[n]]- -ln(U.[0;;1]])n

- -ln(1- -x)=- -ln(1- -u)n1- -x=u1nx=1- -u1n

Was den folgenden Algorithmus ergibt:

x = a
for i in range(n, 0, -1):
    x += (b-x) * (1 - pow(rand(), 1. / i))
    result.append(x) 

Es kann numerische Instabilitäten geben, und das Berechnen powund Teilen für jedes Objekt kann sich als langsamer als das Sortieren herausstellen.

Für ganzzahlige Werte müssen Sie möglicherweise eine andere Verteilung verwenden.

Das Sortieren ist unglaublich billig, verwenden Sie es also einfach

Ö(nLogn)

Hat aufgehört - Anony-Mousse
quelle
1
Es kann Gründe geben, das Sortieren zu vermeiden. Eine ist, wenn Sie eine große Anzahl von Zufallsvariablen generieren möchten, so viele, dass eine Standardsortierroutine diese nicht verarbeiten kann.
whuber
Ich denke, die numerischen Probleme mit Summen unter Verwendung von Gleitkomma-Mathematik werden viel früher zu einem Problem. (Und die Probleme mit zyklischen Mustern in Pseudozufallszahlen!) Es ist ziemlich einfach, den Sortieransatz auf Terabyte und auf verteilte Systeme auf Exabyte zu skalieren.
Hat aufgehört - Anony-Mousse
1012
Ok, sie nicht speichern zu müssen, ist ein Argument. Aber dann brauchen Sie meinen Ansatz, Ihre Variante 3 mit der kumulierten Summe wird nicht funktionieren.
Hat aufgehört - Anony-Mousse
Das ist ein ausgezeichneter Punkt. Jetzt sehe ich die Tugend der zusätzlichen Berechnungen! (+1)
whuber
5

Es hängt auch davon ab, was Sie mit den Zufallszahlen machen. Bei numerischen Integrationsproblemen würde Methode 1 (wenn sie durch Entfernen des Bodenoperators korrigiert wird) eine überlegene Punktmenge erzeugen. Was Sie tun, ist eine Form der geschichteten Probenahme und hat den Vorteil, dass Verklumpungen vermieden werden. Es ist beispielsweise unmöglich, alle Ihre Werte im Bereich von 0- (ba) / n zu erhalten. Für andere Anwendungen kann dies jedoch sehr schlecht sein. Dies hängt davon ab, was Sie damit tun möchten.

user67054
quelle
2
+1 Ich denke, dies ist ein nützlicher Beitrag zu dieser Frage, insbesondere durch die Charakterisierung von Algorithmus 1 im Hinblick auf die Schichtung.
whuber