Skip to content

Commit b3d43b6

Browse files
author
Daniel Kaschek
committed
FIX: Update the vignette. Add jak-stat data as data.frame provided by the dMod package.
1 parent 820d838 commit b3d43b6

File tree

8 files changed

+352
-204
lines changed

8 files changed

+352
-204
lines changed

R/data.R

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -113,6 +113,16 @@ res <- function(data, out, err = NULL) {
113113
}
114114

115115

116+
#' Time-course data for the JAK-STAT cell signaling pathway
117+
#'
118+
#' Phosphorylated Epo receptor (pEpoR), phosphorylated STAT in the
119+
#' cytoplasm (tpSTAT) and total STAT (tSTAT) in the cytoplasmhave been
120+
#' measured at times 0, ..., 60.
121+
#'
122+
#' @name jakstat
123+
#' @docType data
124+
#' @keywords data
125+
NULL
116126

117127

118128

R/toolsMirjam.R

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -11,7 +11,7 @@
1111
#' }
1212
#'
1313
#' @export
14-
lsdMod <- function(classlist = c("odemodel", "parfn", "prdfn", "obsfn", "objfn"), envir = .GlobalEnv){
14+
lsdMod <- function(classlist = c("odemodel", "parfn", "prdfn", "obsfn", "objfn", "datalist"), envir = .GlobalEnv){
1515
glist <- as.list(envir)
1616
out <- list()
1717
for (a in classlist) {

data/jakstat.tab

Lines changed: 48 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,48 @@
1+
"time" "name" "value" "sigma" "condition"
2+
2 "tpSTAT" 0.3315 0.05 "Epo"
3+
4 "tpSTAT" 0.8645 0.066 "Epo"
4+
6 "tpSTAT" 0.9635 0.07 "Epo"
5+
8 "tpSTAT" 0.9279 0.065 "Epo"
6+
10 "tpSTAT" 0.8162 0.051 "Epo"
7+
12 "tpSTAT" 0.7553 0.053 "Epo"
8+
14 "tpSTAT" 0.768 0.051 "Epo"
9+
16 "tpSTAT" 0.8416 0.04 "Epo"
10+
18 "tpSTAT" 0.768 0.04 "Epo"
11+
20 "tpSTAT" 0.801 0.048 "Epo"
12+
25 "tpSTAT" 0.7832 0.052 "Epo"
13+
30 "tpSTAT" 0.8086 0.054 "Epo"
14+
40 "tpSTAT" 0.4888 0.055 "Epo"
15+
50 "tpSTAT" 0.2782 0.044 "Epo"
16+
60 "tpSTAT" 0.2553 0.071 "Epo"
17+
0 "tSTAT" 1 0.084 "Epo"
18+
2 "tSTAT" 0.9275 0.046 "Epo"
19+
4 "tSTAT" 0.7923 0.038 "Epo"
20+
6 "tSTAT" 0.7778 0.032 "Epo"
21+
8 "tSTAT" 0.7053 0.033 "Epo"
22+
10 "tSTAT" 0.6522 0.037 "Epo"
23+
12 "tSTAT" 0.5894 0.039 "Epo"
24+
14 "tSTAT" 0.5894 0.04 "Epo"
25+
16 "tSTAT" 0.6377 0.03 "Epo"
26+
18 "tSTAT" 0.6425 0.028 "Epo"
27+
20 "tSTAT" 0.6908 0.03 "Epo"
28+
25 "tSTAT" 0.6908 0.031 "Epo"
29+
30 "tSTAT" 0.7585 0.032 "Epo"
30+
40 "tSTAT" 0.8068 0.04 "Epo"
31+
50 "tSTAT" 0.9275 0.046 "Epo"
32+
60 "tSTAT" 0.971 0.082 "Epo"
33+
0 "pEpoR" 0.01713 0.05 "Epo"
34+
2 "pEpoR" 0.145 0.05 "Epo"
35+
4 "pEpoR" 0.2442 0.05 "Epo"
36+
6 "pEpoR" 0.7659 0.05 "Epo"
37+
8 "pEpoR" 1 0.05 "Epo"
38+
10 "pEpoR" 0.8605 0.05 "Epo"
39+
12 "pEpoR" 0.7829 0.05 "Epo"
40+
14 "pEpoR" 0.5705 0.05 "Epo"
41+
16 "pEpoR" 0.6217 0.05 "Epo"
42+
18 "pEpoR" 0.331 0.05 "Epo"
43+
20 "pEpoR" 0.3388 0.05 "Epo"
44+
25 "pEpoR" 0.3116 0.05 "Epo"
45+
30 "pEpoR" 0.05062 0.05 "Epo"
46+
40 "pEpoR" 0.02504 0.05 "Epo"
47+
50 "pEpoR" 0.01163 0.05 "Epo"
48+
60 "pEpoR" 0 0.05 "Epo"

vignettes/dMod.Rmd

Lines changed: 21 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -24,7 +24,6 @@ knitr::opts_chunk$set(echo = TRUE, fig.width = 10, fig.height = 8, warning = FAL
2424
## Load framework
2525

2626
```{r}
27-
library(deSolve)
2827
library(dMod)
2928
set.seed(2)
3029
```
@@ -33,34 +32,37 @@ set.seed(2)
3332

3433
### Model definition
3534

36-
Within `dMod`, models are typically defined in a tabular or by a sequence of `addReaction` commands. The model for the JAK-STAT pathway is saved in a file `topology.csv` which is imported into R and converted to an equation list.
35+
Within `dMod`, models are typically defined in a tabular or by a sequence of `addReaction` commands. The model for the JAK-STAT pathway is defined reaction by reaction in the following code chunk.
3736

3837
```{r}
3938
# Generate the ODE model
40-
reactions <- as.eqnlist(read.csv("topology.csv"))
41-
print(reactions)
39+
r <- NULL
40+
r <- addReaction(r, "STAT", "pSTAT", "p1*pEpoR*STAT", "STAT phosphoyrl.")
41+
r <- addReaction(r, "2*pSTAT", "pSTATdimer", "p2*pSTAT^2", "pSTAT dimerization")
42+
r <- addReaction(r, "pSTATdimer", "npSTATdimer", "p3*pSTATdimer", "pSTAT dimer import")
43+
r <- addReaction(r, "npSTATdimer", "2*nSTAT", "p4*npSTATdimer", "dimer dissociation")
44+
r <- addReaction(r, "nSTAT", "STAT", "p5*nSTAT", "nuclear export")
45+
46+
print(r)
4247
```
4348

4449
The first reaction contains the phosphorylated Epo-receptor which is not modeled mechanistically but should be replaced by an algebraic function that can describe the activation of the receptor. Therefore, we replace `pEpoR` in the rates:
4550

4651
```{r}
4752
# Parameterize the receptor phosphorylation
4853
receptor <- "pEpoR*((1 - exp(-time*lambda1))*exp(-time*lambda2))^3"
49-
reactions$rates <- replaceSymbols(
54+
r$rates <- replaceSymbols(
5055
what = "pEpoR",
5156
by = receptor,
52-
x = reactions$rates
57+
x = r$rates
5358
)
5459
```
5560

5661
Next, the equation list representing the ODE is translated into the ODE model `model0`. The function calls `as.eqnvec()` to generate the vector of (ordinary differential) equations, translates them into C code, derives sensitivity equations, etc. Integrating sensitivity equations can increase the integration time. Therefore, it is recommended to tell `odemodel()` which of the parameters will take fixed values and will not be estimated based on experimental data.
5762

5863
```{r}
59-
# Define parameters that are not estimated
60-
fixed.zero <- setdiff(reactions$states, "STAT")
6164
# Generate odemodel
62-
model0 <- odemodel(reactions, fixed = fixed.zero,
63-
modelname = "jak-stat", compile = TRUE)
65+
model0 <- odemodel(r, modelname = "jak-stat", compile = TRUE)
6466
```
6567

6668
### Prediction and observation functions
@@ -92,14 +94,14 @@ observables <- eqnvec(
9294
9395
# Define the observation function. Information about states and dynamic parameters
9496
# is contained in reactions
95-
g <- Y(observables, reactions, modelname = "obsfn", compile = TRUE, attach.input = FALSE)
97+
g <- Y(observables, r, modelname = "obsfn", compile = TRUE, attach.input = FALSE)
9698
```
9799

98100
To predict the observables, observation and prediction function need to be concatenated. What mathematically looks like $(g\circ x)(t, p)$ is translated into code using the `*` operator.
99101

100102
```{r, fig.width = 6, fig.height = 2.5}
101103
# Make a prediction of the observables based on random parameter values
102-
parameters <- union(getParameters(x), getParameters(g))
104+
parameters <- getParameters(x, g)
103105
pars <- structure(runif(length(parameters), 0, 1), names = parameters)
104106
times <- seq(0, 10, len = 100)
105107
prediction <- (g*x)(times, pars)
@@ -111,15 +113,14 @@ plot(prediction)
111113
Parameter transformations can for example be used to fix initial values or perform a log-transform. The basic object is again an equation vector containing the equations which express the inner model parameters as functions of other parameters. The names of the other parameters can coincide with the inner parameters.
112114

113115
```{r}
114-
# Start with the identity transformation
115-
innerpars <- union(attr(x, "parameters"), attr(g, "parameters"))
116-
trafo <- as.eqnvec(innerpars, names = innerpars)
116+
innerpars <- getParameters(x, g)
117117
118+
# Start with the identity transformation
119+
trafo <- repar("x~x", x = innerpars)
118120
# Fix some initial values
119-
trafo[fixed.zero] <- "0"
120-
121+
trafo <- repar("x~0", x = c("pSTAT", "pSTATdimer", "npSTATdimer", "nSTAT"), trafo)
121122
# Log-transform
122-
trafo <- replaceSymbols(innerpars, paste0("exp(", innerpars, ")"), trafo)
123+
trafo <- repar("x~exp(x)", x = innerpars, trafo)
123124
124125
# Generate the parameter transformation function
125126
p <- P(trafo, condition = "Epo")
@@ -158,8 +159,8 @@ Model predictions are organized in lists where the elements correspond to differ
158159
This data frame is converted into a data list that can readily be plotted:
159160

160161
```{r, fig.width = 6, fig.height = 2}
161-
datasheet <- read.csv("pnas_data_original.csv")
162-
data <- as.datalist(datasheet, split.by = "condition")
162+
data(jakstat)
163+
data <- as.datalist(jakstat, split.by = "condition")
163164
plot(data)
164165
```
165166

vignettes/dMod.html

Lines changed: 272 additions & 0 deletions
Large diffs are not rendered by default.

vignettes/pnas_data_original.csv

Lines changed: 0 additions & 48 deletions
This file was deleted.

vignettes/trial-1-02-04-2016-154010/mstrust.log

Lines changed: 0 additions & 135 deletions
This file was deleted.
-24.1 KB
Binary file not shown.

0 commit comments

Comments
 (0)