Jest różnica między dopasowaniem rozkładu gaussowskiego a dopasowaniem krzywej gęstości gaussa. Co normalmixEM robi to pierwsze. To, czego chcesz, to (chyba) drugie.
Dopasowanie dystrybucji to z grubsza to, co zrobisz, jeśli utworzysz histogram Twoich danych i próbowałeś sprawdzić, jaki miał kształt. Zamiast tego po prostu wykreślasz krzywą. Ta krzywa ma garb na środku, tak jak to, co otrzymujesz, wykreślając funkcję gęstości gaussa.
Aby uzyskać to, czego chcesz, może użyć czegoś takiego jak optim, aby dopasować krzywą do danych. Poniższy kod użyje nieliniowych najmniejszych kwadratów, aby znaleźć trzy parametry dające najlepiej dopasowaną krzywą gaussowską: to średnia gaussowska, s to odchylenie standardowe, a k to dowolny parametr skalowania (ponieważ gęstość gaussa jest ograniczona do całkowania do 1, podczas gdy dane nie są „t).
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))
Komentarze
Proponuję użyć do tej analizy nieliniowych metodą najmniejszych kwadratów .
# 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))
Na podstawie danych wyjściowych udało mi się uzyskać następującą dopasowaną „krzywą Gaussa”:
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) )
Dopasowanie nie jest niesamowite … Czy funkcja ta $ x \ mapsto \ sin (x) / x $ nie byłaby lepszy model?
Komentarze