Distância mínima esperada de um ponto com densidade variável

8

Estou vendo como a distância euclidiana mínima esperada entre pontos aleatoriamente uniformes e a origem muda à medida que aumentamos a densidade de pontos aleatórios ( pontos por unidade quadrada ) ao redor da origem. Eu consegui chegar a um relacionamento entre os dois descritos como tal:

Expected Min Distance=12Density

Eu vim com isso executando algumas simulações de Monte Carlo em R e ajustando uma curva manualmente (código abaixo).

Minha pergunta é : eu poderia ter derivado esse resultado teoricamente, e não através de experimentação?

#Stack Overflow example
library(magrittr)
library(ggplot2)


#---------
#FUNCTIONS
#---------
#gen random points within a given radius and given density
gen_circle_points <- function(radius, density) {
  #round radius up then generate points in square with side length = 2*radius
  c_radius <- ceiling(radius)
  coords <- data.frame(
    x = runif((2 * c_radius) ^ 2 * density, -c_radius, c_radius),
    y = runif((2 * c_radius) ^ 2 * density, -c_radius, c_radius)
  )
  return(coords[sqrt(coords$x ^ 2 + coords$y ^ 2) <= radius, ])#filter in circle
}

#Example plot
plot(gen_circle_points(radius = 1,density = 200)) #200 points around origin
points(0,0, col="red",pch=19) #colour origin

insira a descrição da imagem aqui

#return euclidean distances of points generated by gen_circle_points()
calculate_distances <- function(circle_points) {
  return(sqrt(circle_points$x ^ 2 + circle_points$y ^ 2))
}

#find the smallest distance from output of calculate_distances()
calculate_min_value <- function(distances) {
  return(min(distances))
}


#Try a range of values
density_values <- c(1:100)

expected_min_from_density <- sapply(density_values, function(density) {
  #simulate each density value 1000 times and take an average as estimate for
  #expected minimum distance
  sapply(1:1000, function(i) {
    gen_circle_points(radius=1, density=density) %>%
      calculate_distances() %>%
      calculate_min_value()
  }) %>% mean()
})

results <- data.frame(density_values, expected_min_from_density)

#fit based off exploration
theoretical_fit <- data.frame(density = density_values, 
                              fit = 1 / (sqrt(density_values) * 2))

#plot monte carlo (black) and fit (red dashed)
ggplot(results, aes(x = density_values, y = expected_min_from_density)) +
  geom_line() + 
  geom_line(
    data = theoretical_fit,
    aes(x = density, y = fit),
    color = "red",
    linetype = 2
  )

Um gráfico dos valores de densidade em relação ao mínimo esperado, tanto monte carlo quanto teórico

Michael Bird
fonte
A dependência direta (assintótica) da raiz inversa da densidade segue fácil e imediatamente a partir das considerações das unidades de medida; portanto, a única questão diz respeito a por que o múltiplo é1/2.
whuber
@whuber Sim, eu notei que as unidades estavam bem alinhadas e sim, a pergunta é: de onde vieram as duas?
Michael Bird
1
O é a largura do seu quadrado. 2
whuber

Respostas:

8

Considere a distância até a origem de variáveis ​​aleatórias distribuídas independentemente que possuem distribuições uniformes no quadrado( X i , Y i ) [ - 1 , 1 ] 2 .n(Xi,Yi)[1,1]2.

Escrevendo para a distância ao quadrado, a geometria euclidiana nos mostra queRi2=Xi2+Yi2

Pr(Rir1)=14πr2

enquanto (com um pouco mais de trabalho)

Pr(1Rir2)=14(πr2+4r214r2ArcTan(r21)).

Figura 1: Gráfico da função de distribuição

Juntos, estes determinam a função de distribuição comum a todos osFRi.

Como os pontos são independentes, assim como as distâncias onde a função de sobrevivência de énmin ( R i )Ri,min(Ri)

Sn(r)=(1F(r))n,

implicando a menor distância média é

μ(n)=02Sn(r)dr.

Para quase toda a área nesta integral é próxima de portanto podemos aproximar isso comon1,0,

μapprox(n)=01Sn(r)dr=01(1π4r2)ndr.

O erro não é maior que a parte da integral omitida, que por sua vez não é maior que

(21)(1F(1))n=(21)(1π/4)n,

o que obviamente diminui exponencialmente comn.

Por sua vez, podemos aproximar o integrando como

(1π4r2)nexp(12r22/(nπ)).

Até uma constante de normalização, esta é a função de densidade de uma distribuição Normal com média e variância A constante de normalização ausente é0σ2=2/(nπ).

C(n)=12πσ2=12π 2/(nπ)=n2.

Portanto, estendendo a integral de para (que adiciona um erro proporcional a ),1en

μapprox(n)0et2/(2σ2)dt=1C(n)12=1n.

No processo de obtenção dessa aproximação, três erros foram cometidos. Coletivamente, eles estão no máximo na ordem o erro incorrido ao se aproximar de pelo gaussiano.n1,Sn(r)

! [Figura 2: Gráfico de erros de simulação

Esta figura plota vezes a diferença entre e vezes a menor distância média observada em conjuntos de dados simulados separados para cada Como eles diminuem à medida que cresce, isso é evidência de que o erro én1n105n.no(n1/n)=o(n3/2).

Finalmente, o fator da questão deriva do tamanho do quadrado:1/2 a densidade é o número de pontos por unidade de área e o quadrado tem a área , de onden,[1,1]24

2Density=2n/4=n.

Este é o Rcódigo para a simulação:

n.sim <- 1e5  # Size of each simulation
d <- 2        # Dimension
n <- 2^(1:11) # Numbers of points in each simulation
#
# Estimate mean distance to the origin for each `n`.
#
y <- sapply(n, function(n.points) {
  x <- array(runif(d*n.points*n.sim, -1, 1), c(d, n.points, n.sim))
  mean(sqrt(apply(colSums(x^2), 2, min)))
})
#
# Plot the errors (normalized) against `n`.
#
library(ggplot2)
ggplot(data.frame(Log2.n = 1:length(n), Error=sqrt(n)* (1 - y * n^(1/d))),
       aes(Log2.n, Error)) + geom_point() + geom_smooth() 
  ylab("Error * n") + ggtitle("Simulation Means")
whuber
fonte
2
Uau! Que resposta! Muito obrigado, isso é ótimo. Obrigado!
Michael Bird
Olá @whuber, eu estava tentando reproduzir seu e notei que sua equação para não retorna como mostra o gráfico. Quando calculei , obtive que fornece a curva que você forneceu. Você cometeu um erro de digitação? F(r)F(2)1π/4-r(rArcCos(1/r)-Pr(1Rir2)π/4r(rArcCos(1/r)11/r2)
Michael Bird
1
@ Michael Obrigado, há um erro de digitação - mas não é o que você sugere: um dos meus " " deveria ter sido " ". Eu consertei esse. 4r4
whuber