1 Preperation

1.1 Set Working Directory (WD)

R is always pointed at a directory on your computer. By default,

  • R loads files (Plain text files / Data files etc.) from this directory.

  • R outputs files (Plain text files / Data files etc.) to this directory.

  • Set Working Directory (WD)
    • On Windows
      • in R interface choose : Menu - File - Change dir...
      • Choose a local directory (folder)
        • for example I use C:/Users/Admin/Documents/Konstanz_summer_school/Day2 (actrually Day 3) as my working directory
    • On MacOs
      • in R insterface choose : menu - Misc - Change Working Directory...
        • for example I use ~/Documents/Konstanz_summer_school/Day2 as my working directory
    • Verify your WD
      • type getwd() in R console prompt: after >, then press enter key
      • R should print your current WD

1.2 Download dataset and load it to R

  • The dataset is accessible via the following two links Link 1 or Link 2
  • Click on one of the two Links, a new web page will open then you have two options
    • Option 1: press Windows Shortcut : Control + S or MacOS Shortcut : Command ⌘ + S
    • Option 2: right click your mouse (or touchpad) and choose save as
  • save KBrO3_kinetics.txt to your working directory (in my case /Documents/Konstanz_summer_school/Day2)
  • once this is done, type list.files() in R console prompt: after >, then press enter key
  • R should print all files in your working directory. Check if KBrO3_kinetics.txt is listed
  • if not click on

Copy - Paste in R console prompt: after >, then press enter key

tmpData = readLines("https://antoinechn.github.io/Konstanz-SummerSchool/data/KBrO3_kinetics.txt")
writeLines(tmpData, "KBrO3_kinetics.txt")
rm(tmpData)
  • once this is done, type list.files() in R console prompt: after >, then press enter key
list.files()
  • R should print all files in your working directory. Check if KBrO3_kinetics.txt is listed
  • if still not
    • check your internet connection !!!
    • or contact us.

2 KBrO3 in RPTEC cells

2.1 Load data

We will fist model the in vitro kinetics of KBrO3 in RPTEC cells

# read the data from a file
kinetic.data = read.table("KBrO3_kinetics.txt")
exists("kinetic.data")
#> [1] TRUE

Try this

kinetic.data = read.table("https://trello-attachments.s3.amazonaws.com/59552633a56853b5f6d1f1f9/5a92760183331b90626ce166/f99b2b57b957d442cec69446144ff461/KBrO3_kinetics.txt")
exists("kinetic.data")

Try this

kinetic.data = read.table("https://antoinechn.github.io/Konstanz-SummerSchool/data/KBrO3_kinetics.txt")
exists("kinetic.data")

2.3 Plot kinetic.data

  • Plot the data to get a sense of it
# plot them on two panels 
par(mfrow=c(1,2)) # One row, Two columns
plot(kinetic.data[,1], kinetic.data[,2], las=1, type="b",
     xlab="Time (min)", ylab="KBrO3 (mM)",
     log="", pch=16, cex=1.3, cex.lab=1.2, col="orange")
plot(kinetic.data[,1], kinetic.data[,3], las=1, type="b",
     xlab="Time (min)", ylab="",
     log="", pch=16, cex=1.3, cex.lab=1.2, col="red")

2.4 Kinetic Modelling

Now we need a model for those data. They look like decreasing exponentially (quite usual). It is customary to look at them on a log y-scale.

2.4.1 Plot the data on a log y-scale

par(mfrow=c(1,2))
plot(kinetic.data[,1], kinetic.data[,2], las=1, type="b",
     xlab="Time (min)", ylab="KBrO3 (mM)",
     log="y", pch=16, cex=1.3, cex.lab=1.2, col="orange")
plot(kinetic.data[,1], kinetic.data[,3], las=1, type="b",
     xlab="Time (min)", ylab="",
     log="y", pch=16, cex=1.3, cex.lab=1.2, col="red")

2.4.2 Model

It seems that there is only one line (not two or more lines joined with different slopes). In this case a mono-exponential decrease model will probably do. That model would be:

\[KBrO_3 = C_0 \times e^{-k \cdot \text{time}}\]

Question

  • What values should we give to parameters \(C_0\) and \(k\)?
    • for \(C_0\) ?
    • for \(k\), try different values till you match the low dose data, let’s begin with \(k = 0\) :

\(C_0\) is obvious \(0.375\)

C0 = 0.375
k = 0 # Try different values till you match the low dose data

KBrO3 = C0 * exp(-k * kinetic.data[,1])

par(mfrow=c(1,2))
plot(kinetic.data[,1], kinetic.data[,2], las=1, type="b",
     xlab="Time (min)", ylab="KBrO3 (mM)",
     log="y", pch=16, cex=1.3, cex.lab=1.2, col="orange")
lines(kinetic.data[,1], KBrO3, col="blue")

Exercise

  • repeat the above lines, changing k value, till you find a decent fit
  • report the value found: \[k = ?\]

Solution

  • C0 = 0.375
  • k = 0.00209
C0 = 0.375
k = 0.00209

KBrO3 = C0 * exp(-k * kinetic.data[,1])
par(mfrow=c(1,2))
plot(kinetic.data[,1], kinetic.data[,2], las=1, type="b",
     xlab="Time (min)", ylab="KBrO3 (mM)",
     log="y", pch=16, cex=1.3, cex.lab=1.2, col="orange",
     main = "this seems to be a good fit")
lines(kinetic.data[,1], KBrO3, col="blue")

Exercise

  • Does it work at high dose? (\(k\) might depend on the dose…)
  • how would you check it?

Solution

C0 = 6
k = 0.00209

KBrO3 = C0 * exp(-k * kinetic.data[,1])
par(mfrow=c(1,2))
plot(kinetic.data[,1], kinetic.data[,3], las=1, type="b",
     xlab="Time (min)", ylab="KBrO3 (mM)",
     log="y", pch=16, cex=1.3, cex.lab=1.2, col="red")
lines(kinetic.data[,1], KBrO3, col="blue")

3 KBrO3 - GHS - Kinetics

OK, so now we should have the KBrO3 kinetics about OK. let’s reanalyze the KBrO3 - GSH data, taking the kinetics into account.

3.1 Load data_GSH

# re-read the GSH data
data_GSH = read.table("data_GSH.txt", sep="\t", header=TRUE)
data_GSH
# replot the data
plot(data_GSH$KBrO3_mM, data_GSH$pct_control)

Remember that we have fitted a Hill dose-response model to that data

Hill = function (x, EC50, n, from, diff) {
  return(from + diff * (x^n / (x^n + EC50^n)))
}

Discussion

x was the nominal (initial) concentration of KBrO3 in the medium, but now we know that KBrO3 disappears progressively (probably by reduction to bromide).

  • So, what should be x?
    • should we take the concentration of KBrO3 at 1 hour ?
    • or the average concentration over the 1st hour?

Recall The GSH measurements were made after 1 hour.

# If we take the concentration at 1h we should use
KBrO3 = C0 * exp(-k * 1 * 60) # 1 hour is 60 minutes...
# If we take the average we should use
KBrO3 = (-C0*exp(-k * 60)/k + C0*exp(0)/k) / (60 - 0)

Questions

Does it make any difference? why?

Plot the GSH data as a funtion of KBr03 concentration at 60 min

  • Define a sequence of nominal concentrations
Hill = function (x, EC50, n, from, diff) {
  return(from + diff * (x^n / (x^n + EC50^n)))
}
# define an extended Hill dose-response model as a function

C0 = seq(0, 10, 0.1)
# the actual bromate concentration seen by the cells was:
par(mfrow= c(1,2))
KBrO3 = C0 * exp(-k * 1 * 60) # Model 1
y = Hill(x=KBrO3, EC50=10, n=1, from=100, diff=-100)
plot(C0, y, type="b", main = "KBrO3 = C0 * exp(-k * 1 * 60)", ylim = c(50, 100), las = 1)

KBrO3 = (-C0*exp(-k * 60)/k + C0*exp(0)/k) / (60 - 0)
y = Hill(x=KBrO3, EC50=10, n=1, from=100, diff=-100) # Model 2
plot(C0, y, type="b", main = "KBrO3 = (-C0*exp(-k * 60)/k + C0*exp(0)/k) / (60 - 0)", ylim = c(50, 100), las = 1)

Exercise

Now replot the data and overlay with the model

Remark make sure to define the axes limits before doing overlays…

xlims = c(0, 6) # a vector of min and max
ylims = c(0, 120)
KBrO3.experimental.1h = data_GSH$KBrO3_mM * exp(-k * 1 * 60)
plot(KBrO3.experimental.1h, data_GSH$pct_control, xlim=xlims, ylim=ylims, col="red",
     xlab = "KBrO3 concentration at 1h (mM)", ylab = "GSH (% of control)",
     cex.lab = 1.4, las = 1)
par(new=T) # will overlay the next plot
y = Hill(x = KBrO3, EC50 = 1, n = 1, from = 100, diff = -100)
plot(C0, y, type="l", xlim=xlims, ylim=ylims, xlab="", ylab="",
     xaxt = "n", yaxt = "n", col = "red", lwd = 2.5) # turn off labels

  • now play with EC50 and N to fit the data better

we have coupled PK and PD (dose-response)

```

Back to home