data=scan("C:/Documents and Settings/rmil013/My Documents/courses/ML/Data/Cauchy.dat")

#Negative log-likelihood for theta 
cauchy.nll=function(theta,x) sum( log(1+(x-theta)^2) )

#Take a quick look at the log-likelihood
thetapts=seq(0,10,0.01)
loglike=rep(0,length(thetapts))
for(i in 1:length(thetapts)) loglike[i]=-cauchy.nll(thetapts[i],data)
plot(thetapts,loglike,type="l",xlab="theta",ylab="Log-likelihood",las=1)
abline(v=theta,lty=2)
title(paste("Cauchy log-likelihood for sample size 100",", theta=",theta))

#Obtain the MLE and the hessian (2nd deriv of the log-likelihood)                  
nlm(cauchy.nll,0,x=data,hessian=T) 

#OUTPUT
$minimum
[1] 147.0483

$estimate
[1] 3.161341

$gradient
[1] -8.091358e-08

$hessian
        [,1]
[1,] 45.9446

$code
[1] 1

$iterations
[1] 5