Euler method, SEIR model
vishalstar
AnyLogic-Excerpt.pdf
Ilya Grigoryev · AnyLogic 7 in Three Days 101
System Dynamics modeling “System dynamics is a perspective and set of conceptual tools that enable us to understand the structure and dynamics of complex systems. System dynamics is also a rigorous modeling method that enables us to build formal computer simulations of complex systems and use them to design more effective policies and organizations.”
John Sterman, “Business Dynamics: Systems Thinking and Modeling for a Complex World”
The system dynamics method was created in 1950s by MIT Professor Jay Forrester. Drawing on his science and engineering background, Forrester sought to use the
laws of physics, in particular the laws of electrical circuits, to investigate economic
and social systems.
Today, system dynamics is typically used in long-term, strategic models, and it
assumes high levels of object aggregation: SD models represent people, products,
events, and other discrete items by their quantities.
System dynamics is a methodology to study dynamic systems. It suggests you:
�� Model the system as a causally closed structure that defines its own behavior.
� Discover the system's feedback loops (circular causality) balancing or reinforcing. Feedback loops are the heart of system dynamics.
� Identify stocks (accumulations) and flows that affect them.
Stocks are accumulations and characterize the system state. They are the memory
of the system and sources of disequilibrium. The model works only with aggregates
- the stock’s items are indistinguishable. Flows are the rates at which these system
states change.
If you’re having difficulty distinguishing between stocks and flows, consider how
we measure them. Stocks are usually expressed in quantities such as people,
inventory levels, money, or knowledge, while flows are typically measurements of
quantities in a given time period such as clients per month or dollars per year.
102 Ilya Grigoryev · AnyLogic 7 in Three Days
This chapter will to teach you how to use AnyLogic to develop system dynamics
models. If you want more information about the system dynamics approach, we
recommend Business Dynamics: Systems Thinking and Modeling for a Complex World by John Sterman.
Ilya Grigoryev · AnyLogic 7 in Three Days 103
SEIR model We're about to build a model that displays the spread of a contagious disease
among a large population. Our sample model will have a population of 10,000
people – a value we call TotalPopulation – of which one person is infectious.
�� During the infectious phase, a person comes into contact with an average of ContactRateInfectious = 1.25 people each day. If an infectious person comes into contact with a susceptible person, the susceptible person's
probability of infection is Infectivity = 0.6.
� After a susceptible person is infected, the infection latent phase lasts for AverageIncubationTime = 10 days. We'll use the word exposed to describe people who are in the latent phase.
� After the latent phase, infectious phase starts. This phase lasts for AverageIllnessDuration = 15 days.
� Persons who have recovered from the disease are immune to a second infection.
Phase 1. Creating a stock and flow diagram 1. Create a new model by selecting File > New > Model from the menu, and then
name it SEIR.
Let’s start with drawing stock and flow diagram. To model an epidemic's progress,
we need to reduce our population diversity. In this example, we'll consider four
important characteristics:
� Susceptible - people who are not infected by the virus
� Exposed - people who are infected but who can’t infect others
� Infectious - people who are infected and who can infect others
� Recovered – people who have recovered from the virus
SEIR is an acronym that represents the four stages: Susceptible-Exposed-
Infectious-Recovered. The terminology and the overall structure of the problem is
taken from the ("Compartmental models in epidemiology". n.d.) -- namely, from the
SEIR (Susceptible Exposed Infectious Recovered) model.
There are four stocks in our model - one for each stage.
104 Ilya Grigoryev · AnyLogic 7 in Three Days
2. Open the System Dynamics palette. Drag the Stock from the System Dynamics palette on to the diagram. Name it Susceptible.
3. Add three more stocks. Place them as shown in the figure and name them Exposed, Infectious, and Recovered.
Ilya Grigoryev · AnyLogic 7 in Three Days 105
Stocks and flows In System Dynamics, stocks (also known as levels, accumulations, or state variables) represent real-world stocks of material, knowledge, people, money, etc.
Flows define their rate of change - how stock values change and define the system’s dynamics. Here are some examples of stocks and flows:
Stock Inflows Outflows
Population Births Immigration
Deaths
Emigration
Fuel tank Refueling Fuel consumption
Flow may flow out of one stock and flow in another. Such a flow is outflow for the
first stock and inflow for the second one at the same time:
Flow may flow into a stock from nowhere. In this case the cloud (denoting "source")
is drawn at the flow’s starting point.
And symmetrically, flow may flow out from a stock to "nowhere". In this case the
cloud (denoting "sink") is drawn at the flow's end point.
The flow’s arrow shows its direction.
106 Ilya Grigoryev · AnyLogic 7 in Three Days
In our model, susceptible people are exposed to the virus, become infectious, and
then recover. It's a progression that requires our model to use three flows to drive
people from one stock to the next.
4. Add the first flow that flows from the stock Susceptible to Exposed. Double- click the stock where the flow flows out (Susceptible), and then click the stock where it flows in (Exposed).
5. Name the flow ExposedRate.
6. You can look at the formulas of Susceptible and Exposed stocks. Since our ExposedRate flow reduces the value of Susceptible stock and increases Exposed, the formulas should be the same as in the figures below. AnyLogic automatically created these formulas when you added the flow.
5
Ilya Grigoryev · AnyLogic 7 in Three Days 107
Formulas of stocks AnyLogic automatically generates a stock’s formula according to the user’s stock-
and-flow diagram.
Stock value is calculated according to flows flowing in and out from the stock. The
value of inflows – the flows that increase stock value – are added and the value of
outflows – flows that reduce the stock – are subtracted from the stock's current
value:
inflow1 + inflow2 … - outflow1 - outflow2 …
In the classic system, dynamics notation only flows can appear in the formula. The
formula is non-editable and no other elements other than flows flowing in and out
the stock can appear in the formula.
7. Add a flow from Exposed to Infectious, and then name it InfectiousRate.
8. Add a flow from Infectious to Recovered, and then name it RecoveredRate.
9. Rearrange the flow names as shown in the figure below. To do this, select a flow and then drag its name.
7
8
9
108 Ilya Grigoryev · AnyLogic 7 in Three Days
10. Now, let's define the parameters and dependencies. Add five Parameters, name them, and then define their default values according to the information below:
• TotalPopulation = 10 000
• Infectivity = 0.6
• ContactRateInfectious = 1.25
• AverageIncubationTime = 10
• AverageIllnessDuration = 15
11. Define the number of infected people by specifying 1 as the Initial Value of the stock Infectious.
12. Define the Initial Value for the stock Susceptible: TotalPopulation-1.
You may press Ctrl+space (Mac OS: Alt+space) and then select the parameter’s
name from the Code Completion assistant.)
10
Ilya Grigoryev · AnyLogic 7 in Three Days 109
You'll see the red sign to the expression’s left. The reason for the problem is you’ve
defined a dependency between two elements in the stock and flow diagram (the
stock Susceptible's initial value depends on the parameter TotalPopulation), but this dependency is not defined graphically as it should be.
12
110 Ilya Grigoryev · AnyLogic 7 in Three Days
Dependency links Stock and flow diagrams have two types of dependencies:
�� An element (stock, flow, auxiliary, or parameter) is mentioned in a flow or auxiliary's formula. This link type is drawn with a solid line:
� An element is mentioned in the stock's initial value. This link type is drawn with a dotted line:
You should use links to graphically define dependencies among a stock and flow diagram's elements.
If an element A is mentioned in the equation or element B’s initial value, you should first connect these elements with a link from A to B and then type the expression in B’s properties.
13. Draw a dependency link from TotalPopulation to Susceptible:
In the System Dynamics palette, double-click the Link element, click TotalPopulation, and then click the stock Susceptible. You should see the link with small circles drawn on its end points:
Ilya Grigoryev · AnyLogic 7 in Three Days 111
14. Let’s define the formula for the flow ExposedRate.
Click the flow and define the following formula using the Code Completion
assistant:
Infectious*ContactRateInfectious*Infectivity*Susсeptible/TotalPopulation
We need to draw dependency links from the mentioned variables and parameters
to this flow. You may find it tedious to manually draw the links, so we'll show you
how to add links using AnyLogic's link auto-creation mechanism.
15. Right-click ExposedRate flow in the graphical diagram, and choose Fix Variable Links > Create Missing Links from the context menu. Afterward, you should see the links in the stock and flow diagram:
13
14
112 Ilya Grigoryev · AnyLogic 7 in Three Days
16. Define the following formula for InfectiousRate: Exposed/AverageIncubationTime
17. Define the following formula for RecoveredRate: Infectious/AverageIllnessDuration
18. Draw the missing dependency links, and your stock and flow diagram should resemble the following image:
19. Adjust the appearance of dependency links. Modify the links’ bend angles to make the diagram match the figure below. To adjust the link's bend angle,
select it, and then drag the handle in the middle of the link.
Ilya Grigoryev · AnyLogic 7 in Three Days 113
20. Run the model and inspect the dynamics using the variables' inspect windows. To open a variable's inspect window, click the variable to select it. To resize the
window, drag its lower right corner.
21
20
114 Ilya Grigoryev · AnyLogic 7 in Three Days
21. To switch the inspect window to the plot mode, click the leftmost icon in its toolbar.
22. Increase the model execution speed to make the simulation go faster.
Ilya Grigoryev · AnyLogic 7 in Three Days 115
Phase 2. Adding a plot to visualize dynamics
Feedback loops: balancing and reinforcing System dynamics studies causal dependencies in systems. There are two types of
feedback loops: reinforcing and balancing.
Here are some hints how to know the loop type (taken from Wikipedia).
To determine if a causal loop is reinforcing or balancing, start with an assumption
such as "VariableN increases," and then follow the loop.
The loop is:
�� reinforcing if, after going around the loop, you get the same result as the initial assumption.
� balancing if the result contradicts the initial assumption.
You can also use the alternate definition:
� reinforcing loops have an even number of negative links (zero also is even)
� balancing loops have an uneven number of negative links.
We'll add a loop identifier for one loop to show you.
1. Drag the Loop element from the System Dynamics palette on to the diagram, and then place it as shown in the figure.
2. Go to the loop's Properties, change its Type to R (stands for Reinforcing), leave the default Clockwise Direction, and specify the text AnyLogic will display near the loop icon: Contagion.
1
116 Ilya Grigoryev · AnyLogic 7 in Three Days
Loop identifiers Loop is a graphical identifier with a label that briefly describes the loop's meaning and an arrow that shows the loop's direction.
Rather than defining the causal loop, it provides information about how your stock
and flow diagram's variables affect one another. By adding loops, you can help
other users understand the stock and flow diagram's influences and causal
dependencies.
Contagion loop is reinforcing. An increase in Infectious leads to an increase of ExposedRate, which in turn leads to a greater increase of Exposed. All links in this loop are positive.
Please find out what are other loops in this stock and flow diagram. What are their
directions and types?
Let’s add a time chart to plot Susceptible, Exposed, Infectious, and Recovered people.
3. Drag the Time Plot from the Analysis palette on to the diagram, and extend the time plot as shown in the figure below.
4. In the Properties view, go the section Data and click the button to add a new data item.
3
Ilya Grigoryev · AnyLogic 7 in Three Days 117
5. Modify the data item's properties:
�� Title: Susceptible people – title of the data item.
� Value: Susceptible (use Code Completion Master).
6. Add three data items to display the values of stocks Exposed, Infectious, and Recovered in the same way - and don't forget to define the corresponding Titles.
4
118 Ilya Grigoryev · AnyLogic 7 in Three Days
Ilya Grigoryev · AnyLogic 7 in Three Days 119
7. We've finished our last model. Now, run the model and use the chart you added to view its dynamics.
120 Ilya Grigoryev · AnyLogic 7 in Three Days
Phase 3. Parameter variation experiment In this phase, we’ll use AnyLogic’s parameter variation experiment to determine
how different contact rates affect the infection rate.
Parameter variation experiment The parameter variation experiment allows us to create a complex model
simulation that performs a series of single model runs that vary by one or many
parameters. After the experiment is complete, AnyLogic can display each run’s
results on a single diagram to help us better understand how the varying
parameters affected our model’s results.
If we run our experiment with fixed parameter values, we can also assess the effect
of random factors in stochastic models.
1. To add an experiment to the model, right-click the model item (SEIR) in the Projects tree, point to New, and then click Experiment.
2. In the New Experiment wizard’s Name field, type ContactRateVariation. AnyLogic will automatically select the Main agent type as the Top-level agent.
3. In the Experiment Type area, click Parameter Variation, and click Finish.
Ilya Grigoryev · AnyLogic 7 in Three Days 121
After you create the experiment, its diagram and properties will open.
4. In the experiment’s properties, open the Parameters section. The parameters of the experiment’s top-level agent (in our case, Main) display.
By default, all parameters are set as fixed. These values will not change during the parameter variation experiment.
3
122 Ilya Grigoryev · AnyLogic 7 in Three Days
5. To ensure our experiment varies the contact rate, locate the table’s ContactRateInfectious parameter and change its Type to Range.
6. Set the parameter’s minimum and maximum values by setting Min: 0.3 and Max: 2 with a Step of 0.1.
7. In the Properties area, click Create default UI.
The experiment diagram will display the simple user interface.
7
Ilya Grigoryev · AnyLogic 7 in Three Days 123
8. To ensure each run simulates exactly 300 days, we need to limit the model’s lifetime to 300 days. Click ContactRateVariation in the Projects tree to open its properties. In the Properties view, open the Model time section, select Stop at specified time from the Stop list, and type 300 in the Stop time box.
Now, we’ll add a time plot to display our experiment’s results. Our first step is to
gather data about the number of infectious persons.
9. Open the Main diagram, right-click the Infectious stock and then click Create Data Set.
124 Ilya Grigoryev · AnyLogic 7 in Three Days
10. After the InfectiousDS data set displays, navigate to its properties. Since we want to view the infectious disease’s dynamics, leave the Use time as horizontal axis value checkbox selected.
9
Ilya Grigoryev · AnyLogic 7 in Three Days 125
11. Select Update data automatically and leave Recurrence time: 1 to add one data sample to our dataset for each model life day.
12. To obtain data samples for the whole model run, set the dataset to Keep up to 300 latest samples.
We’re ready to add a chart to the ContactRateVariation experiment diagram that will display our results.
13. Open ContactRateVariation diagram, and drag the Time Plot from the Analysis palette.
14. Open the time plot’s properties. In the Scale section, ensure the time plot displays data for 300 model time units by setting the Time window to 300 model time units.
126 Ilya Grigoryev · AnyLogic 7 in Three Days
15. Enlarge the area available for the plot’s legend by dragging the diamond handle toward the top of the screen.
The plot’s curves will each display the results from one model run: the history of
disease spread for a contact rate collected by the InfectiousDS dataset.
16. Click ContactRateVariation in the Projects tree to open its properties, and then add data to the plot by navigating to the Java actions section and typing the following code in the After simulation run field: plot.addDataSet( root.InfectiousDS,
"CR=" + format( root.ContactRateInfectious ) );
While we want our plot to display multiple curves - one for each simulation run – we can’t add the data from each simulation run directly into the plot’s Data properties. Each model run destroys the top-level agent and its data, and since our
model can’t store the accumulated data, we need to manually add each curve to the
chart.
15
16
Ilya Grigoryev · AnyLogic 7 in Three Days 127
At the end of each simulation, AnyLogic stores the data in the Main agent’s InfectiousDS dataset. The experiment’s top-level agent is accessible as root, which means we access the data set as root.InfectiousDS.
We could use addDataset(root.InfectiousDS) to add a dataset to the chart with the default appearance and "Data set" title. However, we want AnyLogic to add a series of legends that will help us identify the plot’s curves, and that’s why we’ll use
another notation of addDataSet() function that has two arguments, addDataSet(DataSet ds, String title) .
We construct a legend for the dataset from the CR= label (where CR stands for contact rate) and the ContactRateInfectious value. Since the model’s top-level (root) agent defines this parameter, we’ll access its value as
root.ContactRateInfectious and use the function format(double value) to control how AnyLogic converts the values to text. For example, we want our chart’s legend
to display calculations such as 0.30001 as rounded values.
17. Open the ContactRateVariation experiment properties’ Advanced section and then clear the Allow parallel evaluations checkbox.
18. We’re ready to run the experiment and use our chart to observe the data we’ve gathered from multiple simulation runs. In the toolbar, on the Run list, select SEIR / ContactRateVariation.
19. In the presentation window, click Run.
128 Ilya Grigoryev · AnyLogic 7 in Three Days
AnyLogic performs a series of runs that all use a different value for the
ContactRateInfectious parameter and then add the simulation output to the chart.
If you review the Parameter variation experiment’s results, you’ll see how
increased contact rates allow the infection to spread more quickly. You’ll also see
that the plot displays 18 iterations – in other words, 18 infection scenarios for
contact rates that range from 0.3 to 2 – that represent the 18 steps in the parameter
variation range that we defined earlier. You can highlight a curve by clicking its
Contact Rate (CR) value in the legend.
19
