2009-09-14 8 views
6

Ich habe eine Dichte Objekt wie folgt erstellt DD:aus einer Dichte Objekt mit R

x1 <- rnorm(1000) 
x2 <- rnorm(1000, 3, 2) 
x <- rbind(x1, x2) 
dd <- density(x) 
plot(dd) 

die eben diese nicht-Gauß-Verteilung erzeugt:

alt text http://www.cerebralmastication.com/wp-content/uploads/2009/09/nongaus.png

würde ich Letztendlich mag es Zufallsabweichungen von dieser Verteilung zu bekommen, ähnlich wie rnorm von einer normalen Verteilung abweicht.

Die Art, wie ich dies zu knacken versuchen ist die CDF meiner Kernel zu bekommen und dann bekommt es mir die Veränderlichen zu sagen, wenn ich ihm eine kumulative Wahrscheinlichkeit (inverse CDF) übergebe. Auf diese Weise kann ich einen Vektor von uniformen Zufallsvariablen in Draws aus der Dichte drehen.

Es scheint so, als ob das, was ich versuche, etwas Grundlegendes sein sollte, was andere vor mir getan haben. Gibt es dafür einen einfachen Weg oder eine einfache Funktion? Ich hasse es, das Rad neu zu erfinden.

FWIW Ich fand this R Help article, aber ich kann nicht was sie tun, und die endgültige Ausgabe scheint nicht zu produzieren, was ich danach bin. Aber es könnte ein Schritt auf dem Weg sein, den ich einfach nicht verstehe.

Ich habe überlegt, nur mit einem Johnson distribution from the suppdists package gehen, aber Johnson wird mir nicht den schönen bimodalen Buckel geben, die meine Daten haben.

+0

mehr eine Statistik Frage als Programmierung ... –

+0

weiß, dass ich die Statistiken. Ich möchte die Statistikmethode in einer bestimmten Sprache implementieren. Das ist Programmieren. –

Antwort

2

Dies ist nur eine Mischung von Normalen. Also warum nicht etwas wie:

rmnorm <- function(n,mean, sd,prob) { 
    nmix <- length(mean) 
    if (length(sd)!=nmix) stop("lengths should be the same.") 
    y <- sample(1:nmix,n,prob=prob, replace=TRUE) 
    mean.mix <- mean[y] 
    sd.mix <- sd[y] 
    rnorm(n,mean.mix,sd.mix) 
} 
plot(density(rmnorm(10000,mean=c(0,3), sd=c(1,2), prob=c(.5,.5)))) 

Dies sollte in Ordnung sein, wenn alles, was Sie brauchen, Proben aus dieser Mischung Verteilung sind.

+0

Ich mag die Idee! Aber mein Beispiel ist zur Veranschaulichung zu stark vereinfacht. In Wirklichkeit kenne ich die beiden Modi nicht und es könnte nur einen Modus und einen langen Schwanz (d. H. Leptokurtose) haben. Aber ich mag dein Beispiel. Ich hätte das nicht so knapp programmieren können. BTW, ich glaube, fehlen ac in: Plot (Dichte (rmnorm (10000, Mittelwert = c (0,3), SD = c (1,2), prob = c (0,5, 0,5)))) –

+0

@JD Long: Danke für die Tippfehler. –

+1

Deshalb wollen Sie Hadleys Antwort - Resampling ist es. Denken Sie daran, dass Ihr Dichtediagramm nur eine Schätzung ist, die auch von Ihrem Glättungsparameter abhängt. –

9

Alternative Ansatz:

sample(x, n, replace = TRUE) 
+1

bootstrapping ftw! –

+0

yeah, ich dachte schon darüber nach. Wenn ich + eine Auslosung von einer normalen Probe mache, sollte ich am Ende meine Schwänze verdicken, genauso wie der Kern, richtig? Unter der Annahme, dass ich meine Norm auf die gleiche Weise wie die Kerning-Methode paramieren würde. –

+2

Ja, fügen Sie normale rvs mit Null Mittelwert und SD = Bandbreite aus Dichte Schätzung: Beispiel (x, n, ersetzen = TRUE) + rnorm (n, 0, SD = 0.4214) Simulieren, wie in Silverman 1986 Buch über diskutiert Dichteschätzung. –

Verwandte Themen