2010-11-18 6 views
9

Ich muss einige Daten über Internet-Sitzungen für eine DSL-Leitung analysieren. Ich wollte sehen, wie die Sitzungsdauern verteilt sind. Ich dachte, ein einfacher Weg, dies zu tun, würde darin bestehen, ein Wahrscheinlichkeitsdichtediagramm der Dauer aller Sitzungen zu erstellen.Getting Probability Dichte von Daten

Ich habe die Daten in R geladen und die density() Funktion verwendet. Also, es war so etwas

plot(density(data$duration), type = "l", col = "blue", main = "Density Plot of Duration", 
    xlab = "duration(h)", ylab = "probability density") 

Ich bin neu in R und diese Art von Analyse. Das war, was ich von Google erfahren habe. Ich bekam eine Verschwörung, aber mir blieben einige Fragen. Ist das die richtige Funktion, um das zu tun, was ich versuche, oder gibt es etwas anderes?

In der Grafik fand ich, dass die Y-Achsenskala von 0 ... 1,5 war. Ich verstehe nicht, wie es 1.5 sein kann, sollte es nicht von 0 ... 1 sein?

Auch möchte ich eine glattere Kurve bekommen. Da der Datensatz sehr groß ist, sind die Zeilen wirklich gezackt. Es wäre schöner, wenn sie geglättet werden, wenn ich das vorstelle. Wie würde ich das machen?

+5

Sie interpretieren die Dichte falsch. Die Dichte von X kann als ein Wert ** gesehen werden, der proportional zu ** der Wahrscheinlichkeit ist, aus der Population eine Zahl zu ziehen, die in der Nähe von X liegt. Jetzt ist definitionsgemäß das Integral der Dichtefunktion gleich 1.Dies bedeutet nicht, dass der maximale Wert der Dichtefunktion 1 sein sollte, er kann leicht größer sein. Für eine F-Verteilung mit df = (1,1) ist der maximale Wert für die Dichte (bei 0) sogar unendlich. –

+0

@Joris Ja ich merke jetzt, dass ich es nicht richtig interpretiert habe. eher vereinfachend nahm ich an, dass es wegen seiner Wahrscheinlichkeitsverteilung weniger als 1 wäre :). – sfactor

Antwort

2

Sie sollten mit dem Bandbreitenparameter (bw) herumspielen, um die Glattheit der Kurve zu ändern. Im Allgemeinen macht R einen guten Job und liefert automatisch eine schöne und glatte Kurve, aber das ist vielleicht nicht der Fall für Ihren spezifischen Datensatz.

Wie für den Anruf Sie verwenden, ja, es ist korrekt, type="l" ist nicht erforderlich, es ist die Standardeinstellung zum Plotten von Dichteobjekten. Die Fläche unter der Kurve (d. H. Das Integral von -Inf bis +Inf Ihrer Dichtefunktion) ist = 1.

Nun, ist eine Dichtekurve die beste Sache in Ihrem Fall? Vielleicht, vielleicht nicht ... es hängt wirklich davon ab, welche Art von Analyse du machen willst. Wahrscheinlich wird die Verwendung von hist ausreichend sein, und vielleicht immer informativer, da Sie bestimmte Bins der Dauer auswählen können (siehe ?hist für weitere Informationen).

+0

danke Ich werde es mir ansehen, aber ich verstehe immer noch nicht, warum die Dichteachse größer als 1 ist. – sfactor

+0

Wie gesagt, es ist die Fläche unter der Kurve (also Summe (dx * y)) = 1 Der tatsächliche Wert der y-Achse variiert in Abhängigkeit von der Bandbreite. Kleinere Bandbreitenwerte erzeugen höhere y-Werte. Versuchen Sie, 'dichte (rnorm (1000), 0.2)' und 'dichte (rnorm (1000), 2)' zu zeichnen, um den Unterschied zu sehen. – nico

+0

Der Hist sieht relativ zur Dichte recht verzerrt aus. ist das wegen der Annahme eines normalen Kerns mit einer poisson-verteilten Variable? –

10

Wie nico sagte, sollten Sie hist überprüfen, aber Sie können auch die zwei von ihnen kombinieren. Dann könnte man stattdessen die Dichte mit lines aufrufen. Beispiel:

duration <- rpois(500, 10) # For duration data I assume Poisson distributed 
hist(duration, 
    probability = TRUE, # In stead of frequency 
    breaks = "FD",  # For more breaks than the default 
    col = "darkslategray4", border = "seashell3") 
lines(density(duration - 0.5), # Add the kernel density estimate (-.5 fix for the bins) 
    col = "firebrick2", lwd = 3) 

Sollten Sie so etwas wie: Histogram of duration

Beachten Sie, dass die Kerndichteschätzung eine Gaußsche Kernel als Standard annimmt. Aber die Bandbreite ist oft der wichtigste Faktor. Wenn Sie density direkt aufrufen meldet es die Standard geschätzte Bandbreite:

> density(duration) 

Call: 
     density.default(x = duration) 

Data: duration (500 obs.);  Bandwidth 'bw' = 0.7752 

     x     y    
Min. : 0.6745 Min. :1.160e-05 
1st Qu.: 7.0872 1st Qu.:1.038e-03 
Median :13.5000 Median :1.932e-02 
Mean :13.5000 Mean :3.895e-02 
3rd Qu.:19.9128 3rd Qu.:7.521e-02 
Max. :26.3255 Max. :1.164e-01 

Hier 0,7752 ist. Überprüfen Sie es auf Ihre Daten und spielen Sie damit herum, wie von nico vorgeschlagen. Vielleicht möchten Sie sich ?bw.nrd ansehen.

+0

sehr gut ~~~~~~~~~~~~~~~~ –

1

Ich wollte dies als Kommentar zur vorherigen Antwort hinzufügen, aber es ist zu groß. Der scheinbare Skew ist auf die Art zurückzuführen, in der die Werte in einem Histogramm zusammengefasst werden. Es ist oft ein Fehler, Histogramme für diskrete Daten zu verwenden. Siehe unten ...

set.seed(1001) 
tmpf <- function() { 
    duration <- rpois(500, 10) # For duration data I assume Poisson distributed 
    hist(duration, 
     probability = TRUE, # In stead of frequency 
     breaks = "FD",  # For more breaks than the default 
     col = "darkslategray4", border = "seashell3", 
     main="",ann=FALSE,axes=FALSE,xlim=c(0,25),ylim=c(0,0.15)) 
    box() 
    lines(density(duration), # Add the kernel density estimate 
     col = "firebrick2", lwd = 3) 
    par(new=TRUE) 
    plot(table(factor(duration,levels=0:25))/length(duration), 
     xlim=c(0,25),ylim=c(0,0.15),col=4,ann=FALSE,axes=FALSE) 
} 

par(mfrow=c(3,3),mar=rep(0,4)) 
replicate(9,tmpf()) 
+0

Ja, das stimmt, die Bins werden immer auf beiden Seiten der ganzen Zahl sein (rechts = WAHR vs. rechts = FALSCH). Ich benutze das meistens nur zur vorherigen Visualisierung von Daten, wenig Schaden dort. Aber es könnte leicht mit einem einfachen -0.5 an die Dichte behoben werden ... – eyjo

+0

@eyjo: das geht davon aus, dass Sie ganzzahlige Pausen verwenden, aber Sie sind nicht darauf beschränkt – nico