Some NaN issues on standard errors of GLLVMs

[rikiherliansyah]

I am currently working on my project for Generalized Linear Latent Variable Models for Poisson-normal model. I can send you the codes if u want to see. The estimated parameters seems correct but when computing the standard errors some parameters yield NaN standard errors (negative element of inv(H)). I have tried to fix this problem but seems nothing wrong. Is there anyone here know what's going on? Thank you.
Warning: Hessian of fixed effects was not positive definite. Maximum gradient component: 109266.6 Warning message: In sqrt(diag(object$cov.fixed)) : NaNs produced

[bbolker]

We probably need a fully reproducible example. Are you using glmmTMB or your own TMB code?

[rikiherliansyah]

I am creating my codes as the following. The model is Poisson with log link function and random effects, u is from standard normal. g(mu) = eta(i, j) = beta(j) + beta(k,j) X(i,k) + lambda(l,j) u(i,l) where i is observational unit, j number of related responses, k number of related explanatory variables, and l is number of latent variables. The estimations seems correct but something wrong with standard errors. I am not sure why. Is that because i define the model incorrectly or any over specified parameters? It would be that great if someone can figure out or ever faced same problem before. Many thanks.

include

include

template Type objective_function::operator() () { DATA_MATRIX(y); DATA_MATRIX(x); PARAMETER_VECTOR(a); PARAMETER_MATRIX(b); PARAMETER_VECTOR(lambda); PARAMETER_MATRIX(u); DATA_INTEGER(num_lv);

int n = y.rows(); int p = y.cols(); //To create lambda as matrix upper triangle matrix newlam(num_lv,p); for (int j=0; j<p; j++){ for (int i=0; i<num_lv; i++){ if (j < i){ newlam(i,j) = 0; } else{ newlam(i,j) = lambda(j); if (i > 0){ newlam(i,j) = lambda(i+j+ip-(i(i-1))/2-2*i); } } } } //std::cout <<newlam<<std::endl;

//multiplication between latent variables and coefficients lambda matrix lam = u*newlam;

//eta function b0 + xb + ulambda matrix eta(n,p); eta = x*b + lam; for (int j=0; j<p; j++){ for (int i=0; i<n; i++){ eta(i,j) = a(j) + eta(i,j); } } //std::cout <<lam<<std::endl; //std::cout <<eta<<std::endl; Type nll = 0.0;

//latent variable model from N(0,1) for (int j=0; j<u.cols();j++){ for (int i=0; i<n; i++) { nll -= dnorm(u(i,j),Type(0),Type(1),true); } } //likelihood poisson model with log link function lambda = (exp(eta)) for (int j=0; j<p;j++){ for (int i=0; i<n; i++) { nll -= dpois(y(i,j), exp(eta(i,j)), true); } } return nll; }

In R library(TMB) compile("linreg.cpp") dyn.load(dynlib("linreg")) library(mvabund) data(antTraits) X <- as.matrix(antTraits$env ) TR <- antTraits$trait y <- as.matrix(antTraits$abund) num.lv = 2 nl = (ncol(birdsy)num.lv)-(num.lv(num.lv-1))/2 lambda = rep(1,nl) u = matrix(1,nrow(birdsy),num.lv) parameters = list(a = rep(0,ncol(birdsy)), b = matrix(0,ncol(birdsx),ncol(birdsy)), lambda = lambda, u = u)

data = list(y = birdsy, x = birdsx, num_lv = num.lv) obj <- MakeADFun(data, parameters, random = "u", DLL = "linreg") system.time(opt <- do.call("optim", obj)) system.time(rep <- sdreport(obj))