Tolkning av Hartigans ' dopptest
On november 18, 2020 by adminJag skulle vilja hitta ett sätt att kvantifiera intensiteten av en del bimodalitet distributioner fick jag empiriskt. Enligt vad jag läser finns det fortfarande en del debatt om sättet att kvantifiera bimodalitet. Jag valde att använda Hartigans ”dip-test som verkar vara det enda tillgängliga på R (originalpapper: http://www.stat.washington.edu/wxs/Stat593-s03/Literature/hartigan85a.pdf Hartigans ”dopptest definieras som: ” Dopptestet mäter multimodalitet i ett prov med den maximala skillnaden över alla provpunkter mellan den empiriska fördelningsfunktionen och den unimodala fördelningsfunktionen som minimerar den maximala skillnaden ”<. / em>.
Jag vill förstå hur jag ska tolka denna statistik innan jag använder den. Jag förväntade mig att dopptestet skulle öka om fördelningen är multimodal (eftersom den definieras som ”den maximala skillnaden från den unimodala fördelningen”). Men : du kan läsa på wikipedia-sidan om multimodal distribution att ”Värden mindre än 0,05 indikerar betydande bimodalitet och värden större än 0,05 men mindre än 0,10 föreslår bimodalitet med marginell betydelse. ”. Ett sådant uttalande kommer från detta papper (fig. 2). Enligt detta dokument är dip-testindexet nära 0 när fördelningen är bimodal. Det förvirrar mig.
För att tolka korrekt Hartigans ”dip test konstruerade jag några distributioner (den ursprungliga koden är från här ) och jag ökade värdet av exp (mu2) (kallas ”Intensitet av bimodularitet” från och med nu – Redigera: Jag skulle ha kallat det ”Intensitet av bimodalitet” ) för att få bimodalitet. I den första grafen kan du se några exempel på distributioner. Sedan uppskattade jag diptestindex (andra graf) och p-värdet (tredje graf) associerade (paket diptest ) till de olika simulerade distributionerna. R-koden som används är i slutet av mitt inlägg.
Vad jag visar här är att dip-testindexet är högt och Pvalue är låg när distributionserna är bimodala. Vilket är i strid med vad du kan läsa på internet.
Jag är ingen expert på statistik, så att jag knappt förstod Hartigans papper. Jag skulle vilja få några kommentarer om hur vi ska tolka Hartigans ”dip-test. Har jag fel någonstans?
Tack alla. Hälsningar,
TA
Exempel på distribution simulerad:
Hartigan ”s dip test index associerad:
Hartigan ”s dip test p.value associerad:
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
Kommentarer
- Mycket noggrann fråga om ett ämne som är arkant av alla statistiker ' s standarder De uppenbara första frågorna innan man ens går in i tolkning är, " Varför behöver du detta test? Vilken information är det avsett att kommunicera? " Kan ge ytterligare kontext för motivationen som har lett dig till den mycket, mycket längre nedströmsfrågan om tolkning av resultaten från " dip test? " Med andra ord, annat än det ' s bekvämlighet med R-programmering, vilken logikväg har lett dig till " dip test " i första hand?
- Tack för att du svarade, Mike. Jag ' arbetar med en teoretisk modell inom evolutionär biologi och jag gör en känslighetsanalys. I synnerhet observerar jag att varierande vissa parametrar ändrar fördelningen av en utgångsvariabel från unimodal till bimodal (vilket faktiskt är väldigt intressant). Det är ' därför jag ' letar efter en enkel statistik för att beskriva multimodulariteten hos en distribution. Det skulle tillåta mig att fokusera känslighetsanalysen på multimodulariteten.
- Jag fick reda på att dopptestet lätt kunde beräknas i R och att det kunde kvantifiera avvikelsen från en unimodal fördelning. Självklart skulle jag vara väldigt intresserad av någon annan statistik som beskriver multimodulariteten hos en distribution.
- Hmmm … passar några ödmjuka polynom kan uppgå till en " fattig man ' s " inställning till att hantera kurvlinjäriteten som du ' observerar och kan lättare distribueras och tolkas än Hartigan ' test. Du ' säger inte om dina frågor inkluderar att hantera tillväxtfunktioner.Till exempel i mänsklig utveckling finns det flera välkända " stötar " i tillväxtbanan vid distinkta punkter i livscykeln. Icke-parametriska modeller har visat sig passa bättre och uppskatta dessa icke-linjärer än parametriska modeller.
- Om de statistiska frågorna: Som sagt tar dip-testet unimodalitet som referens. Jag tror inte ' att avvikelser från det kan tolkas i termer av antalet lägen bara från P-värdet. Jag ' har funnit det oerhört mer användbart att tolka antalet lägen med en kombination av densitetsuppskattning och materiell tolkning.
Svar
Mr. Freeman (författare till papper jag berättade om) berättade för mig att han faktiskt bara tittade på DIP-testets P-värde. Denna förvirring kommer från hans mening:
”HDS-värden sträcker sig från 0 till 1 med värden mindre än .05 som indikerar signifikant bimodalitet och värden större än .05 men mindre än .10 vilket tyder på bimodalitet med marginal betydelse” . HDS-värden motsvarar Pvalue och inte dip-teststatistiken. Det var oklart i tidningen.
Min analys är bra: dip-teststatistiken ökar när fördelningen avviker från en unimodal fördelning.
Bimodalitetstest och Silvermans test kan också beräknas enkelt i R och gör jobbet bra.
Kommentarer
- Registrera & slå samman dina konton. Du kan hitta information om hur du gör det i avsnittet Mitt konto i vårt hjälpcenter .
Lämna ett svar