Cheminformatics/Other/Test ICT 7 Sept.Rmd

181 lines
4.4 KiB
Plaintext

---
title: "Test 7 Sept"
author: "Jonathan Herrewijnen"
date: "September 7, 2018"
output: html_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
## R Markdown
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see <http://rmarkdown.rstudio.com>.
When you click the **Knit** button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
```{r cars}
summary(cars)
```
## Including Plots
You can also embed plots, for example:
```{r pressure, echo=FALSE}
plot(pressure)
```
Note that the `echo = FALSE` parameter was added to the code chunk to prevent printing of the R code that generated the plot.
//Question 1
```{r}
library(Deriv)
t=5.49116
Bp1=exp(3.2)+((5.7*t)/(t+13.0))*exp(-0.026*t)
Bp1
#T VERVANGEN VOOR X
f1=function(x) {exp(3.2)+((5.7*x)/(x+13.0))*exp(-0.026*x)}
#Question 2: use uniroot, or should I? t=5.49116
?uniroot
root=uniroot(function(x), f1(x), c(-1, 1), x=26)$root
root$root
#Werkt niet echt: formule ombouwen
#Bp1={exp(3.2)+((5.7*x)/(x+13.0))*exp(-0.026*x)}=26
#Question 3: Create derivative and put to 0
x=0
df.f1 = Deriv(f1, "x")
df.f1
#FOR QUESTION 3: put in x=0! Run df.f1 for exact formula -> t=0.4384615
function (x)
{
.e1 <- 13 + 0
(5.7 - 0 * (0.1482 + 5.7/.e1)) * exp(-(0.026 * 0))/.e1
}
#Question 4: IF PREVIOUS ANSWER IS CORRECT
#t=100000000 #for Q5
Bp1=exp(3.2)+((5.7*t)/(t+13.0))*exp(-0.026*t)
Bp1 #[1] 24.7164 ??
#Quick Curve
curve(f1(x), -10, 40000, col="red", ylab="f1(x)")
##QUESTION 3+4 LATER?
#Question 5: Just put in a huge number...
#Question 6: FUNCTION + CHANGE TO X
Bp1=exp(3.2)+((5.7*t)/(t+13.0))*exp(-0.026*t)
f1=function(x) {exp(3.2)+((5.7*x)/(x+13.0))*exp(-0.026*x)}
Bp2=exp(3.2)+((6.3*t)/(t+5.0))*exp(-0.06*t)
f2=function(x) {exp(3.2)+((6.3*x)/(x+5.0))*exp(-0.06*x)}
uniroot(function(x) f1(x)-f2(x), c(12, 24))$root
#WAAROM ANDERS ALS c(12, 14)?
x <- exp(10)
#Create curve
curve(f1(x), 13.50909, 13.5091, col="red", ylab="biomarker concentration", xlab="time", ylim=c(26.576903, 26.57695))
curve(f2(x), 13.50909, 13.5091, col="blue", add=TRUE)
#Question 7: Intersect of both lines...
#Done manually need to use uniroot
#Question 8: See new line!
#Question 3: Again!
Bp1=exp(3.2)+((5.7*t)/(t+13.0))*exp(-0.026*t)
f1=function(x) {exp(3.2)+((5.7*x)/(x+13.0))*exp(-0.026*x)}
Bp2=exp(3.2)+((6.3*t)/(t+5.0))*exp(-0.06*t)
f2=function(x) {exp(3.2)+((6.3*x)/(x+5.0))*exp(-0.06*x)}
f1(0)
plot((x), f1, f2)
Deriv(f1, "x")
function (x)
{
.e1 <- 13 + 0
(5.7 - 0 * (0.1482 + 5.7/.e1)) * exp(-(0.026 * 0))/.e1
}
f=expression(exp(3.2)+((5.7*x)/(x+13.0))*exp(-0.026*x))
D(f, 'x')
ThisIs=(5.7/(x + 13) - (5.7 * x)/(x + 13)^2) * exp(-0.026 * x) - ((5.7 *
x)/(x + 13)) * (exp(-0.026 * x) * 0.026)
ThisIs
plot((x),ThisIs)
```
//QUESTION 8+
```{r}
#install.packages("deSolve")
library("deSolve")
#Question 8: new equasion
k1=0.3
k2=0.1
k3=0.9
b=2.0
S=10.0
R=3.0
MR=k1*S*(R-b)+(k2*(R-b)^2)-(k3*(R-b)^3)
#USING ODE
RThisProtein = function(timepoint, state, parameters) {
with(as.list(c(state, parameters)), {
dR = k1*S*(R-b)+(k2*(R-b)^2)-(k3*(R-b)^3)
list(c(dR))
})
}
pars = c(k1=0.3, k2=0.1, k3=0.9, b=2.0, S=30, 30, 30, 30, 30, seq(30, 5, by=1.66667))
inistate = c(R=3.0)
inistate.1 = c(R = 0)
inistate.2 = c(R = 2)
inistate.3 = c(R = 4)
inistate.4 = c(R = 6)
finish = seq(0, 2, by = 0.1)
library("deSolve")
out.1 = ode(y = inistate.1, times = finish, func = RThisProtein, parms = pars)
out.2 = ode(y = inistate.2, times = finish, func = RThisProtein, parms = pars)
out.3 = ode(y = inistate.3, times = finish, func = RThisProtein, parms = pars)
out.4 = ode(y = inistate.4, times = finish, func = RThisProtein, parms = pars)
plot(out.1, main = "time", ylab = "Synthesis of R", col="red", ylim=c(0, 6))
lines(out.2, col="Green")
lines(out.3, col="blue")
lines(out.4, col="yellow")
library("phaseR")
#Question 10: time 0-5
phaseport = phasePortrait(RThisProtein, parameters = pars, ylim = c(0, 20), points = 10,
xlab = "Time", ylab = "dR/dt", state.names = c("R"), main="Q13")
grid() # the function grid adds a grid to the plot
?title
```