Error in optim(par = Tm1, fn = .logLgam, x = x) : object 'Tm1' not found

user_8454718 0
Viewed 43

I am trying to run the MSClaio2008 function of the nsRFA library in R. Following the documentation, I have successfully run the examples given on the web site. However, when I created my own data set to run with the MSClaio2008 function as seen below:

library(nsRFA)
data <- c(91,84,42,30,66,95,65,12,27,61,31,101,48,52,53,55,80,23,87,46,50,33,75,88,54,17,57,39,10,89,59,24,11,43,13,93,105,28,104,18,103,62,22,58,15,34,74,51,97,44,99,76,14,16,109,92,110,40,78,83,60,49,96,36,32,81,68,20,56,25,63,47,37,29,100,71,106,41,90,70,85,38,19,108,73,102,26,82,98,45,77,35,94,79,86,72,107,21,67,69,64)
MSC <- MSClaio2008(data)

it throws this error:

Error in optim(par = Tm1, fn = .logLgam, x = x) : object 'Tm1' not found

I have searched for topics on this platform, I have found this but that did not solve my problem.

statistics r
1 Answers

This looks like a bug in the package, but it is related to the fact that your data are weird. This function is supposedly fitting distributions to observations of hydrological extremes: you have given it an integer sequence from 10 to 110.

As it turns out, when computing the fit for some of the possible distributions (P3 and GEV, I think), the function internally computes the skewness of the data and has a test for what to do if it the skewness is positive or negative. It doesn't consider the possibility that the skewness would be exactly 0!

Inside nsRFA:::ML_estimation:

skx <- sum((x - mean(x))^3)/(length(x) * sd(x)^3)
if (skx > 0) {
    Tm1 <- min(x) - 1
    ...
} else if (skx < 0) {
    Tm1 <- max(x) + 1
    ...
}

If the skewness is exactly zero, then the variable Tm1 doesn't get set and the next step that relies on Tm1 breaks.

If you perturb the data in some way so that the skewness is not exactly zero, this should work, e.g. MSClaio2008(c(8,data)) or MSClaio2008(data+rnorm(length(data), sd=0.001)) ...

It's a weird case, but it would be worth contacting the maintainer (maintainer("nsRFA")) to let them know.

user_1902770
Related