Encontrar extremos locales de una función de densidad usando splines

15

Estoy tratando de encontrar los máximos locales para una función de densidad de probabilidad (encontrada usando R's density método ). No puedo hacer un simple método de "mirar alrededor de los vecinos" (donde uno mira alrededor de un punto para ver si es un máximo local con respecto a sus vecinos) ya que hay un gran volumen de datos. Además, parece más eficiente y genérico usar algo como la interpolación de Spline y luego encontrar las raíces de la primera derivada, en lugar de construir una "mirada a los vecinos" con tolerancia a fallas y otros parámetros.

Entonces, mis preguntas:

  1. Dada una función de splinefun, ¿qué métodos encontrarán los máximos locales?
  2. ¿Hay una manera fácil / estándar de encontrar derivadas de una función devuelta usando splinefun?
  3. ¿Existe una forma mejor / estándar de encontrar los máximos locales de una función de densidad de probabilidad?

Como referencia, a continuación hay una gráfica de mi función de densidad. Otras funciones de densidad con las que estoy trabajando son similares en forma. Debo decir que soy nuevo en R, pero no soy nuevo en programación, por lo que puede haber una biblioteca o paquete estándar para lograr lo que necesito. función de densidad

¡¡Gracias por tu ayuda!!

aaronlevin
fuente
No estoy claro por qué el gran volumen de datos es un problema para el método de "mirar alrededor de los vecinos". density()no estima la densidad para cada dato, estima la densidad en n valores, donde n es un parámetro especificado por el usuario con un valor predeterminado n = 512.
onestop
Mi n para esto es 2 ^ 15 y parece que los datos tienen mucha variación a nivel punto por punto. Intenté escribir un buscador max / min usando algo similar al método de vecindario (vía msExtrema {msProcess}) y solo pude identificar algunos de los máximos, nunca todos, jugando con la configuración de tolerancia.
aaronlevin
2
Mirando el código msExtrema, es un contenedor simple para peaksel splus2Rpaquete, que sería mejor usar directamente si solo desea los máximos locales y no los mínimos locales. No puedo ver por qué usar el valor predeterminado span=3no encontraría todos los máximos locales. Y 2 ^ 15 = 32768 no debería ser lo suficientemente grande como para que la eficiencia sea una gran preocupación.
parada el
La función devuelta por splinefun tiene un argumento "deriv" ​​que es 0 por defecto. Establezca deriv = 1 para la primera derivada.
Cian
1
Hmm, peaksparece tener errores: llama max.colcon la configuración predeterminada de ties.method = "random", que no solo rompe los lazos al azar, sino que también establece una tolerancia relativa de 1e-5 para declarar un empate. Lo primero es confuso, lo último definitivamente no es lo que quieres aquí. peaks()también toma un strictparámetro que está mal documentado y, mirando el código de la función, no hace nada. ¡Ah, las alegrías de las bibliotecas de software aportadas por los usuarios! Es posible que también sea capaz de arreglarlo, aunque, como usted dice que no es nuevo en la programación,
onestop

Respuestas:

14

Lo que quieres hacer se llama detección de picos en quimiometría. Hay varios métodos que puede usar para eso. Aquí solo demuestro un enfoque muy simple.

require(graphics)
#some data
d <- density(faithful$eruptions, bw = "sj")

#make it a time series
ts_y<-ts(d$y)

#calculate turning points (extrema)
require(pastecs)
tp<-turnpoints(ts_y)
#plot
plot(d)
points(d$x[tp$tppos],d$y[tp$tppos],col="red")
Roland
fuente
De todas las soluciones, esto funcionó mejor. 1. Pregunta de seguimiento: ¿hay alguna manera de alternar la tolerancia con los puntos de inflexión? Encontré muchos picos y valles en la porción de cola larga de la función de densidad. 2. Pregunta de seguimiento # 2: ¿cuál es una buena manera de determinar la tolerancia?
aaronlevin
anuncio 1. No lo creo. Está destinado a probar la aleatoriedad de series de tiempo, por lo que la función no necesita eso. Puede intentar evaluar la relevancia / importancia de un pico usted mismo. Por ejemplo, podrías hacer una prueba t contra el vecindario (donde puedes decidir qué tan grande debería ser el vecindario). O puede buscar una función más sofisticada en los paquetes R para evaluar los datos de la espectrometría (de masas) u otros métodos de química analítica.
Roland