[R] Constraint maximum (likelihood) using nlm

[Benjamin Dickgiesser]

Hi,

I'm trying to find the maximum (likelihood) of a function. Therefore,
I'm trying to minimize the negative likelihood function:

# params: vector containing values of mu and sigma
# params[1] - mu, params[2]- sigma
# dat: matrix of data pairs y_i and s_i
# dat[,1] - column of y_i , dat[,2] column of s_i
negll <- function(params,dat,constant=0)
{
	for(i in 1:length(dat[,1]))
	{
		llsum <- log( params[2]^2 + dat[i,2]^2) +
		(( dat[i,1] - params[1])^2/ (params[2]^2 + dat[i,2]^2))
	}
	ll <- -0.5 * llsum  + constant
	return(-ll)
}

Using (find data attached):
data.osl <- read.table("osl.dat",header=TRUE)
data.matrix <- as.matrix(data.frame(data.osl$de,data.osl$se))
nlm(negll,c(0.75,0.5),dat=data.matrix,iterlim=200)

I get estimates for mu and sigma of:
 3.629998e+00 -4.975368e-07

However, sigma obviously has to be >= 0.

Therefore I am trying to transform sigma:

negll.trans <- function(params,...)
{
	params[2] <- log(params[2])
	negll(params,...)

}

where
nlm(negll.trans,c(0.75,0.5),dat=data.matrix)
give me the estimates
3.63 1.00
i.e. mu = 3.63 and sigma = exp(1)

I am not confident that this is the correct answer looking at the graphs:
par(mfrow=c(2,1))
plot(osl$de,osl$se,xlab="equivalent dose",ylab="standard errors",
     main="Equivalent Dose vs. Standard Errors",xlim=c(0,4))
hist(osl$de,xlab="equivalent dose",
     main="Histogram of Equivalent Dose",xlim=c(0,4))

I want to stick with using nlm to minimize the function but can't find
an error in what I am doing.

I'd appreciate your help!

Thank you
Ben