- if (is.null(alpha))
- out <- nlm(function(p) dev(p),
- p = c(0.1, 1, rep(mean(x), ncol(Wstart))),
- hessian = TRUE)
- else
- out <- nlm(function(p) dev(c(alpha, p)),
- p = c(1, rep(mean(x), ncol(Wstart))),
- hessian = TRUE)
- para <- cbind(out$estimate, sqrt(diag(solve(out$hessian))))
+
+ out <- if (is.null(alpha))
+ nlm(function(p) dev(p),
+ p = c(0.1, 1, rep(mean(x), ncol(Wstart))), hessian = TRUE)
+ else nlm(function(p) dev(c(alpha, p)),
+ p = c(1, rep(mean(x), ncol(Wstart))), hessian = TRUE)
+
+ ## if alpha is estimated it may be that the Hessian matrix has the
+ ## corresponding column and row filled with 0, making solve() fail
+ se <- if (is.null(alpha) && all(out$hessian[1, ] == 0))
+ c(NA, sqrt(diag(solve(out$hessian[-1, -1]))))
+ else sqrt(diag(solve(out$hessian)))
+
+ para <- cbind(out$estimate, se)