2013-10-12 9 views
9

Mein erstes Ziel war es, eine Population von einzelnen Punkten zu zeichnen und dann eine konvexe Hülle zu zeichnen, die 80% dieser Population um die Masse der Bevölkerung zentriert.Wie richtig interpretieren ggplot stat_density2d

Nachdem ich eine Reihe von Ideen ausprobiert hatte, war die beste Lösung, die ich mir ausgedacht hatte, ggplot zu verwenden. Während dies für eine qualitative Analyse gut funktioniert, muss ich immer noch eine 80% Grenze angeben. Ich begann nach einem Weg zu suchen, um die 80. Perzentil Bevölkerung Grenze zu skizzieren, aber ich kann stattdessen mit einer 80% Wahrscheinlichkeitsgrenze arbeiten.

Hier ist, wo ich Hilfe suche. Der bin Parameter für kde2d (verwendet von stat_density2d) ist nicht eindeutig dokumentiert. Wenn ich im folgenden Beispiel bin = 4 setze, bin ich korrekt bei der Interpretation der zentralen (grünen) Region, die eine 25% ige Wahrscheinlichkeitsmasse enthält und die kombinierten gelben, roten und grünen Flächen eine 75% ige Wahrscheinlichkeitsmasse darstellen? Wenn dies der Fall ist, würde dann, wenn der Behälter auf = 5 geändert wird, der Bereich, der dann beschrieben wird, gleich einer 80% igen Wahrscheinlichkeits-Masse sein?

set.seed(1) 
n=100 

df <- data.frame(x=rnorm(n, 0, 1), y=rnorm(n, 0, 1)) 

TestData <- ggplot (data = df) + 
    stat_density2d(aes(x = x, y = y, fill = as.factor(..level..)), 
    bins=4, geom = "polygon",) + 
    geom_point(aes(x = x, y = y)) + 
    scale_fill_manual(values = c("yellow","red","green","royalblue", "black")) 

TestData 

enter image description here

wiederholte ich eine Anzahl von Testfällen und manuell die ausgeschlossenen Punkte gezählt [würde gerne einen Weg finden, um sie auf zu zählen, basierend was ..level .. sie innerhalb enthalten waren] aber gegeben die zufällige Art der Daten (sowohl meine realen Daten und die Testdaten) die Anzahl der Punkte außerhalb der stat_density2d Bereich variiert genug, um zu verlangen, um Hilfe zu bitten.

Zusammenfassend, gibt es eine praktische Möglichkeit, ein Polygon um die zentralen 80% der Population von Punkten im Datenrahmen zu zeichnen? Oder bin ich sicher, dass ich stat_density2d verwenden und Bin auf 5 setzen kann, um eine 80% ige Wahrscheinlichkeitsmasse zu erzeugen?


Ausgezeichnete Antwort von Bryan Hanson Zerstreuung der Fuzzy-Vorstellung, dass ich eine nicht dokumentierte bin Parameter in stat_density2d passieren könnte. Die Ergebnisse sahen nahe bei Werten für bin um 4 bis 6 aus, aber wie er sagte, ist die tatsächliche Funktion unbekannt und daher nicht verwendbar.

Ich habe das HDRegionplot wie in der akzeptierten Antwort von DWin verwendet, um mein Problem zu lösen. Dazu fügte ich einen Schwerpunkt (COGravity) hinzu und zeigte in das Polygon (pnt.in.poly) aus dem SDMTools-Paket, um die Analyse abzuschließen.

library(MASS) 
library(coda) 
library(SDMTools) 
library(emdbook) 
library(ggplot2) 


theme_set(theme_bw(16)) 
set.seed(1) 
n=100 

df <- data.frame(x=rnorm(n, 0, 1), y=rnorm(n, 0, 1)) 

HPDregionplot(mcmc(data.matrix(df)), prob=0.8) 
with(df, points(x,y)) 
ContourLines <- as.data.frame(HPDregionplot(mcmc(data.matrix(df)), prob=0.8)) 
df$inpoly <- pnt.in.poly(df, ContourLines[, c("x", "y")])$pip 

dp <- df[df$inpoly == 1,] 
COG100 <- as.data.frame(t(COGravity(df$x, df$y))) 
COG80 <- as.data.frame(t(COGravity(dp$x, dp$y))) 

TestData <- ggplot (data = df) + 
    stat_density2d(aes(x = x, y = y, fill = as.factor(..level..)), 
    bins=5, geom = "polygon",) + 
    geom_point(aes(x = x, y = y, colour = as.factor(inpoly)), alpha = 1) + 
    geom_point(data=COG100, aes(COGx, COGy),colour="white",size=2, shape = 4) + 
    geom_point(data=COG80, aes(COGx, COGy),colour="green",size=4, shape = 3) + 
    geom_polygon(data = ContourLines, aes(x = x, y = y), color = "blue", fill = NA) + 
    scale_fill_manual(values = c("yellow","red","green","royalblue", "brown", "black", "white", "black", "white","black")) + 
    scale_colour_manual(values = c("red", "black")) 
TestData 
nrow(dp)/nrow(df) # actual number of population members inscribed within the 80% probability polgyon 

enter image description here

Antwort

2

HPDregionplot im Paket: emdbook soll das tun. Es verwendet MASS :: kde2d, normalisiert aber das Ergebnis. Es hat den Nachteil, dass es ein mcmc-Objekt erfordert.

library(MASS) 
library(coda) 
HPDregionplot(mcmc(data.matrix(df)), prob=0.8) 
with(df, points(x,y)) 

enter image description here

2

Okay, mir sagen, ich bin nicht ganz sicher dieser Antwort beginnen zu lassen, und es ist nur eine Teilantwort! Es gibt keinen bin Parameter für MASS::kde2d, der die von stat_density2d verwendete Funktion ist.Betrachtet man die Hilfeseite für kde2d und den Code dafür (einfach durch Eingabe des Funktionsnamens in der Konsole), denke ich, dass der bin Parameter h ist (wie diese Funktionen wissen bin zu h ist jedoch nicht klar). Wenn Sie auf die Hilfeseite klicken, sehen Sie, dass h nicht von MASS:bandwidth.nrd berechnet wird. Die Hilfeseite für diese Funktion sagt dies:

# The function is currently defined as 
function(x) 
{ 
    r <- quantile(x, c(0.25, 0.75)) 
    h <- (r[2] - r[1])/1.34 
    4 * 1.06 * min(sqrt(var(x)), h) * length(x)^(-1/5) 
} 

Auf dieser Grundlage habe ich die Antwort auf Ihre letzte Frage denken („Bin ich sicher ...“) ist auf jeden Fall nicht. r in der oben genannten Funktion ist, was Sie für Ihre Annahme sicher sein müssen, aber es ist deutlich verändert, so dass Sie nicht sicher sind. HTH.

Zusätzlicher Gedanke: Haben Sie Beweise dafür, dass Ihr Code Ihr bins Argument verwendet? Ich frage mich, ob es ignoriert wird. Wenn ja, versuchen Sie h anstelle von bins übergeben und sehen, ob es zuhört.