Fortolkning af Hartigans ' dip test
On november 18, 2020 by adminJeg vil gerne finde en måde at kvantificere intensiteten af nogle af bimodaliteten distributioner fik jeg empirisk. Fra det jeg læste, er der stadig en del debat om måden at kvantificere bimodalitet på. Jeg valgte at bruge Hartigans “dip test, som synes at være den eneste tilgængelige på R (originalt papir: http://www.stat.washington.edu/wxs/Stat593-s03/Literature/hartigan85a.pdf Hartigans “dip test defineres som: ” Dip test måler multimodalitet i en prøve med den maksimale forskel over alle prøvepunkter mellem den empiriske fordelingsfunktion og den unimodale fordelingsfunktion, der minimerer den maksimale forskel “<. / em>.
Jeg vil gerne forstå fuldstændigt, hvordan jeg skal fortolke denne statistik, før jeg bruger den. Jeg forventede, at dip-testen ville øges, hvis fordelingen er multimodal (da den er defineret som “den maksimale forskel fra den unimodale fordeling”). Men : du kan læse på wikipedia-siden om multimodal fordeling, at “Værdier mindre end 0,05 indikerer signifikant bimodalitet og værdier større end 0,05, men mindre end 0,10 antyder bimodalitet med marginal betydning. “. En sådan erklæring kommer fra dette papir (fig. 2). Ifølge dette papir er dip testindekset tæt på 0, når fordelingen er bimodal. Det forvirrer mig.
For at fortolke korrekt Hartigans “dip test konstruerede jeg nogle distributioner (den originale kode er fra her ), og jeg steg værdien af exp (mu2) (kaldet “Intensitet af bimodularitet” fra nu af – Rediger: Jeg skulle have kaldt det “Intensitet af bimodalitet” ) for at få bimodalitet. I den første graf kan du se nogle eksempler på distributioner. Derefter estimerede jeg diptestindekset (anden graf) og p-værdien (tredje graf), der var tilknyttet (pakke diptest ) til de forskellige simulerede distributioner. Den anvendte R-kode er i slutningen af mit indlæg.
Det, jeg viser her, er, at dip-testindekset er højt, og Pvalue er lav, når distributionserne er bimodale. Hvilket er i modstrid med hvad du kan læse på internettet.
Jeg er ingen ekspert i statistik, så jeg næppe forstod Hartigans papir. Jeg vil gerne have nogle kommentarer til den rigtige måde, hvorpå vi kan fortolke Hartigans “dip-test. Har jeg forkert et eller andet sted?
Tak til jer alle. Hilsen,
TA
Eksempel på distribution simuleret:
Hartigan “s dip test indeks tilknyttet:
Hartigan “s dip test p.value associeret:
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
- Meget grundig spørgsmålsarbejde om et emne, der er banalt af enhver statistiker ' s standarder De åbenlyse første spørgsmål, inden man overhovedet kommer ind i fortolkning er, " Hvorfor har du brug for denne test? Hvilken information er det beregnet til at kommunikere? " Kan give yderligere kontekst til de motivationer, der har ført dig til det meget, meget længere nedstrøms spørgsmål om fortolkning af resultaterne fra " dip test? " Med andre ord, bortset fra den ' s bekvemmelighed med R-programmering, hvilken logiksti har ført dig til " dip test " i første omgang?
- Tak for dit svar, Mike. Jeg ' Jeg arbejder på en teoretisk model i evolutionær biologi, og jeg udfører en følsomhedsanalyse. Især bemærker jeg, at varierende nogle parametre ændrer fordelingen af en outputvariabel fra unimodal til bimodal (hvilket faktisk er meget interessant). At ' er hvorfor jeg ' leder efter en simpel statistik til at beskrive multimodulariteten af en distribution. Det ville give mig mulighed for at fokusere følsomhedsanalysen på multimodulariteten.
- Jeg fandt ud af, at dip-testen let kunne beregnes i R, og at den kunne kvantificere afvigelsen fra en unimodal fordeling. Selvfølgelig ville jeg være meget interesseret i andre statistikker, der beskriver multimodulariteten af en distribution.
- Hmmm … montering af et par ydmyge polynomer kan udgøre en " stakkels mand ' s " tilgang til at håndtere den krumhed, du ' observerer og kan lettere implementeres og fortolkes end Hartigan ' s test. Du siger ikke ' om dine problemer inkluderer håndtering af vækstfunktioner.For eksempel er der i menneskelig udvikling flere velkendte " bump " i vækstbanen på forskellige punkter i livscyklussen. Ikke-parametriske modeller har vist sig bedre at tilpasse og tilnærme disse ikke-lineariteter end parametriske modeller.
- Om de statistiske spørgsmål: Som sagt tager dip-testen unimodalitet som reference. Jeg tror ikke ' at afvigelser fra det kan fortolkes i form af antallet af tilstande lige fra P-værdien. Jeg ' har fundet det umådeligt mere nyttigt at fortolke antallet af tilstande med en kombination af densitetsestimering og materiel fortolkning.
Svar
Mr. Freeman (forfatter af papiret , jeg fortalte dig om) fortalte mig, at han faktisk kun kiggede på dip-testens Pværdi. Denne forvirring kommer fra hans sætning:
“HDS-værdier spænder fra 0 til 1 med værdier mindre end .05, der indikerer signifikant bimodalitet, og værdier større end .05, men mindre end .10, hvilket antyder bimodalitet med marginal betydning” . HDS-værdier svarer til Pvalue og ikke dip-teststatistikken. Det var uklart i avisen.
Min analyse er god: statistikken for dip-test øges, når fordelingen afviger fra en unimodal fordeling.
Bimodalitetstest og Silvermans test kan også let beregnes i R og gør jobbet godt.
Kommentarer
- Registrer & flet dine konti. Du kan finde oplysninger om, hvordan du gør dette i afsnittet Min konto i vores Hjælp .
Skriv et svar