# fitting con fitdistrplus
install.packages("fitdistrplus")
library(fitdistrplus)

set.seed(1)
n <- 100
x <- rgamma(n = n, shape = 2, rate = 0.5)

plotdist(x)

# max verosimiglianza
fit.mle <- fitdist(data = x, distr = "gamma", method = "mle")
summary(fit.mle)
confint(fit.mle, level = 0.99)
plot(fit.mle)
fit.mle$estimate

# metodo dei momenti
fit.mme <- fitdist(data = x, distr = "gamma", method = "mme")
summary(fit.mme)
plot(fit.mme)

# metodo dei quantili
fit.qme <- fitdist(data = x, distr = "gamma", method = "qme", probs = c(0.5, 0.75))
summary(fit.qme)
plot(fit.qme)

# confronto weibull, lognormale, gamma
fit.mle.w <- fitdist(data = x, distr = "weibull", method = "mle")
fit.mle.l <- fitdist(data = x, distr = "lnorm", method = "mle")

# analisi e confronto
modelli <- list(fit.mle, fit.mle.w, fit.mle.l)
leg <- c("gamma", "weibull", "lognormale")
par(mfrow = c(2, 2))

denscomp(ft = modelli, legendtext = leg)
cdfcomp(ft = modelli, legendtext = leg)
qqcomp(ft = modelli, legendtext = leg)
ppcomp(ft = modelli, legendtext = leg)


gofstat(f = modelli)