# [R] using mle2 for multinomial model optimization

Benedikt Gehr benedikt.gehr at ieu.uzh.ch
Fri Feb 12 14:24:02 CET 2010

```Hi there
I'm trying to find the mle fo a multinomial model ->*L(N,h,S¦x)*. There
is only *N* I want to estimate, which is used in the number of successes
for the last cell probability. These successes are given by:
p^(N-x1-x2-...xi)
All the other parameters (i.e. h and S) I know from somewhere else.
Here is what I've tried to do so far for a imaginary data set:
#######################################################

cohort<-c(50,54,50)     #imaginary harvest numbers of a cohort harvested
over 3 years

l<-length(cohort)+1    #nr of cell probabilities
h<-c(0.2,0.3,1)         #harvest rates for the 3 diferent years
S<-c(0.9,0.8)         #survival rates

mfun <- function(d) {
S<-S            #survival rate
h<-h            #harvest rate
l<-length(d)
p<-rep(0,l)         #Cell probabilities
p<-h

prod<-(1-h)*S
for (i in 2:(l-1)){
p[i]<-prod*h[i]
prod<-prod*(1-h[i])*S[i]

}
p[l]<-1-sum(p[-l])    #last cell probability
-dmultinom(exp(d),prob=p,log=TRUE)    #exp(d)->backtransform the estimates
}

c<-c(cohort,100)            #untransformed initial values for the
optimization,->100=N-x1-x2-x3
nvec<-c(rep("x",l-1),"N")     #names for the initial vector
nrs<-c(1:(l-1),1)            #nrs for the initial vector
svec = log(c)             #transformation of the data to avoid
constraints (x>0)
names(svec) <- paste(nvec,nrs,sep="")
parnames(mfun) <- names(svec)

m1 = mle2(mfun,start=svec,
vecpar=TRUE,fixed=svec[1:(l-1)]) #only estimate
"N-x1-x2-x3",everything else is fixed

coef(m1)
############################################################################

The function "mfun" is working, however the mle2 won't work. How can I
fix this?

Thanks so much for your help

Beni

```