studio - superponer graficas en r

¿Cómo puedo generar puntos que coincidan con un histograma?

Question

Estoy trabajando en un sistema de simulación. Pronto tendré datos experimentales (histogramas) para la distribución de valores en el mundo real para varias entradas de simulación.

Cuando se ejecuta la simulación, me gustaría poder producir valores aleatorios que coincidan con la distribución medida. Prefiero hacer esto sin almacenar los histogramas originales. ¿Cuáles son algunas buenas maneras de

1. ¿Asigna un histograma a un conjunto de parámetros que representan la distribución?
2. ¿Generando valores basados ​​en esos parámetros en tiempo de ejecución?

EDITAR: Los datos de entrada son duraciones de eventos para varios tipos diferentes de eventos. Espero que los diferentes tipos tengan diferentes funciones de distribución.

Answer 1

Más información sobre el problema sería útil. Por ejemplo, ¿sobre qué tipo de valores están los histogramas? ¿Son categóricos (por ejemplo, colores, letras) o continuos (por ejemplo, alturas, tiempo)?

Si los histogramas son datos categóricos, creo que puede ser difícil parametrizar las distribuciones a menos que haya muchas correlaciones entre las categorías.

Si los histogramas son datos continuos, puede intentar ajustar la distribución usando mezclas de gaussianos. Es decir, intente ajustar el histograma usando $ / sum_ {i = 1} ^ n w_i N (m_i, v_i) $ donde m_i y v_i son la media y la varianza. Luego, cuando desee generar datos, primero muestree una i de 1..n con una probabilidad proporcional a los pesos w_i y luego muestree un x ~ n (m_i, v_i) como lo haría con cualquier gaussiano.

De cualquier manera, es posible que desee leer más acerca de los modelos de mezcla .

Answer 2

Para una distribución normal, lo siguiente puede ayudar:

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

Answer 3

Al menos dos opciones:

1. Integre el histograma e invierta numéricamente.
2. Rechazo

Integración numérica

De la Computación en Física Moderna por William R. Gibbs:

Uno siempre puede integrar numéricamente [la función] e invertir el [ cdf ], pero a menudo no es muy satisfactorio, especialmente si el pdf está cambiando rápidamente.

Literalmente construyes una tabla que traduce el rango [0-1) en los rangos apropiados en la distribución objetivo. Luego arroje su PRNG habitual (alta calidad) y traduzca con la tabla. Es engorroso, pero claro, viable y completamente general.

Rechazo:

Normalice el histograma del objetivo, luego

1. Tira los dados para elegir una posición ( x ) a lo largo del rango al azar.
2. Tira de nuevo y selecciona este punto si el nuevo número aleatorio es menor que el histograma normalizado en este contenedor. De lo contrario goto (1).

De nuevo, de mente sencilla pero clara y funcional. Puede ser lento para la distribución con mucha probabilidad muy baja (picos con colas largas).

Con estos dos métodos, puede aproximar los datos con ajustes polinómicos por partes o splines para generar una curva suave si no se desea un histograma de función por pasos, pero deje eso para más adelante ya que puede ser una optimización prematura.

Mejores métodos pueden existir para casos especiales.

Todo esto es bastante estándar y debería aparecer en cualquier libro de texto de Análisis numérico si necesito más detalles.

Answer 4

Para elegir un histograma (original o reducido), el método de alias de Walker es rápido y simple.

Answer 5

Entonces parece que lo que quiero para generar una distribución de probabilidad dada es una Función Cuantil , que es la inversa de la función de distribución acumulativa , como dice @dmckee.

La pregunta es: ¿Cuál es la mejor forma de generar y almacenar una función cuantílica que describa un histograma continuo dado? Tengo la sensación de que la respuesta dependerá en gran medida de la forma de la entrada; si sigue cualquier tipo de patrón, debería haber simplificaciones sobre el caso más general. Actualizaré aquí mientras voy.

Editar:

Tuve una conversación esta semana que me recordó este problema. Si renuncio a describir el histograma como una ecuación, y solo almaceno la tabla, ¿puedo hacer selecciones en O (1) vez? Resulta que puede, sin pérdida de precisión, a costa de O (N lgN) tiempo de construcción.

Crea una matriz de N elementos. Una selección aleatoria uniforme en la matriz encontrará un elemento con probabililidad 1 / N. Para cada elemento, almacene la fracción de visitas para las cuales se debe seleccionar este elemento, y el índice de otro elemento que se seleccionará si este no se selecciona.

Muestreo aleatorio ponderado, implementación 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;
} 

Editar 2: Después de un poco de investigación, descubrí que incluso es posible hacer la construcción en el tiempo O (N). Con un seguimiento cuidadoso, no necesita ordenar la matriz para encontrar los contenedores grandes y pequeños. Implementación actualizada aquí

Answer 6

Si necesita extraer una gran cantidad de muestras con una distribución ponderada de puntos discretos, mire la respuesta a una pregunta similar .

Sin embargo, si necesita aproximar alguna función aleatoria continua utilizando un histograma, entonces su mejor apuesta es probablemente la respuesta de integración numérica de dmckee. Alternativamente, puede usar el aliasing y almacenar el punto a la izquierda, y elegir un número uniforme entre los dos puntos.