[R] Solution to differential equation for nls function

[Ravi Varadhan]

Yes, there is an analytical solution.

Here is how you do it:

bernoulli_anal <- function(t, R0, Rmax, a, b) {
# Note:  a, b, R0, Rmax all have to be positive (strictly > 0)
R1 <- R0 / Rmax
k <- R1^b / (1 - R1^b)
R <- Rmax / (1 + exp(-a*b*t)/k)^(1/b)
return(R)
}

R0 <- 1
Rmax <- 2
a <- 0.5
b <- 1.0

t <- seq(0, 10, length=1000)
y1 <- bernoulli_anal(t, R0, Rmax, a, b)
y2 <- bernoulli_anal(t, R0, Rmax, a, b=0.5)
y3 <- bernoulli_anal(t, R0, Rmax, a, b=1.5)
y4 <- bernoulli_anal(t, R0, Rmax, a, b=2)
plot(t, y1, type="l", ylab="R(t)")
lines(t, y2, col=2)
lines(t, y3, col=3)
lines(t, y4, col=4)

Hope this helps,
Ravi.


-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
Behalf Of mahesh samtani
Sent: Wednesday, June 30, 2010 10:44 AM
To: r-help at r-project.org
Subject: [R] Solution to differential equation for nls function

Hello,
I am trying to find the analytical solution to this differential equation

dR/dt = k1*R (1-(R/Rmax)^k2); R(0) = Ro

k1, k2, and Rmax are parameters that need to fitted, while Ro is the
baseline value (which can be fitted or fixed). The response (R) increases
initially at an exponential rate governed by the rate constant k1.
Response has a S-shaped curve as a function of time and it approaches the
value of Rmax at time approaches infinity.

If there is an analytial solution to this differential equation then it
makes my life easier when trying to perform some non-linear regression.
Kindly provide the integration process so I can learn how to do it myself
for future reference. I believe that this is a bernoulli differential
equation, which could help with integration. I guess the harder way would be
to use integration by parts (I tried hard to find the solution but keep
getting stuck).

Please help,
Mahesh

	[[alternative HTML version deleted]]

______________________________________________
R-help at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.