[R] Best Distribution

[Rolf Turner]

On 01/10/14 08:50, eliza botto wrote:
> Dear useRs,
> I have this following data
>> dput(Prec)
> c(42.2, 45.2, 46, 48, 54, 54.1, 59.4, 61, 62.2, 63.5, 65.024, 71.9, 73.4, 76.6, 76.708, 77.5, 77.724, 78, 81.3, 84.7, 84.836, 85.09, 88.2, 91.4, 94, 95.8, 96, 97.3, 101, 101, 101.5, 102.3, 102.87, 108.7, 109.5, 110.5, 110.7, 112, 114.3, 118.11, 121.412, 128.1, 131, 140, 142, 143.3, 151.4, 153.7, 189.4, 214.3)
> I want to fit gumbel and log-normal distribution on it on the same window to see which distribution fits it the best way.
> Thankyou very much in advance.


If I understand you correctly the attached script should do what you want.

cheers,

Rolf Turner

-- 
Rolf Turner
Technical Editor ANZJS
-------------- next part --------------
#
# Script botto.txt
#

require(MASS)
require(FAdist)

x <- c(42.2, 45.2, 46, 48, 54, 54.1, 59.4, 61, 62.2, 63.5, 65.024,
       71.9, 73.4, 76.6, 76.708, 77.5, 77.724, 78, 81.3, 84.7, 84.836,
       85.09, 88.2, 91.4, 94, 95.8, 96, 97.3, 101, 101, 101.5, 102.3,
       102.87, 108.7, 109.5, 110.5, 110.7, 112, 114.3, 118.11, 121.412,
       128.1, 131, 140, 142, 143.3, 151.4, 153.7, 189.4, 214.3)
dx <- density(x)
s <- sd(x)
xbar <- mean(x)
beta.hat <- 6*s/pi
mu.hat <- xbar - 0.5772*beta.hat
fit.g <- suppressWarnings(fitdistr(x,dgumbel,start=list(location=mu.hat,scale=beta.hat)))
fit.ln <- fitdistr(x,"lognormal")
cg <- fit.g$estimate
cln <- fit.ln$estimate
plot(dx,xlab="x",ylab="density",main="Comparison of fits",ylim=c(0,0.015))
plot(function(x){dgumbel(x,location=cg[1],scale=cg[2])},from=0,to=250,add=TRUE,col="red")
plot(function(x){dlnorm(x,meanlog=cln[1],sdlog=cln[2])},from=0,to=250,add=TRUE,col="blue")
legend("topright",lty=1,col=c("black","red","blue"),
       legend=c("non-parametric density","fitted Gumbel density",
                "fitted log normal density"))