Are there any functions in R for fitting distributions? In particular, I would like to fit the weibull and extreme-value distributions when: i) the data are given and the maximum likelihood estimates are required; ii) the sample moments are available but the complete data record is not (so moment estimates are required). Any suggestions will be gratefully received. _____________________________________ Paul S.P. Cowpertwait IIMS, Massey University, Albany, Private Bag 102 904 North Shore Mail Centre Auckland, NZ Tel (+64) (9) 443 9799 ext 9488 http://www.massey.ac.nz/~pscowper _____________________________________ -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
On Wed, 5 Sep 2001, Cowpertwait, Paul wrote:> > Are there any functions in R for fitting distributions? In particular, I > would like to fit the weibull and extreme-value distributions when: > > i) the data are given and the maximum likelihood estimates are required;See survreg in package survival. It's a trivial case (no regressors, no censoring). survreg(y ~ 1, dist="weibull") survreg(log(y) ~ 1, dist="Weibull") since log-Weibulls are Gumbel. Direct computation of the likelihood is easy: see the recent postings here, or V&R's exercises (e.g 8.10).> ii) the sample moments are available but the complete data record is not > (so moment estimates are required).I think you will need to program that yourself. -- Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272860 (secr) Oxford OX1 3TG, UK Fax: +44 1865 272595 -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
This has come up a couple of times recently. MASS4 will have a
general-purpose function for maximum-likelihood fitting and a set of
wrappers for common distributions.
Here's a preview of the R version (rapidly converted this morning).
You can do things like
fitdist(x, dnorm, list(mean = 0, sd = 1))
fitdist(x, dnorm, list(sd = 1), mean=0.5)
the latter fixing the mean.
This is not much tested (in R and much of it had to be changed for R) but
should give a flavour of what is possible. Eventually there will be
a much more comprehensive interface.
fitdist <- function(x, dist, start, ...)
{
myfn <- function(parm, ...) -sum(log(dens(parm, ...)))
dens1 <- function(parm, x, ...) dist(x, parm, ...)
dens2 <- function(parm, x, ...) dist(x, parm[1], parm[2], ...)
dens3 <- function(parm, x, ...) dist(x, parm[1], parm[2], parm[3], ...)
dens4 <- function(parm, x, ...)
dist(x, parm[1], parm[2], parm[3], , parm[4], ...)
if(missing(x) || length(x) == 0 || mode(x) != "numeric")
stop("`x' must be a non-empty numeric vector")
if(missing(dist) || !is.function(dist))
stop("`dist' must be supplied as a function")
if(missing(start) || !is.list(start))
stop("`start' must be a named list")
nm <- names(start)
## reorder arguments to dist
f <- formals(dist)
args <- names(f)
m <- match(nm, args)
if(any(is.na(m)))
stop("`start' specifies names which are not arguments to
`dist'")
formals(dist) <- c(f[c(1, m)], f[-c(1, m)])
## try to be somewhat efficient
switch(length(nm),
"1" = {dens <- dens1},
"2" = {dens <- dens2},
"3" = {dens <- dens3},
"4" = {dens <- dens3},
stop("only handle 1:4 variable parms"))
res <- optim(start, myfn, x=x, ...)
if(res$convergence > 0) stop("optimization failed")
res$par
}
--
Brian D. Ripley, ripley at stats.ox.ac.uk
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272860 (secr)
Oxford OX1 3TG, UK Fax: +44 1865 272595
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at
stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
Hola! I think the original posters problem was to fit a weibull. this could be done with: library(gnlm)> x <- rweibull(100, 2, 3) > y <- rep(1,100)> fit.dist(x,y,"Weibull")Weibull distribution, n = 100 mean variance alpha.hat mu.hat 2.438586 1.840581 1.890002 6.782797 -log likelihood AIC -296.4451 -294.4451 yi ni pi.hat pi.tilde like.comp resid 1 1.8744688 1 0.01 0.300583890 -3.403142 -5.300156 . . . Kjetil Halvorsen Prof Brian Ripley wrote:> > On Thu, 6 Sep 2001, kjetil halvorsen wrote: > > > Hola! > > > > I think some of this is avaliable as fit.dist in Jim Lindsays > > package gnlm. > > Not according to his help page: > > \name{fit.dist} > \title{Fit Probability Distributions to Frequency Data} > > Not the same problem. > > -- > Brian D. Ripley, ripley at stats.ox.ac.uk > Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ > University of Oxford, Tel: +44 1865 272861 (self) > 1 South Parks Road, +44 1865 272860 (secr) > Oxford OX1 3TG, UK Fax: +44 1865 272595 > > -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- > r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html > Send "info", "help", or "[un]subscribe" > (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch > _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._