Área bajo el "pdf" en la estimación de densidad del núcleo en R

Estoy tratando de usar la función ' densidad ' en R para hacer estimaciones de densidad del núcleo. Tengo algunas dificultades para interpretar los resultados y comparar varios conjuntos de datos, ya que parece que el área bajo la curva no es necesariamente 1. Para cualquier función de densidad de probabilidad (pdf) , necesitamos tener el área . Supongo que la estimación de densidad del núcleo informa el pdf. Estoy usando integrate.xy de sfsmisc para estimar el área bajo la curva. $\phi(x)$ $\int_{-\infty}^\infty \phi(x) dx = 1$

> # generate some data
> xx<-rnorm(10000)
> # get density
> xy <- density(xx)
> # plot it
> plot(xy)

trama de la densidad

> # load the library
> library(sfsmisc)
> integrate.xy(xy$x,xy$y)
[1] 1.000978
> # fair enough, area close to 1
> # use another bw
> xy <- density(xx,bw=.001)
> plot(xy)

densidad con bw = .001

> integrate.xy(xy$x,xy$y)
[1] 6.518703
> xy <- density(xx,bw=1)
> integrate.xy(xy$x,xy$y)
[1] 1.000977
> plot(xy)

densidad con bw = 1

> xy <- density(xx,bw=1e-6)
> integrate.xy(xy$x,xy$y)
[1] 6507.451
> plot(xy)

densidad con bw = 1e-6

¿No debería ser el área bajo la curva siempre 1? Parece que los pequeños anchos de banda son un problema, pero a veces desea mostrar los detalles, etc. en las colas y se necesitan pequeños anchos de banda.

Actualización / respuesta:

$2^{20}$

> xy <- density(xx,n=2^15,bw=.001)
> plot(xy)

densidad con mayor número de puntos para muestrear

> integrate.xy(xy$x,xy$y)
[1] 1.000015
> xy <- density(xx,n=2^20,bw=1e-6)
> integrate.xy(xy$x,xy$y)
[1] 2.812398

r estimation pdf kernel-smoothing auc highBandWidth
fuente

Esto parece una limitación de coma flotante en densidad (): al usar un ancho de banda de 1e-6, está creando (en teoría) una colección de 10,000 picos, cada uno de masa total 1/10000. Esos picos terminan siendo representados principalmente por sus picos, sin que las brechas se caractericen adecuadamente. Simplemente estás empujando la densidad () más allá de sus límites.

whuber

@whuber, por limitación de punto flotante, se refiere a límites de precisión, ya que el uso de flotadores conduciría a una mayor sobreestimación del error en comparación con el uso de dobles. No creo ver cómo sucedería eso, pero me gustaría ver alguna evidencia.

highBandWidth

n

$n$ en el cálculo de densidad.

whuber

1

$1$ ?

HA SALIDO - Anony-Mousse

@ Anony-Mousse, sí, eso es lo que esta pregunta hace. ¿Por qué no está evaluando a 1?

highBandWidth

Respuestas:

Piensa en los integrate.xy()usos de la regla trapezoidal . Para la distribución normal, subestimará el área debajo de la curva en el intervalo (-1,1) donde la densidad es cóncava (y, por lo tanto, la interpolación lineal está por debajo de la densidad real), y la sobreestimará en otra parte (como va la interpolación lineal encima de la verdadera densidad). Como la última región es más grande (en la medida de Lesbegue, si lo desea), la regla trapezoidal tiende a sobreestimar la integral. Ahora, a medida que avanza a anchos de banda más pequeños, casi toda su estimación es convexa por partes, con muchos picos estrechos correspondientes a los puntos de datos y valles entre ellos. Ahí es donde la regla trapezoidal se rompe especialmente mal.

StasK
fuente

eso significa que estamos "sobremuestreando" los picos y "submuestreando" los valles, en cierto sentido ondulado. Dado que la visualización también sigue la regla trapezoidal (interpolación lineal entre muestras), parece que un ancho de banda de kernel demasiado pequeño también es malo para la visualización. Además, si pudiéramos obtener un mayor número de puntos en los que calculamos la densidad, habría menos problemas.

highBandWidth

Esta explicación no retiene el agua. El problema es que la densidad está inadecuadamente discretizada, no que la regla trapezoidal se descomponga gravemente. Integrar () no puede obtener una respuesta correcta porque la densidad () no produce una representación correcta. Para ver esto, solo inspeccione xy $ x: ¡tiene solo 512 valores destinados a representar 10,000 picos estrechos!

whuber

@whuber, eso es lo que decía la respuesta. El punto es que debe usar la regla trapezoidal para un número finito de muestras, y sobreestima el área en comparación con la densidad real en un eje continuo de acuerdo con los núcleos. Mi actualización al final de la pregunta se amplía.

highBandWidth

@high No; La regla trapezoidal está funcionando bien. El problema es que está trabajando con una discretización incorrecta del integrando. ¡No puede tener "muchos picos estrechos correspondientes a los puntos de datos" cuando hay 10,000 puntos de datos y solo 512 valores en la matriz de densidad!

whuber

Mirando estos gráficos, ahora estoy pensando que el problema es con densitymás que con integrate.xy. Con N = 10 000 y BW = 1E6, sería tener para ver un peine con una altura de cada diente de aproximadamente 1E6, y siendo más densa alrededor de 0. En lugar de los dientes, todavía ver una curva en forma de campana reconocible. Entonces, te densityestá engañando, o al menos debería usarse de manera diferente con anchos de banda pequeños: ndebería ser sobre (rango de datos) / (bw) en lugar del valor predeterminado n=512. El integrador debe estar recogiendo uno de estos enormes valores que densityregresa por una infeliz coincidencia.

StasK

-1

Está bien, puedes arreglarlo cambiando y escalando; agregue el número más pequeño de modo que la densidad no sea negativa, luego multiplique todo por una constante tal que el área sea la unidad. Este es el camino fácil.

$L_2$ $c$ $\left[\phi(x)-c\right]^+$

Emre
fuente

Tenga en cuenta que la pregunta es más bien por qué la densityfunción no produce la densidad "adecuada" que se integra a 1, en lugar de cómo solucionarla.

Tim