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
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
|
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/nlme.R
\name{nlme_function}
\alias{nlme_function}
\alias{mean_degparms}
\alias{nlme_data}
\title{Estimation of parameter distributions from mmkin row objects}
\usage{
nlme_function(object)
mean_degparms(object)
nlme_data(object)
}
\arguments{
\item{object}{An mmkin row object containing several fits of the same model to different datasets}
}
\value{
A function that can be used with nlme
A named vector containing mean values of the fitted degradation model parameters
A \code{\link{groupedData}} object
}
\description{
These functions facilitate setting up a nonlinear mixed effects model for
an mmkin row object. An mmkin row object is essentially a list of mkinfit
objects that have been obtained by fitting the same model to a list of
datasets.
}
\examples{
sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120)
m_SFO <- mkinmod(parent = mkinsub("SFO"))
d_SFO_1 <- mkinpredict(m_SFO,
c(k_parent_sink = 0.1),
c(parent = 98), sampling_times)
d_SFO_1_long <- mkin_wide_to_long(d_SFO_1, time = "time")
d_SFO_2 <- mkinpredict(m_SFO,
c(k_parent_sink = 0.05),
c(parent = 102), sampling_times)
d_SFO_2_long <- mkin_wide_to_long(d_SFO_2, time = "time")
d_SFO_3 <- mkinpredict(m_SFO,
c(k_parent_sink = 0.02),
c(parent = 103), sampling_times)
d_SFO_3_long <- mkin_wide_to_long(d_SFO_3, time = "time")
d1 <- add_err(d_SFO_1, function(value) 3, n = 1)
d2 <- add_err(d_SFO_2, function(value) 2, n = 1)
d3 <- add_err(d_SFO_3, function(value) 4, n = 1)
ds <- c(d1 = d1, d2 = d2, d3 = d3)
f <- mmkin("SFO", ds, cores = 1, quiet = TRUE)
mean_dp <- mean_degparms(f)
grouped_data <- nlme_data(f)
nlme_f <- nlme_function(f)
# These assignments are necessary for these objects to be
# visible to nlme and augPred when evaluation is done by
# pkgdown to generated the html docs.
assign("nlme_f", nlme_f, globalenv())
assign("grouped_data", grouped_data, globalenv())
library(nlme)
m_nlme <- nlme(value ~ nlme_f(name, time, parent_0, log_k_parent_sink),
data = grouped_data,
fixed = parent_0 + log_k_parent_sink ~ 1,
random = pdDiag(parent_0 + log_k_parent_sink ~ 1),
start = mean_dp)
summary(m_nlme)
plot(augPred(m_nlme, level = 0:1), layout = c(3, 1))
\dontrun{
# Test on some real data
ds_2 <- lapply(experimental_data_for_UBA_2019[6:10],
function(x) x$data[c("name", "time", "value")])
m_sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
A1 = mkinsub("SFO"), use_of_ff = "min")
m_sfo_sfo_ff <- mkinmod(parent = mkinsub("SFO", "A1"),
A1 = mkinsub("SFO"), use_of_ff = "max")
m_fomc_sfo <- mkinmod(parent = mkinsub("FOMC", "A1"),
A1 = mkinsub("SFO"))
m_dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
A1 = mkinsub("SFO"))
m_sforb_sfo <- mkinmod(parent = mkinsub("SFORB", "A1"),
A1 = mkinsub("SFO"))
f_2 <- mmkin(list("SFO-SFO" = m_sfo_sfo,
"SFO-SFO-ff" = m_sfo_sfo_ff,
"FOMC-SFO" = m_fomc_sfo,
"DFOP-SFO" = m_dfop_sfo,
"SFORB-SFO" = m_sforb_sfo),
ds_2)
grouped_data_2 <- nlme_data(f_2["SFO-SFO", ])
mean_dp_sfo_sfo <- mean_degparms(f_2["SFO-SFO", ])
mean_dp_sfo_sfo_ff <- mean_degparms(f_2["SFO-SFO-ff", ])
mean_dp_fomc_sfo <- mean_degparms(f_2["FOMC-SFO", ])
mean_dp_dfop_sfo <- mean_degparms(f_2["DFOP-SFO", ])
mean_dp_sforb_sfo <- mean_degparms(f_2["SFORB-SFO", ])
nlme_f_sfo_sfo <- nlme_function(f_2["SFO-SFO", ])
nlme_f_sfo_sfo_ff <- nlme_function(f_2["SFO-SFO-ff", ])
nlme_f_fomc_sfo <- nlme_function(f_2["FOMC-SFO", ])
# Allowing for correlations between random effects leads to non-convergence
f_nlme_sfo_sfo <- nlme(value ~ nlme_f_sfo_sfo(name, time,
parent_0, log_k_parent_sink, log_k_parent_A1, log_k_A1_sink),
data = grouped_data_2,
fixed = parent_0 + log_k_parent_sink + log_k_parent_A1 + log_k_A1_sink ~ 1,
random = pdDiag(parent_0 + log_k_parent_sink + log_k_parent_A1 + log_k_A1_sink ~ 1),
start = mean_dp_sfo_sfo)
# The same model fitted with transformed formation fractions does not converge
f_nlme_sfo_sfo_ff <- nlme(value ~ nlme_f_sfo_sfo_ff(name, time,
parent_0, log_k_parent, log_k_A1, f_parent_ilr_1),
data = grouped_data_2,
fixed = parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1,
random = pdDiag(parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1),
start = mean_dp_sfo_sfo_ff)
# It does converge with this version of reduced random effects
f_nlme_sfo_sfo_ff <- nlme(value ~ nlme_f_sfo_sfo_ff(name, time,
parent_0, log_k_parent, log_k_A1, f_parent_ilr_1),
data = grouped_data_2,
fixed = parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1,
random = pdDiag(parent_0 + log_k_parent ~ 1),
start = mean_dp_sfo_sfo_ff)
f_nlme_fomc_sfo <- nlme(value ~ nlme_f_fomc_sfo(name, time,
parent_0, log_alpha, log_beta, log_k_A1, f_parent_ilr_1),
data = grouped_data_2,
fixed = parent_0 + log_alpha + log_beta + log_k_A1 + f_parent_ilr_1 ~ 1,
random = pdDiag(parent_0 + log_alpha + log_beta + log_k_A1 + f_parent_ilr_1 ~ 1),
start = mean_dp_fomc_sfo)
# DFOP-SFO and SFORB-SFO did not converge with full random effects
anova(f_nlme_fomc_sfo, f_nlme_sfo_sfo)
}
}
|
