Kommentare
- r ist eindeutig nicht normal. Seine Verteilung ist recht schief (viele kleine Werte nahe 0, wenige große Werte). Sie sehen dies, indem Sie " hist (r) " eingeben.
- Tipp! Verwenden Sie
dput(r)
, um eine Zeichenfolge zu generieren, die leicht kopiert werden kann. ' n ' kann eingefügt werden. Jetzt müssen wir die Daten manuell in r eingeben … - @RasmusB å å Danke, ich habe nach diesem Befehl gesucht 🙂 Ich habe die Frage bearbeitet.
- Sie ' zeichnen keine Sequenz von Datenwerten, um die zu sehen Verteilung. Stellt $ r $ bereits eine Häufigkeit oder Wahrscheinlichkeitsdichte von Datenwerten dar?
- Ich habe diese Frage in verschiedenen Kontexten mehrmals beantwortet. Eine
R
-Lösung für eine diskrete Variable wie IhreIndex
wird unter stats.stackexchange angezeigt. com / a / 43004/919 ; EineR
-Lösung für eine kontinuierliche Variable befindet sich unter stats.stackexchange.com/questions/70153/… ; Eine Excel-Lösung finden Sie unter stats.stackexchange.com/a/11563/919 .
Antwort
Es gibt einen Unterschied zwischen dem Anpassen einer Gaußschen Verteilung und dem Anpassen einer Gaußschen Dichtekurve . Was normalmixEM
macht das erstere. Was Sie wollen, ist (ich denke) das letztere.
Das Anpassen einer Verteilung ist ungefähr das, was Sie tun würden, wenn Sie ein Histogramm erstellen würden von Ihren Daten und versuchte zu sehen, welche Art von Form es hatte. Stattdessen zeichnen Sie einfach eine Kurve. Diese Kurve hat zufällig einen Buckel in der Mitte, wie Sie ihn durch Zeichnen einer Gaußschen Dichtefunktion erhalten.
Um das zu erhalten, was Sie wollen, Sie Sie können optim
verwenden, um die Kurve an Ihre Daten anzupassen. Der folgende Code verwendet nichtlineare kleinste Quadrate, um die drei Parameter zu finden, die die am besten passende Gauß-Kurve ergeben: m
ist der Gaußsche Mittelwert, s
ist die Standardabweichung und k
ist ein beliebiger Skalierungsparameter (seit Die Gaußsche Dichte muss auf 1 integriert werden, während Ihre Daten nicht „t) sind.
x <- seq_along(r) f <- function(par) { m <- par[1] sd <- par[2] k <- par[3] rhat <- k * exp(-0.5 * ((x - m)/sd)^2) sum((r - rhat)^2) } optim(c(15, 2, 1), f, method="BFGS", control=list(reltol=1e-9))
Kommentare
- Mit dieser Lösung herumgespielt, aber du hast mich geschlagen 🙂 Beim Herumspielen habe ich festgestellt, dass die anfänglichen Startwerte für
optim
sehr wichtig sind Stellen Sie bei Verwendung dieser Methode sicher, dass die Anpassung grafisch überprüft wird.
An swer
Ich schlage vor, nichtlineare kleinste Quadrate für diese Analyse zu verwenden.
# First present the data in a data-frame tab <- data.frame(x=seq_along(r), r=r) #Apply function nls (res <- nls( r ~ k*exp(-1/2*(x-mu)^2/sigma^2), start=c(mu=15,sigma=5,k=1) , data = tab))
Und aus der Ausgabe konnte ich die folgende angepasste „Gaußsche Kurve“ erhalten:
v <- summary(res)$parameters[,"Estimate"] plot(r~x, data=tab) plot(function(x) v[3]*exp(-1/2*(x-v[1])^2/v[2]^2),col=2,add=T,xlim=range(tab$x) )
Die Anpassung ist nicht erstaunlich … Wäre die Funktion „ta $ x \ mapsto \ sin (x) / x $“ nicht a besseres Modell?
Kommentare
- Danke. Ich erhalte die verbleibende Quadratsumme: 0,01997. Ich denke, ich bekomme genau das gleiche mit der Lösung von Hong Ooi oben. Ist der Algo der gleiche? Wie zeichne ich auch das Ergebnis von nls?
- Ja, die Algorithmen sind in dem Sinne dieselben, dass sie, wenn sie funktionieren (', in einigen nicht stecken bleiben lokales Minimum) geben sie die gleiche Antwort. Abhängig vom Wert für
method=
können sie genau gleich sein. - Ich habe zwei Zeilen hinzugefügt, um das Diagramm zu generieren.