It is indeed a negative value for sigma that causes the issue. You can check this by inserting this line
if(sigma <= 0 ) cat("Negative sigma=",sigma,"\n") after the line mu <- x %*% beta in function llk.mar Negative values for sigma can be avoided with the use of a transformation for sigma, forcing it to be always positive. Make optim use log(sigma) as parameter and transform this to sigma with sigma <- exp(parm[l]) in llk.mar. Like this # define the log likelihood to be maximized over llk.mar <- function(parm,y,x){ # parm is the vector of parameters # the last element is sigma # y is the response # x is the design matrix l <- length(parm) beta <- parm[-l] sigma <- exp(parm[l]) # <=== transform x <- as.matrix(x) mu <- x %*% beta if(sigma <= 0 ) cat("Negative sigma=",sigma,"\n") llk <- sum(dnorm(y, mu, sigma,log=TRUE)) return(llk) } # initial values parm <- c(as.vector(coef(fit)),log(summary(fit)$sigma)) # use log(sigma) as independent parameter Caveat: transformations often help in situations like this but can lead to badly scaled problems and are not a universal remedy. /Berend -- View this message in context: http://r.789695.n4.nabble.com/Help-on-simple-problem-with-optim-tp2533420p2533939.html Sent from the R help mailing list archive at Nabble.com. ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.