Message Boards

WOLFRAM COMMUNITY

5636 Views

4 Replies

1 Total Likes

View groups...

Follow this post

Share this post:

GROUPS:

Mathematics Wolfram Language

[?] Implement simple birth death process?

Okkes Dulgerci

Posted 7 years ago

Hi All, I am trying to implement simple birth death process. Why my code does not work? Any suggestion. Thanks. birthDeath[?_, ?_, initialPopulation_, numOfReaction_] := NestList[( ?t1 = RandomVariate[ExponentialDistribution[? #]]; ?t2 = RandomVariate[ExponentialDistribution[? #]]; If[?t1 < ?t2, {# + 1}, {# - 1}]) &, initialPopulation, numOfReaction] birthDeath[3, 1, 10, 2]

Hi All, I am trying to implement simple birth death process. Why my code does not work? Any suggestion. Thanks.

birthDeath[?_, ?_, initialPopulation_, numOfReaction_] :=     
 NestList[(
    ?t1 =  RandomVariate[ExponentialDistribution[?  #]];
    ?t2 = RandomVariate[ExponentialDistribution[? #]]; 
    If[?t1 < ?t2, {# + 1}, {# - 1}]) &,  initialPopulation, numOfReaction]

birthDeath[3, 1, 10, 2]

POSTED BY: Okkes Dulgerci

4 Replies

Sort By:

Gianluca Gorni

Gianluca Gorni, University of Udine

Posted 7 years ago

This may be one way: birthDeath2[\[Lambda]_, \[Mu]_, initialPopulation_, numOfReaction_] := NestList[(\[CapitalDelta]t1 = RandomVariate[ExponentialDistribution[\[Lambda] #[[3]]]]; \[CapitalDelta]t2 = RandomVariate[ExponentialDistribution[\[Mu] #[[3]]]]; {\[CapitalDelta]t1, \[CapitalDelta]t2, If[\[CapitalDelta]t1 < \[CapitalDelta]t2, #[[3]] + 1, #[[3]] - 1]}) &, {0, 0, initialPopulation}, numOfReaction] birthDeath2[3, 1, 10, 4]

This may be one way:

birthDeath2[\[Lambda]_, \[Mu]_, initialPopulation_, numOfReaction_] :=
  NestList[(\[CapitalDelta]t1 = 
     RandomVariate[ExponentialDistribution[\[Lambda] #[[3]]]];
    \[CapitalDelta]t2 = 
     RandomVariate[ExponentialDistribution[\[Mu] #[[3]]]];
    {\[CapitalDelta]t1, \[CapitalDelta]t2, 
     If[\[CapitalDelta]t1 < \[CapitalDelta]t2, #[[3]] + 1, #[[3]] - 
       1]}) &, {0, 0, initialPopulation}, numOfReaction]

birthDeath2[3, 1, 10, 4]

POSTED BY: Gianluca Gorni

Okkes Dulgerci

Posted 7 years ago

I fixed it. exact[\[Lambda]_, \[Mu]_, initialPopulation_] := Module[{}, DSolveValue[{p'[t] == (\[Lambda] - \[Mu]) p[t], p[0] == initialPopulation}, p[t], t]]; birthDeath[\[Lambda]_, \[Mu]_, initialPopulation_, numOfReaction_] := NestList[(\[CapitalDelta]t1 = RandomVariate[ExponentialDistribution[\[Lambda] #[[2]]]]; \[CapitalDelta]t2 = RandomVariate[ExponentialDistribution[\[Mu] #[[2]]]]; \[CapitalDelta]t = Min[\[CapitalDelta]t1, \[CapitalDelta]t2]; {#[[1]] + \[CapitalDelta]t, If[\[CapitalDelta]t1 < \[CapitalDelta]t2, #[[2]] + 1, #[[2]] - 1]}) &, {0, initialPopulation}, numOfReaction] With[{\[Lambda] = 3, \[Mu] = 1, initialPopulation = 10, numOfReaction = 1000, numOfsim = 10}, sim = Table[ birthDeath[\[Lambda], \[Mu], initialPopulation, numOfReaction], numOfsim]; Show[ListStepPlot[sim, PlotLegends -> {"Simulation"}, PlotStyle -> Directive[AbsoluteThickness[0.2]], Frame -> True, PlotTheme -> "Detailed", FrameLabel -> {"Time", "Population"}, ImageSize -> Large], Plot[Evaluate@exact[\[Lambda], \[Mu], initialPopulation], {t, 0, Max@sim[[All, -1, 1]]}, PlotStyle -> {Black, Thick}, PlotLegends -> {"ODE"}], ListLinePlot[Mean@sim, PlotLegends -> {"Avg of Simulation"}, PlotStyle -> {Red, Dashed}]]]

I fixed it.

exact[\[Lambda]_, \[Mu]_, initialPopulation_] := 
 Module[{}, 
  DSolveValue[{p'[t] == (\[Lambda] - \[Mu]) p[t], 
    p[0] == initialPopulation}, p[t], t]]; 
birthDeath[\[Lambda]_, \[Mu]_, initialPopulation_, numOfReaction_] := 
 NestList[(\[CapitalDelta]t1 = 
     RandomVariate[ExponentialDistribution[\[Lambda] #[[2]]]];
    \[CapitalDelta]t2 = 
     RandomVariate[ExponentialDistribution[\[Mu] #[[2]]]];
    \[CapitalDelta]t = Min[\[CapitalDelta]t1, \[CapitalDelta]t2];
    {#[[1]] + \[CapitalDelta]t, 
     If[\[CapitalDelta]t1 < \[CapitalDelta]t2, #[[2]] + 1, #[[2]] - 
       1]}) &, {0, initialPopulation}, numOfReaction]


With[{\[Lambda] = 3, \[Mu] = 1, initialPopulation = 10, 
  numOfReaction = 1000, numOfsim = 10},
 sim = Table[
   birthDeath[\[Lambda], \[Mu], initialPopulation, numOfReaction], 
   numOfsim];
 Show[ListStepPlot[sim, PlotLegends -> {"Simulation"}, 
   PlotStyle -> Directive[AbsoluteThickness[0.2]], Frame -> True, 
   PlotTheme -> "Detailed", FrameLabel -> {"Time", "Population"}, 
   ImageSize -> Large], 
  Plot[Evaluate@exact[\[Lambda], \[Mu], initialPopulation], {t, 0, 
    Max@sim[[All, -1, 1]]}, PlotStyle -> {Black, Thick}, 
   PlotLegends -> {"ODE"}], 
  ListLinePlot[Mean@sim, PlotLegends -> {"Avg of Simulation"}, 
   PlotStyle -> {Red, Dashed}]]]

POSTED BY: Okkes Dulgerci