R optim() function
The issue is almost certainly in the objective function i.e., diagH, since Nelder Mead doesn't use any matrix operations such as Choleski. I think you probably need to adjust the objective function to catch singularities (non-positive definite cases). I do notice that you have two identical parameters, so perhaps the matrix is something with identical rows or columns. You may be able to get more insight by calling optimr() from the optimx package (I'm one of the authors) and using a different but similar optimizer such as nmkb which is all in R so the trace is a bit easier. It also is helpful to try the diagH function with the parameters that seem to be giving trouble in a direct call from R. JN
On 2020-10-28 7:03 p.m., Baker, Kieran via R-help wrote:
Hi R-Help,
I am using R to do functional outlier detection (using PCA to reduce to 2 dimensions - the functional boxplot methodology used in the Rainbow package), and using Hscv.diag function to calculate the bandwidth matrix where this line of code is run:
result <- optim(diag(Hstart), scv.mat.temp, method = "Nelder-Mead", control = list(trace = as.numeric(verbose)))
Within the optim function, there is a call to an external C function:
.External2(C_optim, par, fn1, gr1, method, con, lower, upper)
Where
Par = (0.339, 0.339),
fn1 = function (diagH)
{
H <- diag(diagH) %*% diag(diagH)
if (default.bflag(d = d, n = n))
scvm <- scv.mat(x.star, H, Gs, binned = binned, verbose = FALSE,
bin.par = bin.par.star, deriv.order = r)
else scvm <- scv.mat(x.star, H, Gs, binned = binned, verbose = FALSE,
deriv.order = r)
return(scvm)
},
gr1 is Null,
method = "Nelder-Mead"
con is a list of 18 items (trace=0, fnscale=1, parscale=(1,1), ...),
lower and upper are -Inf and Inf respectively.
When this is run, it returns an error
Error in chol.default(H) :
the leading minor of order 1 is not positive definite
Has anyone any experience with this error, or know where I can locate the code for the C_optim function to continue trying to find the issue and a work around?
Thanks for any help,
Kieran
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see 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.