1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
|
#' Calculate mean degradation parameters for an mmkin row object
#'
#' @return If random is FALSE (default), a named vector containing mean values
#' of the fitted degradation model parameters. If random is TRUE, a list with
#' fixed and random effects, in the format required by the start argument of
#' nlme for the case of a single grouping variable ds.
#' @param object An mmkin row object containing several fits of the same model to different datasets
#' @param random Should a list with fixed and random effects be returned?
#' @param test_log_parms If TRUE, log parameters are only considered in
#' the mean calculations if their untransformed counterparts (most likely
#' rate constants) pass the t-test for significant difference from zero.
#' @param conf.level Possibility to adjust the required confidence level
#' for parameter that are tested if requested by 'test_log_parms'.
#' @param default_log_parms If set to a numeric value, this is used
#' as a default value for the tested log parameters that failed the
#' t-test.
#' @export
mean_degparms <- function(object, random = FALSE, test_log_parms = FALSE, conf.level = 0.6,
default_log_parms = NA)
{
if (nrow(object) > 1) stop("Only row objects allowed")
parm_mat_trans <- sapply(object, parms, transformed = TRUE)
if (test_log_parms) {
parm_mat_dim <- dim(parm_mat_trans)
parm_mat_dimnames <- dimnames(parm_mat_trans)
log_parm_trans_names <- grep("^log_", rownames(parm_mat_trans), value = TRUE)
log_parm_names <- gsub("^log_", "", log_parm_trans_names)
t_test_back_OK <- matrix(
sapply(object, function(o) {
suppressWarnings(summary(o)$bpar[log_parm_names, "Pr(>t)"] < (1 - conf.level))
}), nrow = length(log_parm_names))
rownames(t_test_back_OK) <- log_parm_trans_names
parm_mat_trans_OK <- parm_mat_trans
for (trans_parm in log_parm_trans_names) {
parm_mat_trans_OK[trans_parm, ] <- ifelse(t_test_back_OK[trans_parm, ],
parm_mat_trans[trans_parm, ], default_log_parms)
}
} else {
parm_mat_trans_OK <- parm_mat_trans
}
mean_degparm_names <- setdiff(rownames(parm_mat_trans), names(object[[1]]$errparms))
degparm_mat_trans <- parm_mat_trans[mean_degparm_names, , drop = FALSE]
degparm_mat_trans_OK <- parm_mat_trans_OK[mean_degparm_names, , drop = FALSE]
fixed <- apply(degparm_mat_trans_OK, 1, mean, na.rm = TRUE)
if (random) {
random <- t(apply(degparm_mat_trans[mean_degparm_names, , drop = FALSE], 2, function(column) column - fixed))
# If we only have one parameter, apply returns a vector so we get a single row
if (nrow(degparm_mat_trans) == 1) random <- t(random)
rownames(random) <- levels(nlme_data(object)$ds)
# For nlmixr we can specify starting values for standard deviations eta, and
# we ignore uncertain parameters if test_log_parms is FALSE
eta <- apply(degparm_mat_trans_OK, 1, stats::sd, na.rm = TRUE)
return(list(fixed = fixed, random = list(ds = random), eta = eta))
} else {
return(fixed)
}
}
|
