math - puntos - numeros aleatorios gaussianos en matlab
Genera un punto aleatorio dentro de un círculo(uniformemente) (18)
1) Elige una X al azar entre -1 y 1.
var X:Number = Math.random() * 2 - 1;
2) Usando la fórmula del círculo, calcule los valores máximo y mínimo de Y dado que X y un radio de 1:
var YMin:Number = -Math.sqrt(1 - X * X);
var YMax:Number = Math.sqrt(1 - X * X);
3) Elige una Y aleatoria entre esos extremos:
var Y:Number = Math.random() * (YMax - YMin) + YMin;
4) Incorpore su ubicación y valores de radio en el valor final:
var finalX:Number = X * radius + pos.x;
var finalY:Number = Y * radois + pos.y;
Necesito generar un punto uniformemente aleatorio dentro de un círculo de radio R.
Me doy cuenta de que simplemente eligiendo un ángulo uniformemente aleatorio en el intervalo [0 ... 2π), y un radio uniformemente aleatorio en el intervalo (0 ... R ) terminaría con más puntos hacia el centro, ya que para dos dado radios, los puntos en el radio más pequeño estarán más cerca uno del otro que para los puntos en el radio más grande.
Encontré una entrada de blog sobre esto aquí pero no entiendo su razonamiento. Supongo que es correcto, pero realmente me gustaría entender de dónde obtiene (2 / R 2 ) × r y cómo obtiene la solución final.
Con respecto al muestreo de rechazo: podría generar un punto aleatorio dentro del cuadrado R × R una y otra vez hasta que obtenga uno dentro del círculo. Este enfoque tiene el obvio inconveniente de que no proporciona una garantía para la terminación (aunque es muy poco probable que continúe por mucho tiempo).
Acerquémonos a esto como lo hubiera hecho Arquímedes.
¿Cómo podemos generar un punto uniformemente en un triángulo ABC, donde | AB | = | BC |? Hagamos esto más fácil extendiéndonos a un paralelogramo ABCD. Es fácil generar puntos de manera uniforme en ABCD. Seleccionamos uniformemente un punto aleatorio X en AB e Y en BC y elegimos Z de manera que XBYZ es un paralelogramo. Para obtener un punto elegido uniformemente en el triángulo original, simplemente doblamos cualquier punto que aparezca en ADC de nuevo hacia ABC a lo largo de AC.
Ahora considera un círculo. En el límite podemos pensar en infinitamente muchos triángulos isocelos ABC con B en el origen y A y C en la circunferencia que se aproximan de forma muy cercana. Podemos elegir uno de estos triángulos simplemente eligiendo un ángulo theta. Entonces, necesitamos generar una distancia desde el centro eligiendo un punto en la franja ABC. Nuevamente, extienda a ABCD, donde D ahora tiene el doble del radio del centro del círculo.
Escoger un punto aleatorio en ABCD es fácil usando el método anterior. Elige un punto al azar en AB. Uniformemente escoja un punto al azar en BC. Es decir. escoge un par de números aleatorios xey uniformemente en [0, R] dando distancias desde el centro. Nuestro triángulo es una astilla delgada por lo que AB y BC son esencialmente paralelos. Entonces, el punto Z es simplemente una distancia x + y desde el origen. Si x + y> R, volvemos a doblar.
Aquí está el algoritmo completo para R = 1. Espero que estés de acuerdo, es bastante simple. Utiliza trigonometría, pero puede ofrecer una garantía sobre cuánto tiempo llevará y cuántas llamadas al random() necesita, a diferencia del muestreo de rechazo.
t = 2*pi*random()
u = random()+random()
r = if u>1 then 2-u else u
[r*cos(t), r*sin(t)]
Aquí está en Mathematica.
f[] := Block[{u, t, r},
u = Random[] + Random[];
t = Random[] 2 Pi;
r = If[u > 1, 2 - u, u];
{r Cos[t], r Sin[t]}
]
ListPlot[Table[f[], {10000}], AspectRatio -> Automatic]
Aquí está mi código de Python para generar num aleatorios de un círculo de radio rad :
import matplotlib.pyplot as plt
import numpy as np
rad = 10
num = 1000
t = np.random.uniform(0.0, 2.0*np.pi, num)
r = rad * np.sqrt(np.random.uniform(0.0, 1.0, num))
x = r * np.cos(t)
y = r * np.sin(t)
plt.plot(x, y, "ro", ms=1)
plt.axis([-15, 15, -15, 15])
plt.show()
Aquí hay una solución rápida y simple.
Elija dos números aleatorios en el rango (0, 1), a saber, a y b . Si b < a , cambiarlos. Su punto es (b*R*cos(2*pi*a/b), b*R*sin(2*pi*a/b)) .
Puedes pensar en esta solución de la siguiente manera. Si tomó el círculo, lo cortó, luego lo enderezó, obtendría un triángulo rectángulo. Escala ese triángulo hacia abajo, y tendrás un triángulo de (0, 0) a (1, 0) a (1, 1) y de nuevo a (0, 0) . Todas estas transformaciones cambian la densidad de manera uniforme. Lo que ha hecho se selecciona de manera uniforme como un punto aleatorio en el triángulo e invierte el proceso para obtener un punto en el círculo.
Creo que en este caso el uso de coordenadas polares es una forma de complicar el problema, sería mucho más fácil si selecciona puntos aleatorios en un cuadrado de 2R lados y luego selecciona los puntos (x, y) tales que x ^ 2 + y ^ 2 <= R ^ 2.
Deje ρ (radio) y φ (acimut) ser dos variables aleatorias correspondientes a las coordenadas polares de un punto arbitrario dentro del círculo. Si los puntos están distribuidos uniformemente, ¿cuál es la función de distribución de ρ y φ?
Para cualquier r: 0 <r <R la probabilidad de coordenada de radio ρ sea menor que r es
P [ρ <r] = P [el punto está dentro de un círculo de radio r] = S1 / S0 = (r / R) 2
Donde S1 y S0 son las áreas del círculo de radio r y R, respectivamente. Entonces la CDF se puede dar como:
0 if r<=0
CDF = (r/R)**2 if 0 < r <= R
1 if r > R
Y PDF:
PDF = d/dr(CDF) = 2 * (r/R**2) (0 < r <= R).
Tenga en cuenta que para R = 1 la variable aleatoria sqrt (X) donde X es uniforme en [0, 1) tiene esta FCD exacta (porque P [sqrt (X) <y] = P [x <y ** 2] = y * * 2 para 0 <y <= 1).
La distribución de φ es obviamente uniforme de 0 a 2 * π. Ahora puedes crear coordenadas polares al azar y convertirlas a cartesianas usando ecuaciones trigonométricas:
x = ρ * cos(φ)
y = ρ * sin(φ)
No puedo resistirme a publicar el código python para R = 1.
from matplotlib import pyplot as plt
import numpy as np
rho = np.sqrt(np.random.uniform(0, 1, 5000))
phi = np.random.uniform(0, 2*np.pi, 5000)
x = rho * np.cos(phi)
y = rho * np.sin(phi)
plt.scatter(x, y, s = 4)
Conseguirás
El elemento de área en un círculo es dA = rdr * dphi. Ese factor adicional destruyó tu idea de elegir aleatoriamente ar y phi. Mientras phi se distribuye plano, r no lo es, pero plano en 1 / r (es decir, es más probable que llegue al límite que "la diana").
Entonces, para generar puntos distribuidos uniformemente sobre el círculo, seleccione phi de una distribución plana yr de una distribución 1 / r.
Alternativamente, use el método Monte Carlo propuesto por Mehrdad.
EDITAR
Para elegir un plano aleatorio en 1 / r, puede elegir un x aleatorio del intervalo [1 / R, infinito] y calcular r = 1 / x. r se distribuye plano en 1 / r.
Para calcular una phi aleatoria, elija una x aleatoria del intervalo [0, 1] y calcule phi = 2 * pi * x.
Existe una relación lineal entre el radio y el número de puntos "cercanos" a ese radio, por lo que debe usar una distribución de radio que también hace que el número de puntos de datos cerca de un radio r proporcional a r .
La razón por la cual la solución ingenua no funciona es porque da una mayor densidad de probabilidad a los puntos más cercanos al centro del círculo. En otras palabras, el círculo que tiene radio r / 2 tiene probabilidad r / 2 de obtener un punto seleccionado en él, pero tiene área (número de puntos) pi * r ^ 2/4.
Por lo tanto, queremos que una densidad de probabilidad de radio tenga la siguiente propiedad:
La probabilidad de elegir un radio menor o igual a un r dado tiene que ser proporcional al área del círculo con radio r. (porque queremos tener una distribución uniforme en los puntos y áreas más grandes significa más puntos)
En otras palabras, queremos que la probabilidad de elegir un radio entre [0, r] sea igual a su parte del área total del círculo. El área total del círculo es pi * R ^ 2, y el área del círculo con radio r es pi * r ^ 2. Por lo tanto, nos gustaría que la probabilidad de elegir un radio entre [0, r] sea (pi * r ^ 2) / (pi * R ^ 2) = r ^ 2 / R ^ 2.
Ahora viene la matemática:
La probabilidad de elegir un radio entre [0, r] es la integral de p (r) dr de 0 a r (eso es solo porque agregamos todas las probabilidades de los radios más pequeños). Por lo tanto, queremos integral (p (r) dr) = r ^ 2 / R ^ 2. Podemos ver claramente que R ^ 2 es una constante, por lo que todo lo que tenemos que hacer es descubrir qué p (r), cuando se integra, nos daría algo como r ^ 2. La respuesta es claramente r * constante. integral (r * constante dr) = r ^ 2/2 * constante. Esto tiene que ser igual a r ^ 2 / R ^ 2, por lo tanto constante = 2 / R ^ 2. Por lo tanto, tiene la distribución de probabilidad p (r) = r * 2 / R ^ 2
Nota: Otra forma más intuitiva de pensar sobre el problema es imaginar que estás tratando de dar a cada círculo de radio con una densidad de probabilidad igual a la proporción del número de puntos que tiene en su circunferencia. Por lo tanto, un círculo que tiene radio r tendrá 2 * pi * r "puntos" en su circunferencia. El número total de puntos es pi * R ^ 2. Por lo tanto, debe dar al círculo una probabilidad de ra igual a (2 * pi * r) / (pi * R ^ 2) = 2 * r / R ^ 2. Esto es mucho más fácil de entender y más intuitivo, pero no es tan matemáticamente sólido.
No sé si esta pregunta aún está abierta para una nueva solución con todas las respuestas ya dadas, pero casualmente me he enfrentado a la misma pregunta. Traté de "razonar" conmigo mismo para encontrar una solución, y encontré una. Podría ser lo mismo que algunos ya han sugerido aquí, pero de todos modos aquí está:
para que dos elementos de la superficie del círculo sean iguales, suponiendo dr iguales, debemos tener dtheta1 / dtheta2 = r2 / r1. Escribir la expresión de la probabilidad para ese elemento como P (r, theta) = P {r1 <r <r1 + dr, theta1 <theta <theta + dtheta1} = f (r, theta) * dr * dtheta1, y establecer los dos las probabilidades (para r1 y r2) son iguales, llegamos a (suponiendo que r y theta son independientes) f (r1) / r1 = f (r2) / r2 = constante, lo que da f (r) = c * r. Y el resto, determinar la constante c se deriva de la condición en que f (r) es un PDF.
Piensa en ello de esta manera. Si tiene un rectángulo donde un eje es un radio y uno es un ángulo, y toma los puntos dentro de este rectángulo que están cerca del radio 0. Estos caerán muy cerca del origen (que está muy cerca el uno del otro). Sin embargo, los puntos cerca del radio R, todos caerán cerca del borde del círculo (es decir, muy separados el uno del otro).
Esto podría darle una idea de por qué está obteniendo este comportamiento.
El factor que se deriva de ese enlace le indica la cantidad de área correspondiente en el rectángulo que necesita ajustarse para que no dependa del radio una vez que esté asignado al círculo.
Editar: Entonces, lo que escribe en el enlace que comparte es: "Eso es bastante fácil de hacer calculando el inverso de la distribución acumulativa, y obtenemos para r:".
La premisa básica es que aquí puede crear una variable con una distribución deseada a partir de un uniforme mapeando el uniforme mediante la función inversa de la función de distribución acumulativa de la función de densidad de probabilidad deseada. ¿Por qué? Simplemente dalo por hecho por ahora, pero esto es un hecho.
Aquí está mi explicación intuitiva de las matemáticas. La función de densidad f (r) con respecto a r tiene que ser proporcional a la r misma. Comprender este hecho es parte de cualquier libro de cálculo básico. Ver secciones sobre elementos de área polar. Algunos otros carteles han mencionado esto.
Entonces lo llamaremos f (r) = C * r;
Esto resulta ser la mayor parte del trabajo. Ahora, como f (r) debe ser una densidad de probabilidad, puede ver fácilmente que al integrar f (r) sobre el intervalo (0, R) se obtiene que C = 2 / R ^ 2 (esto es un ejercicio para el lector) .)
Por lo tanto, f (r) = 2 * r / R ^ 2
OK, así es como se obtiene la fórmula en el enlace.
Entonces, la parte final va desde la variable aleatoria uniforme u en (0,1) debe mapear por la función inversa de la función de distribución acumulada desde esta densidad deseada f (r). Para entender por qué este es el caso, necesita encontrar un texto de probabilidad avanzada como Papoulis probablemente (o derivarlo usted mismo).
Integrando f (r) obtienes F (r) = r ^ 2 / R ^ 2
Para encontrar la función inversa de esto, establece u = r ^ 2 / R ^ 2 y luego resuelve para r, lo que le da r = R * sqrt (u)
Esto también tiene sentido intuitivamente, u = 0 debería corresponder a r = 0. Además, u = 1 shoudl map a r = R. Además, va por la función de raíz cuadrada, que tiene sentido y coincide con el enlace.
Primero generamos un cdf [x] que es
La probabilidad de que un punto esté a menos de la distancia x del centro del círculo. Suponga que el círculo tiene un radio de R.
obviamente si x es cero entonces cdf [0] = 0
obviamente si x es R entonces el cdf [R] = 1
obviamente si x = r entonces el cdf [r] = (Pi r ^ 2) / (Pi R ^ 2)
Esto se debe a que cada "área pequeña" en el círculo tiene la misma probabilidad de ser elegido, por lo que la probabilidad es proporcional al área en cuestión. Y el área dada una distancia x desde el centro del círculo es Pi r ^ 2
entonces cdf [x] = x ^ 2 / R ^ 2 porque los Pi se cancelan unos a otros
tenemos cdf [x] = x ^ 2 / R ^ 2 donde x va de 0 a R
Entonces resolvemos para x
R^2 cdf[x] = x^2
x = R Sqrt[ cdf[x] ]
Ahora podemos reemplazar cdf con un número aleatorio de 0 a 1
x = R Sqrt[ RandomReal[{0,1}] ]
Finalmente
r = R Sqrt[ RandomReal[{0,1}] ];
theta = 360 deg * RandomReal[{0,1}];
{r,theta}
obtenemos las coordenadas polares {0.601168 R, 311.915 grados}
Realmente depende de lo que quieres decir con "uniformemente aleatorio". Este es un punto sutil y puede leer más sobre él en la página wiki aquí: http://en.wikipedia.org/wiki/Bertrand_paradox_%28probability%29 , donde el mismo problema, dando diferentes interpretaciones a ''uniformemente aleatorio'' da diferentes respuestas!
Dependiendo de cómo elija los puntos, la distribución puede variar, aunque sean uniformemente aleatorios en algún sentido.
Parece que la entrada del blog está tratando de hacerlo uniformemente aleatorio en el siguiente sentido: si tomas un subcirculo del círculo, con el mismo centro, entonces la probabilidad de que el punto caiga en esa región es proporcional al área de la región. Eso, creo, está intentando seguir la interpretación ahora estándar de ''uniformemente aleatorio'' para regiones 2D con áreas definidas en ellas : la probabilidad de que un punto caiga en cualquier región (con un área bien definida) es proporcional al área de esa región.
Solución en Java y el ejemplo de distribución (2000 puntos)
public void getRandomPointInCircle() {
double t = 2 * Math.PI * Math.random();
double r = Math.sqrt(Math.random());
double x = r * Math.cos(t);
double y = r * Math.sin(t);
System.out.println(x);
System.out.println(y);
}
basado en la solución previa https://.com/a/5838055/5224246 de @sigfpe
Tenga en cuenta la densidad del punto en proporcional al cuadrado inverso del radio, por lo tanto, en lugar de elegir r de [0, r_max] , elija de [0, r_max^2] , luego calcule sus coordenadas como:
x = sqrt(r) * cos(angle)
y = sqrt(r) * sin(angle)
Esto le dará distribución uniforme de puntos en un disco.
Todavía no estoy seguro del ''(2 / R2) × r'' exacto, pero lo que es evidente es el número de puntos requeridos para ser distribuidos en la unidad dada ''dr'', es decir, el aumento en r será proporcional a r2 y no a r.
compruebe de esta manera ... el número de puntos en algún ángulo theta y entre r (0.1r a 0.2r) es decir, la fracción de ry el número de puntos entre r (0.6r a 0.7r) sería igual si usa generación estándar, ya que la diferencia es solo 0.1r entre dos intervalos. pero dado que el área cubierta entre los puntos (0.6r a 0.7r) será mucho más grande que el área cubierta entre 0.1r a 0.2r, la misma cantidad de puntos estará espaciadamente espaciada en un área más grande, esto supongo que usted ya lo sabe, entonces la función generar los puntos aleatorios no debe ser lineal sino cuadrático (dado que el número de puntos que se debe distribuir en la unidad dada ''dr'', es decir, aumentar en r será proporcional a r2 y no r), en este caso será inverso cuadrático, ya que el delta que tenemos (0.1r) en ambos intervalos debe ser cuadrado de alguna función para que pueda actuar como valor inicial para la generación lineal de puntos (ya que afterwords, esta semilla se usa linealmente en las funciones sen y cos), entonces Saber, dr debe ser un valor cuadrático y para hacer que esta semilla sea cuadrática, necesitamos originar estos valores a partir de la raíz cuadrada de r not r en sí, espero que esto lo haga un poco más claro.
Una solución de programador:
- Crea un mapa de bits (una matriz de valores booleanos). Puede ser tan grande como quieras
- Dibuja un círculo en ese mapa de bits.
- Crea una tabla de búsqueda de los puntos del círculo.
- Elija un índice aleatorio en esta tabla de búsqueda.
const int RADIUS = 64;
const int MATRIX_SIZE = RADIUS * 2;
bool matrix[MATRIX_SIZE][MATRIX_SIZE] = {0};
struct Point { int x; int y; };
Point lookupTable[MATRIX_SIZE * MATRIX_SIZE];
void init()
{
int numberOfOnBits = 0;
for (int x = 0 ; x < MATRIX_SIZE ; ++x)
{
for (int y = 0 ; y < MATRIX_SIZE ; ++y)
{
if (x * x + y * y < RADIUS * RADIUS)
{
matrix[x][y] = true;
loopUpTable[numberOfOnBits].x = x;
loopUpTable[numberOfOnBits].y = y;
++numberOfOnBits;
} // if
} // for
} // for
} // ()
Point choose()
{
int randomIndex = randomInt(numberOfBits);
return loopUpTable[randomIndex];
} // ()
El mapa de bits solo es necesario para la explicación de la lógica. Este es el código sin el mapa de bits:
const int RADIUS = 64;
const int MATRIX_SIZE = RADIUS * 2;
struct Point { int x; int y; };
Point lookupTable[MATRIX_SIZE * MATRIX_SIZE];
void init()
{
int numberOfOnBits = 0;
for (int x = 0 ; x < MATRIX_SIZE ; ++x)
{
for (int y = 0 ; y < MATRIX_SIZE ; ++y)
{
if (x * x + y * y < RADIUS * RADIUS)
{
loopUpTable[numberOfOnBits].x = x;
loopUpTable[numberOfOnBits].y = y;
++numberOfOnBits;
} // if
} // for
} // for
} // ()
Point choose()
{
int randomIndex = randomInt(numberOfBits);
return loopUpTable[randomIndex];
} // ()
Utilicé una vez este método: esto puede ser totalmente no optimizado (es decir, usa una matriz de puntos por lo que no se puede usar en círculos grandes) pero da una distribución aleatoria suficiente. Puede omitir la creación de la matriz y dibujar directamente si lo desea. El método es aleatorizar todos los puntos en un rectángulo que caen dentro del círculo.
bool[,] getMatrix(System.Drawing.Rectangle r) {
bool[,] matrix = new bool[r.Width, r.Height];
return matrix;
}
void fillMatrix(ref bool[,] matrix, Vector center) {
double radius = center.X;
Random r = new Random();
for (int y = 0; y < matrix.GetLength(0); y++) {
for (int x = 0; x < matrix.GetLength(1); x++)
{
double distance = (center - new Vector(x, y)).Length;
if (distance < radius) {
matrix[x, y] = r.NextDouble() > 0.5;
}
}
}
}
private void drawMatrix(Vector centerPoint, double radius, bool[,] matrix) {
var g = this.CreateGraphics();
Bitmap pixel = new Bitmap(1,1);
pixel.SetPixel(0, 0, Color.Black);
for (int y = 0; y < matrix.GetLength(0); y++)
{
for (int x = 0; x < matrix.GetLength(1); x++)
{
if (matrix[x, y]) {
g.DrawImage(pixel, new PointF((float)(centerPoint.X - radius + x), (float)(centerPoint.Y - radius + y)));
}
}
}
g.Dispose();
}
private void button1_Click(object sender, EventArgs e)
{
System.Drawing.Rectangle r = new System.Drawing.Rectangle(100,100,200,200);
double radius = r.Width / 2;
Vector center = new Vector(r.Left + radius, r.Top + radius);
Vector normalizedCenter = new Vector(radius, radius);
bool[,] matrix = getMatrix(r);
fillMatrix(ref matrix, normalizedCenter);
drawMatrix(center, radius, matrix);
}