This function sets up a nonlinear mixed effects model for an mmkin row object for use with the saemix package. 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.
saemix_model(object) saemix_data(object, ...)
Arguments
| object | An mmkin row object containing several fits of the same model to different datasets |
|---|---|
| ... | Further parameters passed to saemix::saemixData |
Value
An saemix::SaemixModel object.
An saemix::SaemixData object.
Examples
ds <- lapply(experimental_data_for_UBA_2019[6:10], function(x) subset(x$data[c("name", "time", "value")])) names(ds) <- paste("Dataset", 6:10) sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"), A1 = mkinsub("SFO"))#>f_mmkin <- mmkin(list("SFO-SFO" = sfo_sfo), ds, quiet = TRUE, cores = 5) # \dontrun{ if (require(saemix)) { m_saemix <- saemix_model(f_mmkin) d_saemix <- saemix_data(f_mmkin) saemix_options <- list(seed = 123456, save = FALSE, save.graphs = FALSE) saemix(m_saemix, d_saemix, saemix_options) }#>#> #>#> #> #> The following SaemixModel object was successfully created: #> #> Nonlinear mixed-effects model #> Model function: Mixed model generated from mmkin object Model type: structural #> function (psi, id, xidep) #> { #> uid <- unique(id) #> res_list <- parallel::mclapply(uid, function(i) { #> transparms_optim <- psi[i, ] #> names(transparms_optim) <- names(degparms_optim) #> odeini_optim <- transparms_optim[odeini_optim_parm_names] #> names(odeini_optim) <- gsub("_0$", "", odeini_optim_parm_names) #> odeini <- c(odeini_optim, odeini_fixed)[names(mkin_model$diffs)] #> ode_transparms_optim_names <- setdiff(names(transparms_optim), #> odeini_optim_parm_names) #> odeparms_optim <- backtransform_odeparms(transparms_optim[ode_transparms_optim_names], #> mkin_model, transform_rates = object[[1]]$transform_rates, #> transform_fractions = object[[1]]$transform_fractions) #> odeparms <- c(odeparms_optim, odeparms_fixed) #> xidep_i <- subset(xidep, id == i) #> if (analytical) { #> out_values <- mkin_model$deg_func(xidep_i, odeini, #> odeparms) #> } #> else { #> i_time <- xidep_i$time #> i_name <- xidep_i$name #> out_wide <- mkinpredict(mkin_model, odeparms = odeparms, #> odeini = odeini, solution_type = object[[1]]$solution_type, #> outtimes = sort(unique(i_time))) #> out_index <- cbind(as.character(i_time), as.character(i_name)) #> out_values <- out_wide[out_index] #> } #> return(out_values) #> }, mc.cores = 15) #> res <- unlist(res_list) #> return(res) #> } #> <bytecode: 0x555559875398> #> <environment: 0x55555973a248> #> Nb of parameters: 4 #> parameter names: parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> distribution: #> Parameter Distribution Estimated #> [1,] parent_0 normal Estimated #> [2,] log_k_parent normal Estimated #> [3,] log_k_A1 normal Estimated #> [4,] f_parent_ilr_1 normal Estimated #> Variance-covariance matrix: #> parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> parent_0 1 0 0 0 #> log_k_parent 0 1 0 0 #> log_k_A1 0 0 1 0 #> f_parent_ilr_1 0 0 0 1 #> Error model: constant , initial values: a.1=1 #> No covariate in the model. #> Initial values #> parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> Pop.CondInit 86.53449 -3.207005 -3.060308 -1.920449 #> #> #> The following SaemixData object was successfully created: #> #> Object of class SaemixData #> longitudinal data for use with the SAEM algorithm #> Dataset ds_saemix #> Structured data: value ~ time + name | ds #> X variable for graphs: time () #> Running main SAEM algorithm #> [1] "Mon May 25 12:48:51 2020" #> .#> .#> .#> .#> #> Minimisation finished #> [1] "Mon May 25 12:56:39 2020"#> Nonlinear mixed-effects model fit by the SAEM algorithm #> ----------------------------------- #> ---- Data ---- #> ----------------------------------- #> Object of class SaemixData #> longitudinal data for use with the SAEM algorithm #> Dataset ds_saemix #> Structured data: value ~ time + name | ds #> X variable for graphs: time () #> Dataset characteristics: #> number of subjects: 5 #> number of observations: 170 #> average/min/max nb obs: 34.00 / 30 / 38 #> First 10 lines of data: #> ds time name value mdv cens occ ytype #> 1 Dataset 6 0 parent 97.2 0 0 1 1 #> 2 Dataset 6 0 parent 96.4 0 0 1 1 #> 3 Dataset 6 3 parent 71.1 0 0 1 1 #> 4 Dataset 6 3 parent 69.2 0 0 1 1 #> 5 Dataset 6 6 parent 58.1 0 0 1 1 #> 6 Dataset 6 6 parent 56.6 0 0 1 1 #> 7 Dataset 6 10 parent 44.4 0 0 1 1 #> 8 Dataset 6 10 parent 43.4 0 0 1 1 #> 9 Dataset 6 20 parent 33.3 0 0 1 1 #> 10 Dataset 6 20 parent 29.2 0 0 1 1 #> ----------------------------------- #> ---- Model ---- #> ----------------------------------- #> Nonlinear mixed-effects model #> Model function: Mixed model generated from mmkin object Model type: structural #> function (psi, id, xidep) #> { #> uid <- unique(id) #> res_list <- parallel::mclapply(uid, function(i) { #> transparms_optim <- psi[i, ] #> names(transparms_optim) <- names(degparms_optim) #> odeini_optim <- transparms_optim[odeini_optim_parm_names] #> names(odeini_optim) <- gsub("_0$", "", odeini_optim_parm_names) #> odeini <- c(odeini_optim, odeini_fixed)[names(mkin_model$diffs)] #> ode_transparms_optim_names <- setdiff(names(transparms_optim), #> odeini_optim_parm_names) #> odeparms_optim <- backtransform_odeparms(transparms_optim[ode_transparms_optim_names], #> mkin_model, transform_rates = object[[1]]$transform_rates, #> transform_fractions = object[[1]]$transform_fractions) #> odeparms <- c(odeparms_optim, odeparms_fixed) #> xidep_i <- subset(xidep, id == i) #> if (analytical) { #> out_values <- mkin_model$deg_func(xidep_i, odeini, #> odeparms) #> } #> else { #> i_time <- xidep_i$time #> i_name <- xidep_i$name #> out_wide <- mkinpredict(mkin_model, odeparms = odeparms, #> odeini = odeini, solution_type = object[[1]]$solution_type, #> outtimes = sort(unique(i_time))) #> out_index <- cbind(as.character(i_time), as.character(i_name)) #> out_values <- out_wide[out_index] #> } #> return(out_values) #> }, mc.cores = 15) #> res <- unlist(res_list) #> return(res) #> } #> <bytecode: 0x555559875398> #> <environment: 0x55555973a248> #> Nb of parameters: 4 #> parameter names: parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> distribution: #> Parameter Distribution Estimated #> [1,] parent_0 normal Estimated #> [2,] log_k_parent normal Estimated #> [3,] log_k_A1 normal Estimated #> [4,] f_parent_ilr_1 normal Estimated #> Variance-covariance matrix: #> parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> parent_0 1 0 0 0 #> log_k_parent 0 1 0 0 #> log_k_A1 0 0 1 0 #> f_parent_ilr_1 0 0 0 1 #> Error model: constant , initial values: a.1=1 #> No covariate in the model. #> Initial values #> parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> Pop.CondInit 86.53449 -3.207005 -3.060308 -1.920449 #> ----------------------------------- #> ---- Key algorithm options ---- #> ----------------------------------- #> Estimation of individual parameters (MAP) #> Estimation of standard errors and linearised log-likelihood #> Estimation of log-likelihood by importance sampling #> Number of iterations: K1=300, K2=100 #> Number of chains: 10 #> Seed: 123456 #> Number of MCMC iterations for IS: 5000 #> Simulations: #> nb of simulated datasets used for npde: 1000 #> nb of simulated datasets used for VPC: 100 #> Input/output #> save the results to a file: FALSE #> save the graphs to files: FALSE #> ---------------------------------------------------- #> ---- Results ---- #> ---------------------------------------------------- #> ----------------- Fixed effects ------------------ #> ---------------------------------------------------- #> Parameter Estimate SE CV(%) #> [1,] parent_0 86.21 1.51 1.7 #> [2,] log_k_parent -3.21 0.59 18.5 #> [3,] log_k_A1 -4.64 0.29 6.3 #> [4,] f_parent_ilr_1 -0.32 0.30 93.2 #> [5,] a.1 4.69 0.27 5.8 #> ---------------------------------------------------- #> ----------- Variance of random effects ----------- #> ---------------------------------------------------- #> Parameter Estimate SE CV(%) #> parent_0 omega2.parent_0 6.07 7.08 117 #> log_k_parent omega2.log_k_parent 1.75 1.11 63 #> log_k_A1 omega2.log_k_A1 0.26 0.26 101 #> f_parent_ilr_1 omega2.f_parent_ilr_1 0.38 0.27 71 #> ---------------------------------------------------- #> ------ Correlation matrix of random effects ------ #> ---------------------------------------------------- #> omega2.parent_0 omega2.log_k_parent omega2.log_k_A1 #> omega2.parent_0 1 0 0 #> omega2.log_k_parent 0 1 0 #> omega2.log_k_A1 0 0 1 #> omega2.f_parent_ilr_1 0 0 0 #> omega2.f_parent_ilr_1 #> omega2.parent_0 0 #> omega2.log_k_parent 0 #> omega2.log_k_A1 0 #> omega2.f_parent_ilr_1 1 #> ---------------------------------------------------- #> --------------- Statistical criteria ------------- #> ---------------------------------------------------- #> Likelihood computed by linearisation #> -2LL= 1064.397 #> AIC = 1082.397 #> BIC = 1078.882 #> #> Likelihood computed by importance sampling #> -2LL= 1063.161 #> AIC = 1081.161 #> BIC = 1077.646 #> ----------------------------------------------------#> Nonlinear mixed-effects model fit by the SAEM algorithm #> ----------------------------------------- #> ---- Data and Model ---- #> ----------------------------------------- #> Data #> Dataset ds_saemix #> Longitudinal data: value ~ time + name | ds #> #> Model: #> Mixed model generated from mmkin object #> 4 parameters: parent_0 log_k_parent log_k_A1 f_parent_ilr_1 #> error model: constant #> No covariate #> #> Key options #> Estimation of individual parameters (MAP) #> Estimation of standard errors and linearised log-likelihood #> Estimation of log-likelihood by importance sampling #> Number of iterations: K1=300, K2=100 #> Number of chains: 10 #> Seed: 123456 #> Number of MCMC iterations for IS: 5000 #> Input/output #> results not saved #> no graphs #> ---------------------------------------------------- #> ---- Results ---- #> Fixed effects #> Parameter Estimate SE CV(%) #> parent_0 86.214 1.506 1.75 #> log_k_parent -3.210 0.593 18.47 #> log_k_A1 -4.643 0.294 6.34 #> f_parent_ilr_1 -0.322 0.300 93.24 #> a.1 4.689 0.270 5.76 #> #> Variance of random effects #> Parameter Estimate SE CV(%) #> omega2.parent_0 6.068 7.078 116.7 #> omega2.log_k_parent 1.752 1.111 63.4 #> omega2.log_k_A1 0.256 0.257 100.5 #> omega2.f_parent_ilr_1 0.385 0.273 70.8 #> #> Statistical criteria #> Likelihood computed by linearisation #> -2LL= 1064.397 #> AIC= 1082.397 #> BIC= 1078.882 #> Likelihood computed by importance sampling #> -2LL= 1063.161 #> AIC= 1081.161 #> BIC= 1077.646# }
