aboutsummaryrefslogtreecommitdiff
path: root/R/mkinfit.R
diff options
context:
space:
mode:
authorJohannes Ranke <jranke@uni-bremen.de>2016-04-21 17:34:38 +0200
committerJohannes Ranke <jranke@uni-bremen.de>2016-04-21 17:37:54 +0200
commit0b36fb666b84321814eabbc1a25687ee70d789e6 (patch)
treec2cecc577993eb641504b6c2d2f1c82b59ece744 /R/mkinfit.R
parentae5e2bf30393a1a97dcc2da4f8e8fd7be7028117 (diff)
Avoid unnecessary stop in edge case, improve warnings
mkinfit: Do not error out in cases where the fit converges, but the Jacobian for the untransformed model cost can not be estimated. Give a warning instead and return NA for the t-test results. summary.mkinfit: Give a warning message when the covariance matrix can not be obtained. A test has been added to containing such an edge case to check that the warnings are correctly issued and the fit does not terminate.
Diffstat (limited to 'R/mkinfit.R')
-rw-r--r--R/mkinfit.R20
1 files changed, 14 insertions, 6 deletions
diff --git a/R/mkinfit.R b/R/mkinfit.R
index 929c73f0..2302544b 100644
--- a/R/mkinfit.R
+++ b/R/mkinfit.R
@@ -483,9 +483,15 @@ mkinfit <- function(mkinmod, observed,
# Estimate the Hessian for the model cost without parameter transformations
# to make it possible to obtain the usual t-test
# Code ported from FME::modFit
- Jac_notrans <- gradient(function(p, ...) cost_notrans(p)$residuals$res,
- bparms.optim, centered = TRUE)
- fit$hessian_notrans <- 2 * t(Jac_notrans) %*% Jac_notrans
+ Jac_notrans <- try(gradient(function(p, ...) cost_notrans(p)$residuals$res,
+ bparms.optim, centered = TRUE), silent = TRUE)
+ if (inherits(Jac_notrans, "try-error")) {
+ warning("Calculation of the Jacobian failed for the cost function of the untransformed model.\n",
+ "No t-test results will be available")
+ fit$hessian_notrans <- NA
+ } else {
+ fit$hessian_notrans <- 2 * t(Jac_notrans) %*% Jac_notrans
+ }
# Collect initial parameter values in three dataframes
fit$start <- data.frame(value = c(state.ini.optim,
@@ -548,7 +554,7 @@ summary.mkinfit <- function(object, data = TRUE, distimes = TRUE, alpha = 0.05,
covar_notrans <- try(solve(0.5*object$hessian_notrans), silent = TRUE) # unscaled covariance
rdf <- object$df.residual
resvar <- object$ssr / rdf
- if (!is.numeric(covar)) {
+ if (!is.numeric(covar) | is.na(covar[1])) {
covar <- NULL
se <- lci <- uci <- rep(NA, p)
} else {
@@ -558,7 +564,7 @@ summary.mkinfit <- function(object, data = TRUE, distimes = TRUE, alpha = 0.05,
uci <- param + qt(1-alpha/2, rdf) * se
}
- if (!is.numeric(covar_notrans)) {
+ if (!is.numeric(covar_notrans) | is.na(covar_notrans[1])) {
covar_notrans <- NULL
se_notrans <- tval <- pval <- rep(NA, p)
} else {
@@ -690,7 +696,9 @@ print.summary.mkinfit <- function(x, digits = max(3, getOption("digits") - 3), .
rownames(Corr) <- colnames(Corr) <- rownames(x$par)
print(Corr, digits = digits, ...)
} else {
- cat("Could not estimate covariance matrix; singular system:\n")
+ msg <- "Could not estimate covariance matrix; singular system:\n"
+ warning(msg)
+ cat(msg)
}
}

Contact - Imprint