Interpretación de Hartigans ' dip test
On noviembre 18, 2020 by adminMe gustaría encontrar una manera de cuantificar la intensidad de la bimodalidad de algunos distribuciones que obtuve empíricamente. Por lo que leí, todavía hay cierto debate sobre la forma de cuantificar la bimodalidad. Elegí usar la prueba de inmersión de Hartigans, que parece ser la única disponible en R (documento original: http://www.stat.washington.edu/wxs/Stat593-s03/Literature/hartigan85a.pdf ). La prueba de inmersión de Hartigans se define como: «La prueba de inmersión mide la multimodalidad en una muestra por la diferencia máxima, sobre todos los puntos de la muestra, entre la función de distribución empírica y la función de distribución unimodal que minimiza esa diferencia máxima» .
Me gustaría entender completamente cómo debo interpretar estas estadísticas antes de usarlas. Esperaba que la prueba de inmersión aumentara si la distribución es multimodal (como se define como «la diferencia máxima de la distribución unimodal»). Pero : puedes leer en la página de wikipedia acerca de la distribución multimodal que «Los valores inferiores a 0.05 indican valores y bimodalidad significativos mayor que 0.05 pero menor que 0.10 sugiere bimodalidad con importancia marginal «. . Dicha declaración proviene de este artículo (Fig. 2). Según este artículo, el índice de la prueba de inmersión es cercano a 0 cuando la distribución es bimodal. Me confunde.
Para interpretar correctamente la prueba de inmersión de Hartigans, construí algunas distribuciones (el código original es de aquí ) y aumenté el valor de exp (mu2) (llamado «Intensidad de bimodularidad» de ahora en adelante – Editar: debería haberlo llamado «Intensidad de bimodalidad» ) para obtener bimodalidad. En el primer gráfico, puedes ver algún ejemplo de distribuciones. Luego calculé el índice diptest (segundo gráfico) y el valor p (tercer gráfico) asociado (paquete diptest ) a esas diferentes distribuciones simuladas. El código R utilizado está al final de mi publicación.
Lo que muestro aquí es que el índice de prueba de inmersión es alto y el valor P es bajo cuando las distribuciones son bimodales. Lo cual es contrario a lo que se puede leer en Internet.
No soy un experto en estadísticas, por lo que apenas entendí el artículo de Hartigans. Me gustaría recibir algunos comentarios sobre la forma correcta en que debemos interpretar la prueba de inmersión de Hartigans. ¿Me equivoco en alguna parte?
Gracias a todos. Saludos,
TA
Ejemplo de distribución simulada:
Índice de prueba de caída de Hartigan asociado:
P.value de prueba de inmersión de Hartigan asociado:
library(diptest) library(ggplot2) # CONSTANT PARAMETERS sig1 <- log(3) sig2 <- log(3) cpct <- 0.5 N=1000 #CREATING BIMOD DISTRIBUTION bimodalDistFunc <- function (n,cpct, mu1, mu2, sig1, sig2) { y0 <- rlnorm(n,mean=mu1, sd = sig1) y1 <- rlnorm(n,mean=mu2, sd = sig2) flag <- rbinom(n,size=1,prob=cpct) y <- y0*(1 - flag) + y1*flag } #DIP TEST DIP_TEST <- function(bimodalData) { TEST <- dip.test(bimodalData) return(TEST$statistic[[1]]) # return(TEST$p.value[[1]]) to get the p value } DIP_TEST(bimodalData) # SIMULATION exp_mu1 = 1 max_exp_mu2 = 100 intervStep = 100 repPerInt = 10 # single distibutions expMu2Value <- c() bimodalData <- c() mu1 <- log(exp_mu1) mu2 <- log(exp_mu1) bimodalData <- c(bimodalData,log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2))) expMu2Value <- c(expMu2Value,rep(exp_mu1,length(log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2))))) mu1 <- log(exp_mu1) mu2 <- log(max_exp_mu2) bimodalData <- c(bimodalData,log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2))) expMu2Value <- c(expMu2Value,rep(max_exp_mu2,length(log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2))))) mu1 <- log(exp_mu1) mu2 <- log(trunc((max_exp_mu2-exp_mu1)/2+1)) bimodalData <- c(bimodalData,log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2))) expMu2Value <- c(expMu2Value,rep(trunc((max_exp_mu2-exp_mu1)/2+1),length(log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2))))) tableExamples <- data.frame(expMu2Value,bimodalData) tableExamples$expMu2Value <- as.factor(tableExamples$expMu2Value) ExamplePlot <- ggplot(tableExamples)+ geom_histogram(aes(bimodalData),color="white")+ ylab("Count")+ xlab("")+ facet_wrap(~expMu2Value)+ ggtitle("Intensity of bimodularity") # calculation of the dip test index exp_mu2Int = seq(from=exp_mu1,to=max_exp_mu2,length.out=intervStep) expmu2Vec = c() dipStat = c() testDone = c() for(exp_mu2 in exp_mu2Int){ mu1 <- log(exp_mu1) mu2 <- log(exp_mu2) for(rep in 1:repPerInt){ bimodalData <- log(bimodalDistFunc(n=N,cpct,mu1,mu2, sig1,sig2)) diptestone = DIP_TEST(bimodalData) expmu2Vec = c(expmu2Vec,exp_mu2) dipStat = c(dipStat,diptestone) testDone = c(testDone,"diptest") } } table = data.frame(expmu2Vec,dipStat,testDone) IndexPlot <- ggplot(table)+ geom_point(aes(expmu2Vec,dipStat,color=testDone))+ ylab("Index")+ xlab("Intensity of Bimodularity")+ scale_color_discrete(name="Test") ExamplePlot IndexPlot
Comentarios
- Una pregunta muy completa sobre un tema que es arcano por cualquier ' estadístico s estándares . Las primeras preguntas obvias, antes incluso de entrar en la interpretación es, " ¿Por qué necesita esta prueba? ¿Qué información se pretende comunicar? " Puede proporcionar un contexto adicional para las motivaciones que lo han llevado a la cuestión mucho, mucho más posterior de la interpretación de los resultados de la " prueba de inmersión? " En otras palabras, que no sea ' s conveniencia de la programación en R, ¿qué ruta de lógica lo ha llevado a la " prueba de inmersión " en primer lugar?
- Gracias por tu respuesta, Mike. Estoy ' trabajando en un modelo teórico en biología evolutiva y estoy realizando un análisis de sensibilidad. En particular, observo que la variación de algunos parámetros modifica la distribución de una variable de salida de unimodal a bimodal (lo que en realidad es muy interesante). Esa ' es la razón por la que ' busco estadísticas simples para describir la multimodularidad de una distribución. Me permitiría enfocar el análisis de sensibilidad en la multimodularidad.
- Descubrí que la prueba de inmersión se podía calcular fácilmente en R y que podía cuantificar la desviación de una distribución unimodal. Por supuesto, me interesaría mucho cualquier otra estadística que describa la multimodularidad de una distribución.
- Hmmm … ajustar unos pocos polinomios humildes podría equivaler a " Pobre ' s " enfoque para lidiar con la curvilinealidad que ' estás observando y podría implementarse e interpretarse más fácilmente que la prueba de Hartigan ' s. No ' no dice si sus problemas incluyen el manejo de funciones de crecimiento.Por ejemplo, en el desarrollo humano, hay varios " baches " bien conocidos en la trayectoria de crecimiento en distintos puntos del ciclo de vida. Se ha encontrado que los modelos no paramétricos se ajustan mejor y se aproximan mejor a estas no linealidades que los modelos paramétricos.
- Sobre las cuestiones estadísticas: Como se dijo, la prueba de inmersión toma la unimodalidad como referencia. No ' no creo que las desviaciones se puedan interpretar en términos del número de modos solo a partir del valor P. ' he encontrado inmensamente más útil interpretar el número de modos con una combinación de estimación de densidad e interpretación sustantiva.
Respuesta
Sr. Freeman (autor del artículo del que les hablé) me dijo que en realidad solo estaba mirando el valor P de la prueba de inmersión. Esta confusión proviene de su oración:
«Los valores de HDS van de 0 a 1 con valores menores que .05 indican bimodalidad significativa, y valores mayores que .05 pero menores que .10 sugieren bimodalidad con significancia marginal» . Los valores HDS corresponden al valor P y no a las estadísticas de la prueba de inmersión. No estaba claro en el periódico.
Mi análisis es bueno: las estadísticas de la prueba de inmersión aumentan cuando la distribución se desvía de una distribución unimodal.
La prueba de bimodalidad y la prueba de Silverman también se pueden calcular fácilmente en R y funcionan bien.
Comentarios
- Regístrese & y combine sus cuentas. Puede encontrar información sobre cómo hacerlo en la sección Mi cuenta de nuestra centro de ayuda .
Deja una respuesta