dat <- read.table("bh.dat", header=TRUE)
par <- list(loga=0, logb=0, logsigma=0)

nll<-function(theta){
  par <- relist(theta, skeleton=par)  
  attach(c(par,dat))
  sigma <- exp(logsigma)
  pred <- loga+log(ssb)-log(1+exp(logb)*ssb)
  nll <- -sum(dnorm(logR,pred,sigma,TRUE))
  detach(c(par,dat))
  return(nll)
}

opt <- nlminb(unlist(par), nll)
print(opt$par)