typical R-User here, applying a bunch of packages to my data, hoping for a convincing result although I understand only half of what I do at most (and then getting rejected by the reviewers, surprise). But bear with me, I am here to learn.
I have a time-series of monthly precipitation data, for which I'd like to find a suitable distribution. The Cullen and Frey graph from the
fitdistrplus package indicates that a lognormal distribution should outperform a gamma distribution:
and the results of fitting a gamma, Weibull and lognormal distribution seem to say the same:
Goodness-of-fit statistics Weibull lognormal gamma Kolmogorov-Smirnov statistic 0.0263855 0.03245692 0.07115254 Cramer-von Mises statistic 0.1427163 0.11374430 1.13149886 Anderson-Darling statistic 1.0934641 0.62867455 6.73258465 Goodness-of-fit criteria Weibull lognormal gamma Akaike's Information Criterion 7845.461 7842.426 7954.480 Bayesian Information Criterion 7854.940 7851.905 7963.959
This is of particular interest to me, as in my field of research a gamma distribution is often used for precipitation data, without considering that there might be better fitting distributions.
I was quite satisfied with my results until I stumbled upon the GAMLSS package, which offers an all-in-one solution for finding a suiting distribution to your data. So I applied the easy to use `fitDist' function to see what's the outcome:
> fit <- fitDist(df$prcp, k = round(log(length(df$prcp))), type = "realplus", trace = F , try.gamlss = TRUE) There were 11 warnings (use warnings() to see them) > summary(fit) ******************************************************************* Family: c("BCCG", "Box-Cox-Cole-Green") Call: gamlssML(formula = y, family = DIST[i], data = sys.parent()) Fitting method: "nlminb" Coefficient(s): Estimate Std. Error t value Pr(>|t|) eta.mu 45.1741563 0.9560784 47.24943 < 2.22e-16 *** eta.sigma -0.5522153 0.0265535 -20.79636 < 2.22e-16 *** eta.nu 0.4250894 0.0439835 9.66474 < 2.22e-16 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Degrees of Freedom for the fit: 3 Residual Deg. of Freedom 842 Global Deviance: 7827.23 AIC: 7833.23 SBC: 7847.45
First of all, the function offers me a distribution I've never heard of, but that's fair since (you might've already guessed it) I am quite new to this statistical topic it makes sense that I don't know "all" distributions. But, as I investigated the results of the other distributions, it got quite interesting (values are the BIC):
> fit$fits BCCG BCCGo GA GG BCT BCTo BCPE BCPEo WEI2 WEI WEI3 GB2 GIG exGAUS LNO 7848.230 7848.230 7852.426 7853.005 7854.763 7854.763 7854.989 7854.989 7855.461 7855.461 7855.461 7857.338 7859.426 7871.134 7964.480 LOGNO2 LOGNO IG EXP PARETO2o PARETO2 GP IGAMMA 7964.480 7964.480 8260.617 8290.105 8297.106 8297.106 8297.112 8521.495
Right behind the mysterious BCCG-distributions follows the GA (gamma) distribution, quite closely followed by Weibull, and then, lognormal ranked even lower.
So here are my questions:
- What's the reason for these contrary results? Is it due to different fitting methods?
- What can I further do to find the most fitting distribution?
I know in general that all models are wrong and that it probably won't make a huge difference whether I pick the result of
GAMLSS for my further analysis (calculation of SPI). But still, I am wondering about the comparability of research:
- In a publication, I would probably state the package that I used for the fitting and present some goodness-of-fit values (such as the BIC e.g.). But, obviously, there's more to that and the results depend on the package used. But in the publications that I investigated, there's rarely more information offered. Doesn't that hinder the comparability of research?