CODEWITHSUNDEEP
r
user67275
user67275

Reputation: 2180

optim() function error in parameter estimation

The following code is to find the 3 parameters of the generalized beta distribution.

beta1 = function(a,b,t) { beta(a+(1/t),b) }
beta2 = function(a,b,t) { beta(a+(2/t),b) }

geb11 = function(a,b,t) { beta2(a,b,t)/beta(a,b) }
geb12 = function(a,b,t) { (beta1(a,b,t)-beta2(a,b,t))/beta(a,b) }
geb22 = function(a,b,t) { 1 + (beta2(a,b,t)-2*beta1(a,b,t))/beta(a,b) }

gbetloglik = function(v) {
   a = v[1]; b = v[2]; t = v[3]
   loglik = n1*log(geb11(a,b,t)) + n3*log(geb12(a,b,t)) + n2*log(geb22(a,b,t))
   return(-loglik)
}

n1 = 127; n2 = 111; n3 = 262
abt = optim(c(5,5,1),gbetloglik,lower=c(0.01,0.01,0.1),method="L-BFGS-B")$par

Then it gives an error as follows with many warnings.

Error in optim(c(2, 8, 1), gbetloglik, lower = c(0.01, 0.01, 0.1), method = "L-BFGS-B") : 
  L-BFGS-B needs finite values of 'fn'

The warning messages are 

1: In beta(a + (2/t), b) : underflow occurred in 'beta'
2: In beta(a, b) : underflow occurred in 'beta'
3: In beta(a + (1/t), b) : underflow occurred in 'beta'
...

n1, n2 and n3 are generated by some random number generation scenarios, and (n1,n2,n3) decides whether an error occurs. (For some (n1,n2,n3) the error does not occur)

I don't understand what the reason of this error is and what "underflow in 'beta' means.

Upvotes: 0

Views: 486

Answers (1)

G. Grothendieck
G. Grothendieck

Reputation: 270195

As the warning message indicates there is underflow when computing beta.

Try this (where we have added the ## line):

gbetloglik = function(v) {
   a = v[1]; b = v[2]; t = v[3]
   if (beta(a, b) == 0) return(1e6) ##
   loglik = n1*log(geb11(a,b,t)) + n3*log(geb12(a,b,t)) + n2*log(geb22(a,b,t))
   return(-loglik)
}

Upvotes: 1

Related Questions