aboutsummaryrefslogtreecommitdiff
path: root/R
diff options
context:
space:
mode:
authorJohannes Ranke <jranke@uni-bremen.de>2020-05-29 16:05:11 +0200
committerJohannes Ranke <jranke@uni-bremen.de>2020-05-29 16:05:11 +0200
commite6f9e9ca89e35e610d9895b979f1351a47451db0 (patch)
treedd9d389c05e35db7a86abd578751199cd2c6a1be /R
parent510436646b1bdd5b8cfab70be29334bd3cc9c828 (diff)
Improve handling of warnings, reorganize tests
Diffstat (limited to 'R')
-rw-r--r--R/mkinfit.R27
-rw-r--r--R/summary.mkinfit.R16
-rw-r--r--R/summary_DFOP_FOCUS_C.txt82
3 files changed, 111 insertions, 14 deletions
diff --git a/R/mkinfit.R b/R/mkinfit.R
index d7b1b7f4..ec2d3412 100644
--- a/R/mkinfit.R
+++ b/R/mkinfit.R
@@ -233,7 +233,7 @@ if(getRversion() >= '2.15.1') utils::globalVariables(c("name", "time", "value"))
#' fit.FOMC = mkinfit("FOMC", FOCUS_2006_D, quiet = TRUE, error_model = "tc")
#' fit.FOMC_SFO <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE,
#' parms.ini = fit.FOMC$bparms.ode, error_model = "tc")
-#'
+#' }
#' @export
mkinfit <- function(mkinmod, observed,
parms.ini = "auto",
@@ -258,6 +258,8 @@ mkinfit <- function(mkinmod, observed,
{
call <- match.call()
+ summary_warnings <- character()
+
# Derive the name used for the model
if (is.character(mkinmod)) mkinmod_name <- mkinmod
else mkinmod_name <- deparse(substitute(mkinmod))
@@ -289,7 +291,9 @@ mkinfit <- function(mkinmod, observed,
# Also remove zero values to avoid instabilities (e.g. of the 'tc' error model)
if (any(observed$value == 0)) {
- warning("Observations with value of zero were removed from the data")
+ zero_warning <- "Observations with value of zero were removed from the data"
+ summary_warnings <- c(summary_warnings, zero_warning)
+ warning(zero_warning)
observed <- subset(observed, value != 0)
}
@@ -848,8 +852,9 @@ mkinfit <- function(mkinmod, observed,
fit$error_model_algorithm <- error_model_algorithm
if (fit$convergence != 0) {
- fit$warning = paste0("Optimisation did not converge:\n", fit$message)
- warning(fit$warning)
+ convergence_warning = paste0("Optimisation did not converge:\n", fit$message)
+ summary_warnings <- c(warnings, convergence_warning)
+ warning(convergence_warning)
} else {
if(!quiet) message("Optimisation successfully terminated.\n")
}
@@ -918,14 +923,22 @@ mkinfit <- function(mkinmod, observed,
fit$errparms <- errparms
fit$df.residual <- n_observed - length(c(degparms, errparms))
+ # Assign the class here so method dispatch works for residuals
+ class(fit) <- c("mkinfit")
+
# Check for normal distribution of residuals
- fit$shapiro.p <- shapiro.test(residuals.mkinfit(fit, standardized = TRUE))$p.value
- if (fit$shapiro.p < 0.05) warning("The p-value for the Shapiro-Wilk test of normality on standardized residuals is < 0.05")
+ fit$shapiro.p <- shapiro.test(residuals(fit, standardized = TRUE))$p.value
+ if (fit$shapiro.p < 0.05) {
+ shapiro_warning <- paste("Shapiro-Wilk test for standardized residuals: p = ", signif(fit$shapiro.p, 3))
+ warning(shapiro_warning)
+ summary_warnings <- c(summary_warnings, shapiro_warning)
+ }
+
+ fit$summary_warnings <- summary_warnings
fit$date <- date()
fit$version <- as.character(utils::packageVersion("mkin"))
fit$Rversion <- paste(R.version$major, R.version$minor, sep=".")
- class(fit) <- c("mkinfit")
return(fit)
}
diff --git a/R/summary.mkinfit.R b/R/summary.mkinfit.R
index 2c291ffd..f9858c32 100644
--- a/R/summary.mkinfit.R
+++ b/R/summary.mkinfit.R
@@ -122,7 +122,7 @@ summary.mkinfit <- function(object, data = TRUE, distimes = TRUE, alpha = 0.05,
date.fit = object$date,
date.summary = date(),
solution_type = object$solution_type,
- warning = object$warning,
+ warnings = object$summary_warnings,
use_of_ff = object$mkinmod$use_of_ff,
error_model_algorithm = object$error_model_algorithm,
df = c(p, rdf),
@@ -190,8 +190,6 @@ print.summary.mkinfit <- function(x, digits = max(3, getOption("digits") - 3), .
cat("Date of fit: ", x$date.fit, "\n")
cat("Date of summary:", x$date.summary, "\n")
- if (!is.null(x$warning)) cat("\n\nWarning:", x$warning, "\n\n")
-
cat("\nEquations:\n")
nice_diffs <- gsub("^(d.*) =", "\\1/dt =", x[["diffs"]])
writeLines(strwrap(nice_diffs, exdent = 11))
@@ -223,16 +221,20 @@ print.summary.mkinfit <- function(x, digits = max(3, getOption("digits") - 3), .
if(length(x$fixed$value) == 0) cat("None\n")
else print(x$fixed)
+ # We used to only have one warning - show this for summarising old objects
+ if (!is.null(x[["warning"]])) cat("\n\nWarning:", x$warning, "\n\n")
+
+ if (length(x$warnings) > 0) {
+ cat("\n\nWarning(s):", "\n")
+ cat(x$warnings, sep = "\n")
+ }
+
if (!is.null(x$AIC)) {
cat("\nResults:\n\n")
print(data.frame(AIC = x$AIC, BIC = x$BIC, logLik = x$logLik,
row.names = " "))
}
- if (!is.null(x$shapiro.p) && x$shapiro.p < 0.05) {
- warning("The p-value for the Shapiro-Wilk test of normality on standardized residuals is < 0.05")
- }
-
cat("\nOptimised, transformed parameters with symmetric confidence intervals:\n")
print(signif(x$par, digits = digits))
diff --git a/R/summary_DFOP_FOCUS_C.txt b/R/summary_DFOP_FOCUS_C.txt
new file mode 100644
index 00000000..ab64a588
--- /dev/null
+++ b/R/summary_DFOP_FOCUS_C.txt
@@ -0,0 +1,82 @@
+mkin version used for fitting: Dummy 0.0 for testing
+R version used for fitting: Dummy R version for testing
+Date of fit: Dummy date for testing
+Date of summary: Dummy date for testing
+
+Equations:
+d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
+ time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time)))
+ * parent
+
+Model predictions using solution type analytical
+
+Fitted using test 0 model solutions performed in test time 0 s
+
+Error model: Constant variance
+
+Error model algorithm: OLS
+
+Starting values for parameters to be optimised:
+ value type
+parent_0 85.10 state
+k1 0.10 deparm
+k2 0.01 deparm
+g 0.50 deparm
+
+Starting values for the transformed parameters actually optimised:
+ value lower upper
+parent_0 85.100000 -Inf Inf
+log_k1 -2.302585 -Inf Inf
+log_k2 -4.605170 -Inf Inf
+g_ilr 0.000000 -Inf Inf
+
+Fixed parameter values:
+None
+
+Results:
+
+ AIC BIC logLik
+ 29.02372 30.00984 -9.511861
+
+Optimised, transformed parameters with symmetric confidence intervals:
+ Estimate Std. Error Lower Upper
+parent_0 85.0000 0.66620 83.1500 86.8500
+log_k1 -0.7775 0.03380 -0.8713 -0.6836
+log_k2 -4.0260 0.13100 -4.3890 -3.6620
+g_ilr 1.2490 0.05811 1.0870 1.4100
+sigma 0.6962 0.16410 0.2406 1.1520
+
+Parameter correlation:
+[1] "Correlation matrix is platform dependent, not tested"
+
+Backtransformed parameters:
+Confidence intervals for internally transformed parameters are asymmetric.
+t-test (unrealistically) based on the assumption of normal distribution
+for estimators of untransformed parameters.
+ Estimate t value Pr(>t) Lower Upper
+parent_0 85.00000 127.600 1.131e-08 83.15000 86.85000
+k1 0.45960 29.580 3.887e-06 0.41840 0.50480
+k2 0.01785 7.636 7.901e-04 0.01241 0.02568
+g 0.85390 83.310 6.221e-08 0.82310 0.88020
+sigma 0.69620 4.243 6.618e-03 0.24060 1.15200
+
+FOCUS Chi2 error levels in percent:
+ err.min n.optim df
+All data 2.661 4 5
+parent 2.661 4 5
+
+Estimated disappearance times:
+ DT50 DT90 DT50_k1 DT50_k2
+parent 1.887 21.25 1.508 38.83
+
+Data:
+ time variable observed predicted residual
+ 0 parent 85.1 85.003 0.09726
+ 1 parent 57.9 58.039 -0.13912
+ 3 parent 29.9 30.054 -0.15351
+ 7 parent 14.6 13.866 0.73388
+ 14 parent 9.7 9.787 -0.08657
+ 28 parent 6.6 7.532 -0.93205
+ 63 parent 4.0 4.033 -0.03269
+ 91 parent 3.9 2.447 1.45348
+ 119 parent 0.6 1.484 -0.88424

Contact - Imprint