0

I am new in R and I am having difficulty in figuring out how to fit Log Pearson Type III using method = mme to the data.

YearThai <- c(2000,2001,2002,2003,2004,2005,2006,2007,2008,2009,2010,2011,2012,2013,2014,2015,2016,2017,2018)
FlowThai <- c(2009.6,2071.8,3534.5,1813.6,1674.17,1892,4020,2572.4,2303.2,2299.8,3242,3720,1939.9,2221.6,1127,706.8,2109.7,2571.4,880.8)
DataThai <- as.data.frame(cbind(YearThai,FlowThai))

LFlowThai <- log(FlowThai)
m <- mean(LFlowThai)
v <- var(LFlowThai)
s <- sd(LFlowThai)
g <- e1071::skewness(LFlowThai, type=1)

n <- length(LFlowThai)
g <- g*(sqrt(n*(n-1))/(n-2))*(1+8.5/n)

my.shape <- 4/(g)^2
my.scale <- abs(g)/(g*(s/sqrt(my.shape)))
my.location <- m-(my.shape*my.scale)

my.param <- list(shape=my.shape, scale=my.scale, location=my.location)

dPIII<-function(x, shape, location, scale) PearsonDS::dpearsonIII(x, shape, location, scale, log=FALSE)
pPIII<-function(q, shape, location, scale) PearsonDS::ppearsonIII(q, shape, location, scale, lower.tail = TRUE, log.p = FALSE)
qPIII<-function(p, shape, location, scale) PearsonDS::qpearsonIII(p, shape, location, scale, lower.tail = TRUE, log.p = FALSE)

fitdistrplus::fitdist(LFlowThai, distr = "PIII", method = "mme", order = c(1,2,3), start = my.param)

The error I get is this:

Error in mmedist(data, distname, start = arg_startfix$start.arg, fix.arg = arg_startfix$fix.arg,  : 
the empirical moment function must be defined

There is no problem when using method = mle, but our research requires the use of method of moments for the parameter estimation.

Thank you in advance to those who can help.

JOHN WILLIAM VEDEJA
  • 11
  • 2

1 Answers1

1

I'm having the exact same problem as you. What i've found googling about it is that you need to provide the empirical moments function and the theoretical moments function, for Pearson III distribution, when using fitdist from fitdistrplus.

If you use the fasstr package, you'll see that it fits Pearson type III to your data by matching moments. That's because these package provides, to fitdist, the "memp" and "mPearsonIII" functions needed to do so.

memp function:

function(x, order) mean(x^order)

mPearson III function:

function (order, shape, location, scale) 
{
  if (order == 1) 
    return(location + shape * scale)
  if (order == 2) 
    return(scale * scale * shape)
  if (order == 3) 
    return(2/sqrt(shape) * sign(scale))
}

Using compute_annual_frequencies() from fasstr gives you an object named Freq_Fitting with the estimated parameters for Pearson type III using MOM (object class is fitdistrplus). Sadly, is has some issues that i'm trying to fix (quantiles and scale issues when plotting).

Skytled
  • 11
  • 3