R Under development (unstable) (2023-12-07 r85661 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2023 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(bbmle)
Loading required package: stats4
> 
> funcresp <-
+ structure(list(Initial = as.integer(c(5, 5, 10, 10, 15, 15, 20, 
+ 20, 30, 30, 50, 50, 75, 75, 100, 100)), Killed = as.integer(c(1, 
+ 2, 5, 6, 10, 9, 7, 10, 11, 15, 5, 21, 32, 18, 25, 35))), .Names = c("Initial", 
+ "Killed"), class = "data.frame", row.names = c("1", "2", "3", 
+ "4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", "15", 
+ "16"))
> 
> binomNLL2 = function(p) {
+   a = p[1]
+   h = p[2]
+   ##  cat(a,h,"\n")
+   p = a/(1+a*h*N)
+   -sum(dbinom(k,prob=p,size=N,log=TRUE))
+ }
> 
> N=0; k=0
> parnames(binomNLL2) = c("a","h")
> m2a = mle2(binomNLL2,start=c(a=0.5,h=0.0125),
+       data=with(funcresp,list(N=Initial,k=Killed)))
> p1a = profile(m2a)
There were 50 or more warnings (use warnings() to see the first 50)
> c2a = print(confint(p1a),digits=3)
    2.5 % 97.5 %
a 0.40250 0.6825
h 0.00699 0.0264
> 
> binomNLL2b = function(p,N,k) {
+   a = p[1]
+   h = p[2]
+   ##  cat(a,h,"\n")
+   p = a/(1+a*h*N)
+   -sum(dbinom(k,prob=p,size=N,log=TRUE))
+ }
> parnames(binomNLL2b) = c("a","h")
> m2b = mle2(binomNLL2,start=c(a=0.5,h=0.0125),
+   data=with(funcresp,list(N=Initial,k=Killed)))
> c2b = confint(m2b)
There were 50 or more warnings (use warnings() to see the first 50)
> 
> N=funcresp$Initial; k=funcresp$Killed
> m2c = mle2(binomNLL2,start=c(a=0.5,h=0.0125))
> c2c = confint(m2c)
There were 50 or more warnings (use warnings() to see the first 50)
> print(c2c,digits=3)
    2.5 % 97.5 %
a 0.40250 0.6825
h 0.00699 0.0264
> 
> 
> proc.time()
   user  system elapsed 
   2.46    0.25    2.68