Come posso generare punti che corrispondono a un istogramma?

StackOverflow https://stackoverflow.com/questions/423006

Domanda

Sto lavorando a un sistema di simulazione. Presto avrò dati sperimentali (istogrammi) per la distribuzione dei valori nel mondo reale per diversi input di simulazione.

Quando viene eseguita la simulazione, vorrei essere in grado di produrre valori casuali che corrispondano alla distribuzione misurata. Preferirei farlo senza memorizzare gli istogrammi originali. Quali sono alcuni buoni modi di

  1. Mappare un istogramma su un insieme di parametri che rappresentano la distribuzione?
  2. Generazione di valori basati su tali parametri in fase di esecuzione?

EDIT: i dati di input sono durate degli eventi per diversi tipi di eventi. Mi aspetto che tipi diversi abbiano funzioni di distribuzione diverse.

È stato utile?

Soluzione

Almeno due opzioni:

  1. Integra l'istogramma e inverti numericamente.
  2. Rifiuto

Integrazione numerica

Da Calcolo della fisica moderna di William R. Gibbs:

  

Si può sempre integrare numericamente [la funzione] e invertire il [ cdf ]   ma questo spesso non è molto soddisfacente, specialmente se il pdf sta cambiando   rapidamente.

Si crea letteralmente una tabella che traduce l'intervallo [0-1) in intervalli appropriati nella distribuzione di destinazione. Quindi lancia il tuo solito PRNG (alta qualità) e traduci con il tavolo. È ingombrante, ma chiaro, praticabile e completamente generale.

Rifiuto:

Normalizza l'istogramma di destinazione, quindi

  1. Lancia i dadi per scegliere una posizione ( x ) lungo il raggio in modo casuale.
  2. Lancia di nuovo e seleziona questo punto se il nuovo numero casuale è inferiore all'istogramma normalizzato in questo cestino. Altrimenti vai a (1).

Ancora una volta, semplice ma chiaro e funzionante. Può essere lento per la distribuzione con molta probabilità molto bassa (picchi con code lunghe).

Con entrambi questi metodi, puoi approssimare i dati con adattamenti polinomiali a tratti o spline per generare una curva uniforme se non si desidera un istogramma a funzione di passaggio --- ma lasciarlo per dopo come potrebbe essere un'ottimizzazione prematura.

Potrebbero esistere metodi migliori per casi speciali.

Tutto questo è piuttosto standard e dovrebbe apparire in qualsiasi libro di testo di Analisi numerica se ho bisogno di maggiori dettagli.

Altri suggerimenti

Ulteriori informazioni sul problema sarebbero utili. Ad esempio, che tipo di valori sono gli istogrammi? Sono categoriali (ad es. Colori, lettere) o continui (ad es. Altezze, tempo)?

Se gli istogrammi si trovano su dati categorici, penso che potrebbe essere difficile parametrizzare le distribuzioni a meno che non ci siano molte correlazioni tra le categorie.

Se gli istogrammi si trovano su dati continui, potresti provare ad adattare la distribuzione usando miscele di gaussiani. Ossia, prova ad adattare l'istogramma usando $ \ sum_ {i = 1} ^ n w_i N (m_i, v_i) $ dove m_i e v_i sono la media e la varianza. Quindi, quando si desidera generare dati, prima si campiona un i da 1..n con probabilità proporzionale ai pesi w_i e quindi si campiona un x ~ n (m_i, v_i) come si farebbe da qualsiasi gaussiano.

Ad ogni modo, potresti voler leggere di più sui modelli di miscele .

Quindi sembra che ciò che voglio per generare una data distribuzione di probabilità sia una Funzione quantile , che è l'inverso di funzione di distribuzione cumulativa , come dice @dmckee.

La domanda diventa: qual è il modo migliore per generare e memorizzare una funzione quantile che descrive un dato istogramma continuo? Ho la sensazione che la risposta dipenderà molto dalla forma dell'input: se segue qualsiasi tipo di modello ci dovrebbero essere semplificazioni nel caso più generale. Aggiornerò qui mentre vado.

Modifica:

Questa settimana ho avuto una conversazione che mi ha ricordato questo problema. Se rinuncio a descrivere l'istogramma come un'equazione e memorizzo semplicemente la tabella, posso fare selezioni in O (1) tempo? Si scopre che è possibile, senza alcuna perdita di precisione, a scapito dei tempi di costruzione di O (N lgN).

Crea una matrice di N elementi. Una selezione casuale uniforme nell'array troverà un oggetto con probabilità 1 / N. Per ogni oggetto, memorizza la frazione dei risultati per i quali questo oggetto dovrebbe essere effettivamente selezionato e l'indice di un altro oggetto che verrà selezionato in caso contrario.

Campionamento casuale ponderato, implementazione C: 

//data structure
typedef struct wrs_data {
  double share; 
  int pair;
  int idx;
} wrs_t;


//sort helper
int wrs_sharecmp(const void* a, const void* b) {
  double delta = ((wrs_t*)a)->share - ((wrs_t*)b)->share;
  return (delta<0) ? -1 : (delta>0);
}


//Initialize the data structure
wrs_t* wrs_create(int* weights, size_t N) {
  wrs_t* data = malloc(sizeof(wrs_t));
  double sum = 0;
  int i;
  for (i=0;i<N;i++) { sum+=weights[i]; }
  for (i=0;i<N;i++) {
    //what percent of the ideal distribution is in this bucket?
    data[i].share = weights[i]/(sum/N); 
    data[i].pair = N;
    data[i].idx = i;
  }
  //sort ascending by size
  qsort(data,N, sizeof(wrs_t),wrs_sharecmp);

  int j=N-1; //the biggest bucket
  for (i=0;i<j;i++) {
    int check = i;
    double excess = 1.0 - data[check].share;
    while (excess>0 && i<j) {
      //If this bucket has less samples than a flat distribution,
      //it will be hit more frequently than it should be.  
      //So send excess hits to a bucket which has too many samples.
      data[check].pair=j; 
      // Account for the fact that the paired bucket will be hit more often,
      data[j].share -= excess;  
      excess = 1.0 - data[j].share;
      // If paired bucket now has excess hits, send to new largest bucket at j-1
      if (excess >= 0) { check=j--;} 
    }
  }
  return data;
}


int wrs_pick(wrs_t* collection, size_t N)
//O(1) weighted random sampling (after preparing the collection).
//Randomly select a bucket, and a percentage.
//If the percentage is greater than that bucket's share of hits, 
// use it's paired bucket.
{
  int idx = rand_in_range(0,N);
  double pct = rand_percent();
  if (pct > collection[idx].share) { idx = collection[idx].pair; }
  return collection[idx].idx;
}

Modifica 2:   Dopo un po 'di ricerche, ho scoperto che è persino possibile eseguire la costruzione in O (N). Con un attento monitoraggio, non è necessario ordinare l'array per trovare i contenitori grandi e piccoli. Implementazione aggiornata qui

Se devi prelevare un gran numero di campioni con una distribuzione ponderata di punti discreti, guarda una risposta a una domanda simile .

Tuttavia, se hai bisogno di approssimare alcune funzioni casuali continue usando un istogramma, allora la tua scommessa migliore è probabilmente la risposta di integrazione numerica di dmckee. In alternativa, puoi utilizzare l'aliasing, memorizzare il punto a sinistra e scegliere un numero uniforme tra i due punti.

Per scegliere da un istogramma (originale o ridotto), Metodo alias di Walker è veloce e semplice.

Per una distribuzione normale, potrebbe essere utile quanto segue:

http://en.wikipedia.org/wiki/Normal_distribution#Generating_valuesavaria_al_aria

Non affiliato a StackOverflow
scroll top