[R] output hessian matrix in constrOptim

[Spencer Graves]

	  The documentation for "constrOptim" says it has a "..." argument, 
which can contain "Other arguments passed to 'optim', which will pass 
them to 'f' and 'grad' if it does not used them."  I routinely get an 
estimated "hessian" from 'optim' just by specifying the argument 
'hessian=TRUE'.  This led me to try the following

# an example from 'constrOptim' that works:
      fQP <- function(b) {-sum(c(0,5,0)*b)+0.5*sum(b*b)}
      Amat       <- matrix(c(-4,-3,0,2,1,0,0,-2,1),3,3)
      bvec       <- c(-8,2,0)
      constrOptim(c(2,-1,-1), fQP, NULL, ui=t(Amat),ci=bvec)

# The following generated an error:
constrOptim(c(2,-1,-1), fQP, NULL, ui=t(Amat),ci=bvec,
                     hessian=TRUE)
Error in f(theta, ...) : unused argument(s) (hessian ...)

	  If it were my problem, I might first try 'fit <- constrOptim(...)'. 
Then if the optimum was NOT at a constraint, I'd use the outpt from 
'constrOptim' as an input to "optim" with 'hessian=TRUE'.  If it was at 
a constraint, I'd have to think very carefully about what I wanted in 
terms of a "hessian", because that does not seem obvious to me at the 
moment.  Then I'd program something to get that.  If I wanted more than 
that, I'd list 'constrOptim' and / or 'optim', and walk through them 
line by line using "debug" until I figured out what I wanted to do next.

	  If you'd like more help from this listserve, please submit another 
post -- but first, PLEASE do read the posting guide! 
http://www.R-project.org/posting-guide.html.

	  hope this helps,
	  spencer graves

liu wrote:

> Hi, 
> 
> Is there any way to get the hessian matrix from the "constrOptim"  
function without supplying gradient function?   Thanks.
> 
> 
> 
> 		
> ---------------------------------
> Bring words and photos together (easily) with
> 
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html