2012-06-12 30 views
9

Spero di creare 3 numeri (quasi non negativi) quasi casuali che sommati a uno, e ripetere più e più volte.Genera 3 numeri casuali che sommano a 1 in R

Fondamentalmente sto provando a suddividere qualcosa in tre parti casuali su molte prove.

Mentre io sono a conoscenza di

a = runif (3,0,1)

stavo pensando che potrei usare 1-a come il massimo nel prossimo periodo se, ma sembra disordinato .

Ma questi naturalmente non si sommano a uno. Qualche idea, oh saggio stackoverflow?

+2

È un'opzione per rinormalizzare i numeri casuali dopo la generazione? –

+0

Che ne dici di generare 2 numeri casuali aeb? Quindi a + b + c = 1 => c = 1 - (a + b) –

+0

e se a e b somma a maggiore di 1? – mmann1123

risposta

9

solo casuali 2 cifre da (0, 1) e se assumere la sua a e b poi avete ottenuto:

rand1 = min(a, b) 
rand2 = abs(a - b) 
rand3 = 1 - max(a, b) 
+0

Inoltre, devi ripetere la generazione del secondo numero se a == b ... (dovrebbe essere MOLTO caso raro) – ddzialak

+0

@user a = 0,85 , b = 0.99 quindi hai i numeri: 0,85, 0,14, 0,01 (per quanto mi riguarda questi sono 3 numeri casuali molto buoni da 0..1) – ddzialak

+3

La distribuzione risultante sembra non essere esattamente banale: http: //www.jstor. org/discover/10.2307/2983572? uid = 2129 & uid = 2 & uid = 70 & uid = 4 & sid = 21100849643501 e un documento successivo a cui è possibile accedere liberamente http://doc.utwente.nl/70657/1/Sleutel67random.pdf – Christian

4

Immagino che dipende da quale distribuzione si vuole sui numeri, ma qui è un modo :

diff(c(0, sort(runif(2)), 1)) 

Usa replicate per ottenere il maggior numero di set che si desidera:

> x <- replicate(5, diff(c(0, sort(runif(2)), 1))) 
> x 
      [,1]  [,2]  [,3]  [,4]  [,5] 
[1,] 0.66855903 0.01338052 0.3722026 0.4299087 0.67537181 
[2,] 0.32130979 0.69666871 0.2670380 0.3359640 0.25860581 
[3,] 0.01013117 0.28995078 0.3607594 0.2341273 0.06602238 
> colSums(x) 
[1] 1 1 1 1 1 
11

Questa domanda riguarda questioni più sottili di quanto potrebbe essere a prima vista evidente. Dopo aver guardato il seguente, si consiglia di riflettere attentamente sul processo che si sta utilizzando questi numeri per rappresentare:

## My initial idea (and commenter Anders Gustafsson's): 
## Sample 3 random numbers from [0,1], sum them, and normalize 
jobFun <- function(n) { 
    m <- matrix(runif(3*n,0,1), ncol=3) 
    m<- sweep(m, 1, rowSums(m), FUN="/") 
    m 
} 

## Andrie's solution. Sample 1 number from [0,1], then break upper 
## interval in two. (aka "Broken stick" distribution). 
andFun <- function(n){ 
    x1 <- runif(n) 
    x2 <- runif(n)*(1-x1) 
    matrix(c(x1, x2, 1-(x1+x2)), ncol=3) 
} 

## ddzialak's solution (vectorized by me) 
ddzFun <- function(n) { 
    a <- runif(n, 0, 1) 
    b <- runif(n, 0, 1) 
    rand1 = pmin(a, b) 
    rand2 = abs(a - b) 
    rand3 = 1 - pmax(a, b) 
    cbind(rand1, rand2, rand3) 
} 

## Simulate 10k triplets using each of the functions above 
JOB <- jobFun(10000) 
AND <- andFun(10000) 
DDZ <- ddzFun(10000) 

## Plot the distributions of values 
par(mfcol=c(2,2)) 
hist(JOB, main="JOB") 
hist(AND, main="AND") 
hist(DDZ, main="DDZ") 

enter image description here

+0

Nice, I was pensando di tracciare i risultati ma lo hai già fatto. È interessante vedere che apparentemente nessuna delle soluzioni fa realmente ciò che si vorrebbe intuitivamente. È anche interessante il fatto che in questi diagrammi non si possa davvero vedere che DDZ fa la cosa giusta secondo i mezzi mentre AND non lo fa nemmeno. – Christian

6

Quando si desidera generare casualmente i numeri che si aggiungono a 1 (o qualche altro valore) quindi si dovrebbe guardare il Dirichlet Distribution.

c'è una funzione rdirichlet nel pacchetto gtools e correre RSiteSearch('Dirichlet') porta in primo piano alcuni colpi bel che potrebbe facilmente portare a strumenti per fare questo (e non è difficile da codice a mano sia per semplici distribuzioni di Dirichlet).

2

Questo problema e le diverse soluzioni proposte mi hanno incuriosito. Ho fatto un piccolo test dei tre algoritmi di base suggeriti e quali valori medi avrebbero prodotto per i numeri generati.

choose_one_and_divide_rest 
means:    [ 0.49999212 0.24982403 0.25018384] 
standard deviations: [ 0.28849948 0.22032758 0.22049302] 
time needed to fill array of size 1000000 was 26.874945879 seconds 

choose_two_points_and_use_intervals 
means:    [ 0.33301421 0.33392816 0.33305763] 
standard deviations: [ 0.23565652 0.23579615 0.23554689] 
time needed to fill array of size 1000000 was 28.8600130081 seconds 

choose_three_and_normalize 
means:    [ 0.33334531 0.33336692 0.33328777] 
standard deviations: [ 0.17964206 0.17974085 0.17968462] 
time needed to fill array of size 1000000 was 27.4301018715 seconds 

Il tempo misure sono da prendere con un grano di sale in quanto potrebbero essere più influenzati dalla gestione della memoria Python che con lo stesso algoritmo. Sono troppo pigro per farlo correttamente con timeit. L'ho fatto su Atom a 1 GHz in modo che questo spieghi perché ci è voluto così tanto tempo.

In ogni caso, scegliere_uno_e_divide_rest è l'algoritmo suggerito da Andrie e il poster della domanda stessa (AND): si sceglie un valore a in [0,1], quindi uno in [a, 1] e poi si guarda quello che hai lasciato. Si aggiunge a uno, ma questo è tutto, la prima divisione è due volte più grande degli altri due. Uno potrebbe averlo indovinato ...

choose_two_points_and_use_intervals è la risposta accettata da ddzialak (DDZ). Prende due punti nell'intervallo [0,1] e usa la dimensione dei tre sub-intervalli creati da questi punti come i tre numeri. Funziona come un fascino e i mezzi sono tutti 1/3.

choose_three_and_normalize è la soluzione di Anders Gustafsson e Josh O'Brien (JOB). Genera solo tre numeri in [0,1] e li normalizza di nuovo a una somma di 1. Funziona altrettanto bene e sorprendentemente un po 'più veloce nella mia implementazione Python. La varianza è leggermente inferiore rispetto alla seconda soluzione.

Là ce l'hai. Non ho idea di quale sia la distribuzione beta di queste soluzioni o quale insieme di parametri nel documento corrispondente a cui ho fatto riferimento in un commento, ma forse qualcun altro può capirlo.

Problemi correlati