Setup

First, load the libraries needed for the workflow

library(tidyverse)
library(CoRC)

Input

Put in all variables needed to customise the workflow

Put in the path to the model (SBML), as url or local path (as string) and load the model. Use loadModel() for .cps models.

modelPathSBML <- "./Models/Schaber2012SUP_model.xml"
loadSBML(modelPathSBML)
#> # A COPASI model reference:
#> Model name: "Schaber2012SUP"
#> Number of compartments: 1
#> Number of species: 3
#> Number of reactions: 2

Load the data into the variable data.

data <- read.csv("./Data/Schaber2012SUP_data2", sep = ",")

Define experiment(s) for Parameter Estimation. For multiple experiments use list of defineExperiments().

experiment <- defineExperiments(
  data          = data,
  type          = c("time", "dependent", "ignore"),
  mapping       = c(NA, "{Values[TpFit]}", NA),
  weight_method = "mean_square"
)

Set Parameter Estimation Parameters. For more parameters, copy the defineParameterEstimationParameter() function and append the list.

addParameterEstimationParameter(
    ref = "{(T -> Tp; S).k}",
    start_value = 0.1,
    lower_bound = 0,
    upper_bound = 100)


addParameterEstimationParameter(
    ref = "{(Tp -> T).k1}",
    start_value = 0.1,
    lower_bound = 0,
    upper_bound = 100)

Set Parameter Estimation Settings.

setParameterEstimationSettings(
  method = "LevenbergMarquardt"
)

Set Sigma Point Method Settings. Choose alpha, beta, kappa, and the measurement error (var)

alpha <- 0.5
beta <- 2
kappa <- 3
measurement_error <- 0.1

Do you want to save the result as a RDS object? Also give the file, where you want to save it.

save <- TRUE
file <- "docs/SP_Schaber"

Workflow

This should generally not need to be changed, only if the workflow is to be changed more generally.

# Run Sigma Point Method
result <-
  runSigmaPoint(
    alpha             = alpha,
    beta              = beta,
    kappa             = kappa,
    var               = measurement_error,
    experiments       = experiment,
    mean_fit_as_basis = TRUE
  )

# Save Result
if (save){
  saveRDS(result, 
          file = file)
}

Output

This workflow produces the following output:

result$sp_mean
#> (T -> Tp; S).k   (Tp -> T).k1 
#>       2.443987       1.247621
result$sp_cov_matrix
#>                (T -> Tp; S).k (Tp -> T).k1
#> (T -> Tp; S).k       6.830238     2.904945
#> (Tp -> T).k1         2.904945     1.297842