R help - Apr 2018 - Event-triggered change in value with a time-delay

Hello,

I am solving a set of ODEs using deSolve and have run into a problem I
would appreciate some advice on. One of the parameters (m) in the ODEs
changes between two states when one of the variables (D) crosses a
threshold (D_T) for the first time in either direction. Additionally, when
the variable crosses the threshold (either by increasing or decreasing),
there is a time delay (delay) before the parameter shifts values.
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Toy model:
# define parameters
P <- c(m = 10,
       R_1 = 0.5,
       R_2 = 0.1,
       D_0 = 50,
       D_T = 5,
       delay = 3
)
# set initial values
y_0 <- c(D = 0, Y = 1)

# set timesteps
times <- seq(0, 100, length = 101)

# define function
bit <- function(t,y,parms){
  with(as.list(c(parms,y)),{


    # How I would implement the change if there was no time delay
    # m <- ifelse(D>=D_T, m*.68, m)

    # Option 1: nested if statements with time delay
    # if (D >= D_T) {
    #   t_start1 <- t + delay
    #   if (t >= t_start1){
    #     m <- 6.8
    #   }
    # }
    # if (D < D_T) {
    #   t_start2<- t + delay
    #   if (t >= t_start2){
    #     m <- 10
    #   }
    # }


    dD <- D_0 * (R_1 / (R_1 - R_2)) * (-R_2 * exp(-R_2 * t) + R_1 *
exp(-R_1 * t))
    dY <- Y + m;
    res <- c(dD, dY)
    list(res, m = m)
  })
}


# Solve model
library(deSolve)
out <- ode(y = y_0, times = times, func = bit, parms = P)
out.df <- as.data.frame(out)
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
I have come up with a number of ideas, but I either haven't been able to
get them to work or don't know if R has a command for it:

1. Nested "if" statements with something to track time (as shown
above).
One issue here is* t-start *is redefined with each time step in the ode
solver so *t* never reaches *t_start +delay*.

2. Event triggered changes. From the examples I have found (link
<https://cran.r-project.org/web/packages/deSolve/vignettes/deSolve.pdf> p.
24-25), the time at which the event happens is either pre-set or occurs at
a root of the ODE.  Also- would this strategy require multiple events? One
event in which the variable crossing the threshold triggers a timer,
another for the timer triggering the parameter change?

3. Lagvalue. This might take care of the time lag component, but I am still
struggling to understand (1) how lagvalue works and (2) how to integrate
both the threshold and timer components.

   if ( D < D_T )
     m<- ifelse(lagvalue(t-delay)<=0, 6.8, 10)
   end

4. Sigmoidal function. Someone suggested this to me earlier, but I am not
sure how this would work. Is it possible to have a sigmoidal function that
switches repeatedly between two states (as opposed to the characteristic S
curve that transitions once from one state to another)?

5. Mathmematica has a WhenEvent function- is there something similar in R?


Thank you for your time and any tips you might have!

	[[alternative HTML version deleted]]

> On Apr 18, 2018, at 1:04 AM, Hannah Meredith <hrmeredith12 at
gmail.com> wrote:
> 
> Hello,
> 
> I am solving a set of ODEs using deSolve and have run into a problem I
> would appreciate some advice on. One of the parameters (m) in the ODEs
> changes between two states when one of the variables (D) crosses a
> threshold (D_T) for the first time in either direction. Additionally, when
> the variable crosses the threshold (either by increasing or decreasing),
> there is a time delay (delay) before the parameter shifts values.
I think you migh not have fully reviewed the available functions in that
package. Doing an sos::findFn-search on ODE's with time delay brings up help
pages for the `dede` function in deSolve.

For the conditional value of the within regime-fixed parameter `m`, couldn't
you just use this:

    ... + m*(6.8+ 3.2*(D < D_T))

Caveat: very lightly tested in your toy example, which did generate a result
where D cross that threshold without blowing up.

-- 

David.

> ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> Toy model:
> # define parameters
> P <- c(m = 10,
>       R_1 = 0.5,
>       R_2 = 0.1,
>       D_0 = 50,
>       D_T = 5,
>       delay = 3
> )
> # set initial values
> y_0 <- c(D = 0, Y = 1)
> 
> # set timesteps
> times <- seq(0, 100, length = 101)
> 
> # define function
> bit <- function(t,y,parms){
>  with(as.list(c(parms,y)),{
> 
> 
>    # How I would implement the change if there was no time delay
>    # m <- ifelse(D>=D_T, m*.68, m)
> 
>    # Option 1: nested if statements with time delay
>    # if (D >= D_T) {
>    #   t_start1 <- t + delay
>    #   if (t >= t_start1){
>    #     m <- 6.8
>    #   }
>    # }
>    # if (D < D_T) {
>    #   t_start2<- t + delay
>    #   if (t >= t_start2){
>    #     m <- 10
>    #   }
>    # }
> 
> 
>    dD <- D_0 * (R_1 / (R_1 - R_2)) * (-R_2 * exp(-R_2 * t) + R_1 *
> exp(-R_1 * t))
>    dY <- Y + m;
>    res <- c(dD, dY)
>    list(res, m = m)
>  })
> }
> 
> 
> # Solve model
> library(deSolve)
> out <- ode(y = y_0, times = times, func = bit, parms = P)
> out.df <- as.data.frame(out)
> ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> I have come up with a number of ideas, but I either haven't been able
to
> get them to work or don't know if R has a command for it:
> 
> 1. Nested "if" statements with something to track time (as shown
above).
> One issue here is* t-start *is redefined with each time step in the ode
> solver so *t* never reaches *t_start +delay*.
> 
> 2. Event triggered changes. From the examples I have found (link
>
<https://cran.r-project.org/web/packages/deSolve/vignettes/deSolve.pdf> p.
> 24-25), the time at which the event happens is either pre-set or occurs at
> a root of the ODE.  Also- would this strategy require multiple events? One
> event in which the variable crossing the threshold triggers a timer,
> another for the timer triggering the parameter change?
> 
> 3. Lagvalue. This might take care of the time lag component, but I am still
> struggling to understand (1) how lagvalue works and (2) how to integrate
> both the threshold and timer components.
> 
>   if ( D < D_T )
>     m<- ifelse(lagvalue(t-delay)<=0, 6.8, 10)
>   end
> 
> 4. Sigmoidal function. Someone suggested this to me earlier, but I am not
> sure how this would work. Is it possible to have a sigmoidal function that
> switches repeatedly between two states (as opposed to the characteristic S
> curve that transitions once from one state to another)?
> 
> 5. Mathmematica has a WhenEvent function- is there something similar in R?
> 
> 
> Thank you for your time and any tips you might have!
> 
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.
David Winsemius
Alameda, CA, USA

'Any technology distinguishable from magic is insufficiently advanced.' 
-Gehm's Corollary to Clarke's Third Law

Hi David,

Thank you for your suggestion- I will look into it some more.

Thanks again,
Hannah

On Wed, Apr 18, 2018 at 5:09 PM, David Winsemius <dwinsemius at
comcast.net>
wrote:
>
> > On Apr 18, 2018, at 1:04 AM, Hannah Meredith <hrmeredith12 at
gmail.com>
> wrote:
> >
> > Hello,
> >
> > I am solving a set of ODEs using deSolve and have run into a problem I
> > would appreciate some advice on. One of the parameters (m) in the ODEs
> > changes between two states when one of the variables (D) crosses a
> > threshold (D_T) for the first time in either direction. Additionally,
> when
> > the variable crosses the threshold (either by increasing or
decreasing),
> > there is a time delay (delay) before the parameter shifts values.
>
> I think you migh not have fully reviewed the available functions in that
> package. Doing an sos::findFn-search on ODE's with time delay brings up
> help pages for the `dede` function in deSolve.
>
> For the conditional value of the within regime-fixed parameter `m`,
> couldn't you just use this:
>
>     ... + m*(6.8+ 3.2*(D < D_T))
>
> Caveat: very lightly tested in your toy example, which did generate a
> result where D cross that threshold without blowing up.
>
> --
>
> David.
>
>
> > ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> > Toy model:
> > # define parameters
> > P <- c(m = 10,
> >       R_1 = 0.5,
> >       R_2 = 0.1,
> >       D_0 = 50,
> >       D_T = 5,
> >       delay = 3
> > )
> > # set initial values
> > y_0 <- c(D = 0, Y = 1)
> >
> > # set timesteps
> > times <- seq(0, 100, length = 101)
> >
> > # define function
> > bit <- function(t,y,parms){
> >  with(as.list(c(parms,y)),{
> >
> >
> >    # How I would implement the change if there was no time delay
> >    # m <- ifelse(D>=D_T, m*.68, m)
> >
> >    # Option 1: nested if statements with time delay
> >    # if (D >= D_T) {
> >    #   t_start1 <- t + delay
> >    #   if (t >= t_start1){
> >    #     m <- 6.8
> >    #   }
> >    # }
> >    # if (D < D_T) {
> >    #   t_start2<- t + delay
> >    #   if (t >= t_start2){
> >    #     m <- 10
> >    #   }
> >    # }
> >
> >
> >    dD <- D_0 * (R_1 / (R_1 - R_2)) * (-R_2 * exp(-R_2 * t) + R_1 *
> > exp(-R_1 * t))
> >    dY <- Y + m;
> >    res <- c(dD, dY)
> >    list(res, m = m)
> >  })
> > }
> >
> >
> > # Solve model
> > library(deSolve)
> > out <- ode(y = y_0, times = times, func = bit, parms = P)
> > out.df <- as.data.frame(out)
> > ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> > I have come up with a number of ideas, but I either haven't been
able to
> > get them to work or don't know if R has a command for it:
> >
> > 1. Nested "if" statements with something to track time (as
shown above).
> > One issue here is* t-start *is redefined with each time step in the
ode
> > solver so *t* never reaches *t_start +delay*.
> >
> > 2. Event triggered changes. From the examples I have found (link
> >
<https://cran.r-project.org/web/packages/deSolve/vignettes/deSolve.pdf>
> p.
> > 24-25), the time at which the event happens is either pre-set or
occurs
> at
> > a root of the ODE.  Also- would this strategy require multiple events?
> One
> > event in which the variable crossing the threshold triggers a timer,
> > another for the timer triggering the parameter change?
> >
> > 3. Lagvalue. This might take care of the time lag component, but I am
> still
> > struggling to understand (1) how lagvalue works and (2) how to
integrate
> > both the threshold and timer components.
> >
> >   if ( D < D_T )
> >     m<- ifelse(lagvalue(t-delay)<=0, 6.8, 10)
> >   end
> >
> > 4. Sigmoidal function. Someone suggested this to me earlier, but I am
not
> > sure how this would work. Is it possible to have a sigmoidal function
> that
> > switches repeatedly between two states (as opposed to the
characteristic
> S
> > curve that transitions once from one state to another)?
> >
> > 5. Mathmematica has a WhenEvent function- is there something similar
in
> R?
> >
> >
> > Thank you for your time and any tips you might have!
> >
> >       [[alternative HTML version deleted]]
> >
> > ______________________________________________
> > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> > 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.
>
> David Winsemius
> Alameda, CA, USA
>
> 'Any technology distinguishable from magic is insufficiently
advanced.'
>  -Gehm's Corollary to Clarke's Third Law
>
>
>
>
>
>

-- 
Hannah Meredith | PhD
Whitaker Scholar
London School of Hygiene and Tropical Medicine

	[[alternative HTML version deleted]]