Re: [R] rzinb (VGAM) and dnbinom in optim

From: Katharine Mullen <kate_at_few.vu.nl>
Date: Fri, 18 Apr 2008 18:21:33 +0200 (CEST)

I'm not going to look into what's happening in-depth but it looks like the bounds for your parameters need to be set with care; the below (with slight re-def. of your residual function) gives results that seem to match vglm approximately, with 1e-10 standing in for a bound of 0 and 1-e-10 standing in for a bound of 1.

library(VGAM)
phi <- 0.09
size <- 0.7
munb <- 72
x <- rzinb(200, phi=phi, size=size, munb=munb)

fit <- vglm(x~1, zinegbinomial,trace=TRUE)

fncZiNB <- function(xx, x){

   phi <- xx[1]
   size <- xx[2]
   munb <- xx[3]
   -sum(dzinb(x=x,phi=phi, size=size, munb=munb,log=TRUE)) }

 solut <- optim(fn=fncZiNB,par=c(0.09, 72, 0.7),x=x,method =
"L-BFGS-B",  lower=c(1e-10,1e-10,1e-10),upper=c(1-1e-10,Inf,Inf),
control=list(trace=TRUE),hessian=TRUE)

##optim
solut$par
##vglm
Coef(fit)

On Fri, 18 Apr 2008, Tim Howard wrote:

> Dear R-help gurus (and T.Yee, the VGAM maintainer) -
> I've been banging my head against the keyboard for too long now, hopefully someone can pick up on the errors of my ways...
>
> I am trying to use optim to fit a zero-inflated negative binomial distribution. No matter what I try I can't get optim to recognize my initial parameters. I think the problem is that dnbinom allows either 'prob' or 'mu' and I am trying to pass it 'mu'. I've tried a wrapper function but it still hangs. An example:
>
> #set up dummy data,load VGAM, which is where I am getting the zero-inflated neg binom
> #I get the same problem without VGAM, using dnbinom in the wrapper instead.
> library(VGAM)
> phi <- 0.09
> size <- 0.7
> munb <- 72
> x <- rzinb(200, phi=phi, size=size, munb=munb)
>
> #VGAM can fit using vglm... but I'd like to try optim
> > fit <- vglm(x~1, zinegbinomial,trace=TRUE)
> VGLM linear loop 1 : loglikelihood = -964.1915
> VGLM linear loop 2 : loglikelihood = -964.1392
> VGLM linear loop 3 : loglikelihood = -964.1392
> > Coef(fit)
> phi munb k
> 0.1200317 62.8882874 0.8179183
>
> > #build my wrapper function for dzinb
> > fncZiNB <- function(phi, size, munb){
> + -sum(dzinb(x=x,phi=phi, size=size, munb=munb,log=TRUE))}
> >
> > #test the wrapper, seems to work.
> > fncZiNB(phi=0.09, size=0.7, munb=72)
> [1] 969.1435
>
> #try it in optim
> > solut <- optim(fn=fncZiNB,par=c(phi=0.09, munb=72, size=0.7),method = "L-BFGS-B", hessian=TRUE)
> Error in dnbinom(x = x, size = size, prob = prob, log = log.arg) :
> argument "size" is missing, with no default
>
> I have the same problem with dnbinom.
> I'm sure I'm doing something obvious.... any suggestions?
> Session info below. Thanks in advance.
> Tim Howard
> New York Natural Heritage Program
>
> > sessionInfo()
> R version 2.6.2 (2008-02-08)
> i386-pc-mingw32
>
> locale:
> LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
>
> attached base packages:
> [1] stats4 splines stats graphics grDevices utils datasets methods
> [9] base
>
> other attached packages:
> [1] VGAM_0.7-6 randomForest_4.5-23
>
> ______________________________________________
> R-help_at_r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Fri 18 Apr 2008 - 16:26:30 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Sat 19 Apr 2008 - 16:30:30 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive