R : Copyright 2001, The R Development Core Team
Version 1.3.1  (2001-08-31)

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.

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

> library(MASS)
> data(Boston)
> crim <- as.integer(Boston$crim > quantile(Boston$crim,0.8))
> med.value <- Boston$medv
> logitreg <- function(beta,x,y)
+ {   eta <- exp(beta[1] + beta[2]*x)
+     p <- eta/(1+eta)
+     return(-sum(ifelse(y,log(p),log(1-p)))) }
> mle <- optim(c(0,0),logitreg,x=med.value,y=crim,hessian=TRUE)
> se <- sqrt(diag(solve(mle$hessian)))
> cbind(lower=mle$par - 1.96*se, upper=mle$par + 1.96*se)
          lower      upper
[1,]  2.1449125  4.0137381
[2,] -0.2896639 -0.1842742
> a <- seq(-10,15,l=50)
> b <- seq(-1,0.2,l=50)
> L <- matrix(NA,50,50)
> for(i in 1:50)
+     for(j in 1:50)
+         L[i,j] <- -logitreg(c(a[i],b[j]), med.value, crim)
> par(mfrow=c(1,2))
> contour(a,b,L)
> persp(a,b,L)
>