`panel6pt1.fn` <-
function(){
# R script to fit closed population model with heterogeneous detection
# probabilities ("Model Mh") to the horned-lizard data. This is 
# the logit-normal version of Model Mh.

expit<-function(x){
1/(1+exp(-x))
}


nx<-c(34, 16, 10, 4, 2, 2,0,0,0,0,0,0,0,0)
nind<-sum(nx)
J<-14


Mhlik<-function(parms){
mu<-parms[1]
sigma<-exp(parms[2])
n0<-exp(parms[3])

il<-rep(NA,J+1)
for(k in 0:J){
il[k+1]<-integrate( 
function(x){
 dbinom(k,J,expit(x))*dnorm(x,mu,sigma) 
},
lower=-Inf,upper=Inf)$value
}

-1*(    lgamma(n0+nind+1) - lgamma(n0+1) + sum(c(n0,nx)*log(il)))
}

print(tmp<-nlm(Mhlik,c(-1,0,log(10) ),hessian=TRUE))
print(sqrt(diag(solve(tmp$hessian))))



lik0<-function(parms){
p<-  expit(parms[1])
n0<- exp(parms[2])
N<-nind + n0
cpvec<- dbinom(0:J,J,p)
-1*(lgamma(N+1) - lgamma(n0+1) + sum(c(n0,nx)*log(cpvec) ))

}

cat("Model M0",fill=TRUE)
 print(nlm(lik0,c(0,2),hessian=TRUE))


}