7.6.1 Numerical simulations of the SIR model
In the following section the data fit without direct influence of rain data is called Simu-lation 1 and SimuSimu-lation 2 includes the rain data. The numerical and graphical results of these simulations are represented in Table 7.5 and 7.7. The magnitude of the parameters β0, γ and ϕ is the same whereby Simulation 2 shows a higher average transmission rate forβ0. To compare the degree of periodic variation we examine the size of the term
β1
Z t−η1
52
t−η522
pc(τ)dτ < β1·η2−η1
52 ·c= 45.50·1.8
52 ·16 = 25.20. (7.14) In Simulation 2 this term fluctuates between 0 and 25.20 whereas in Simulation 1 we obtain β1 = 7.53. The cut–off c and the interval limits [t−η2/52;t−η1/52] fit to the assumptions in Section 7.3 since 15mm ≤ c ≤ 20mm and η2 = 7.80 are sensible values in terms of the cluster and cross–correlation between the moving average of rainfall and dengue data. The hospitalization rate γ is in the expected range of 25%. Due to the optimization of the phase shiftϕ, the timing of the model fits well to the dengue peaks in the data set. When comparing the respective graphs, it is noticeable that a more realistic dynamic can be recognized by adding the rain data. This is particularly noticeable because certain fluctuations during the individual periods are reflected, such as in 2009. It is also striking that between 2015 and 2016 Simulation 2 does not reproduce the very low values in contrast to Simulation 1. The rain data indicate that the longest period with very low precipitation was before 2015 which is reflected in the dengue data but not in the model.
On the one hand sideSIRmodels are only suitable for short periods of time and the start time of the data fit should be set later. But this would be also true for Simulation 1. A possible explanation is that 2016 was an El Nino year which could lead to outliers since then the functional form for the seasonality could be wrong. Further deviations can be found in various patches in 2010, 2013 and 2014 because in these cases the timing of the peaks differs from other years. Overall, it should be noted that the different intensities of the dengue peaks and their relation to each other are well represented by the SIR model. The results also allow conclusions about the behavior of the mosquito population.
Provided the size of the vector population can be modeled byMi(t) =Mi0eωχsin(ω(t+ϕ/52)), the term
eχωsin(ω(t+52ϕ))≈1 + χ ωsin
ω t+ ϕ
52
describes the seasonal variation of the mosquito population. As shown in Section 7.4 in Simulation 1 the fraction ββ1
0 approximately corresponds to the size of the expression χω. Since
χωsin (ω(t+ϕ/52))
≤0.16 applies for ωχ = ββ1
0 we can conclude that the size of the mosquito population varies by a maximum of approximately 20%. In Simulation 2 the
termββ1
0
Rt−η1/52
t−η2/52 pc(τ)dτcorresponds to a time–dependent expression χ(t)ω withχ:D→R+0. Using the result from (7.14) we get
χ(t)
ω sin (ω(t+ϕ/52))
≤ 0.44 and consequently a variation of at most about 50%. The graph of the transmission rate β(t) in the second simulation reflects the dynamic processes within the model. For example, the effects of the cut–off are visible in the years 2014–2015, whenβ(t) assumes the valueβ0 over longer periods of time. In contrast, the very dry time at the end of 2015 is clearly shown inβ(t) in the same way.
Table 7.6: Results of the parameter fit andβ(t) in Simulation 1 and 2.
Parameters β0 β1 c η2 ϕ γ α µ ω η1
Simulation 1 47.31 7.53 / / 10.48 0.24 26 691 2π /
2009 2010 2011 2012 2013 2014 2015 2016 2017 40
45 50 55
(t)
Transmission Rate
β(t) =β0+β1sin (ω(t+ϕ/52)) β(t)∈[39.74; 54.84]
Parameters β0 β1 c η2 ϕ γ α µ ω η1
Simulation 2 57.43 45.50 16 7.80 9.00 0.25 26 691 2π 6
2009 2010 2011 2012 2013 2014 2015 2016 2017 40
50 60 70 80
(t)
Transmission Rate
β(t) =β0+β1
Rt−η1/52
t−η2/52 pc(τ)dτsin (ω(t+ϕ/52)) β(t)∈[43.83; 75.28]
7.6.2 Prediction quality of the SIR model
The prediction quality of the model is the basis for further application with regard to various control methods. For this reason, we test the presented model using the available data sets. In the following, the parameters inu= (β0, β1, c, η2, ϕ, γ, S0, R0)T are fitted up to a time ˆt∈[t0;t1]. Then the model continues with these findings and makes a prediction to the end time 2017. In Simulation 3 we use the available rain data and in Simulation 4 the average rainfall data of previous years within the prognosis interval (see Tables 7.10 and 7.11). To give more weight to current than to past data we introduce a weight function H:D→R+ with
H(t) =w·exp − t−ˆt2
2σ2
! +z .
The parameters w = 50, σ = 524 and z = 1 are selected so that the period of the last four weeks before ˆtis weighted considerably more strongly. Consequently, we obtain the
Table 7.7: Results of the parameter fit in Simulation 1 and 2.
Simulation 1 Simulation 2
Si0 Ii0 Ri0 Ni Si0 Ii0 Ri0 Ni
20090 2010 2011 2012 2013 2014 2015 2016 2017 100
200 300 400
Dengue Cases
South Jakarta
Id(t) I(t) (fitted)
20090 2010 2011 2012 2013 2014 2015 2016 2017 100
200 300 400
Dengue Cases
South Jakarta
Id(t) I(t) (fitted)
1166192 612 1018908 2185711 935017 658 1250036 2185711
20090 2010 2011 2012 2013 2014 2015 2016 2017 200
400 600
Dengue Cases
East Jakarta
Id(t) I(t) (fitted)
20090 2010 2011 2012 2013 2014 2015 2016 2017 200
400 600
Dengue Cases
East Jakarta
Id(t) I(t) (fitted)
1664665 632 1178518 2843816 1343231 627 1499958 2843816
20090 2010 2011 2012 2013 2014 2015 2016 2017 50
100 150 200
Dengue Cases
Central Jakarta
Id(t) I(t) (fitted)
20090 2010 2011 2012 2013 2014 2015 2016 2017 50
100 150 200
Dengue Cases
Central Jakarta
Id(t) I(t) (fitted)
404974 313 508895 914182 313890 340 599952 914182
20090 2010 2011 2012 2013 2014 2015 2016 2017 100
200 300 400
Dengue Cases
West Jakarta
Id(t) I(t) (fitted)
20090 2010 2011 2012 2013 2014 2015 2016 2017 100
200 300 400
Dengue Cases
West Jakarta
Id(t) I(t) (fitted)
1279767 112 1183681 2463560 1013347 158 1450055 2463560
20090 2010 2011 2012 2013 2014 2015 2016 2017 100
200 300
Dengue Cases
North Jakarta
Id(t) I(t) (fitted)
20090 2010 2011 2012 2013 2014 2015 2016 2017 100
200 300
Dengue Cases
North Jakarta
Id(t) I(t) (fitted)
922332 479 824504 1747315 746904 466 999945 1747315
minimization problem
minu J(u) = min
u n
X
i=1
Z tˆ t0
H(t) γIi(t)−Iid(t)2
maxt∈DIid(t) dt+ kuk
N 2
.
In practice, the optimization up to ˆt should be constantly updated. In our simulations we have optimized in most cases until eight weeks after the turn of the year. The reason for this is that in Jakarta the number of dengue cases increases significantly every year during this period and the model is to be tested for its prediction quality for the following season. Additionally we calculate based on theL1 norm
E1 = Z ˆt+1
ˆt
γIi(t)−Iid(t)
maxt∈[t;ˆˆt+1]Iid(t)dt and E2= 1 t1−ˆt
Z t1
ˆt
γIi(t)−Iid(t) maxt∈[t;tˆ 1]Iid(t)dt .
These values are used to determine the deviation of the model in relation to the corre-sponding maximum value within the data. Although the forecast for the coming season is in the foreground, the model also reveals tendencies in the following years. For example, the relatively high increase of registered cases in 2016 is indicated not only in the direct forecast of the corresponding year, but also in the long–term forecasts. In most cases, the model provides information about the expected size of dengue cases in the coming season.
It is noticeable that in most cases the relation of the forthcoming peak to the previous peaks is quite realistic. With regard to the short–term forecast, it makes no significant difference whether the model continues to calculate with the actual or average rain data.
Large deviations between Simulation 3 and 4 arise only with the long–term forecast until 2017.
Remark 7.6.3. Note that, although the model is very homogeneous among the regions, due to the same hospitality rate and the same ratio of mosquito to humans and other parameters, the fitting results for the model are in good agreement with the data. This is also an effect of the mobility matrix, which introduces a spatial–inhomogeneity into the system. In a forthcoming study this particular role of the mobility matrix will be considered.
Table 7.8: Quantitative comparision of the prediction error for South Jakarta.
Simulation 3: The parameter fit is executed in the interval
2009; ˆt . The forecast based on this is carried out with the actual rain data for the period t; 2017ˆ
.
Simulation 4: The forecast in ˆt; 2017
is done with the average rain data from the period
2009; ˆt .
ˆt= 2016 + 8weeks tˆ= 2015 + 8weeks ˆt= 2014 + 8weeks E1 E2 E1 E2 E1 E2
Simulation 3 0.12 0.12 0.21 0.12 0.3 0.41
Simulation 4 0.26 0.26 0.17 0.15 0.39 0.67
ˆt= 2013 + 8weeks tˆ= 2012 + 8weeks ˆt= 2015 + 24weeks E1 E2 E1 E2 E1 E2
Simulation 3 0.22 0.18 0.33 0.36 0.17 0.13
Simulation 4 0.19 0.68 0.30 0.56 0.14 0.17
From Table 7.8, and the numerical findings in Appendix 7.C we see that the Simulation 3 has always smaller prediction errors. Of course this is to expect, since more detailed information about the rainfall in the forecasted period is incorporated. Although the error is slightly higher, it is remarkable, that peaks of the prediction of Simulation 4 are in a very good agreement with the peaks in the dengue data, although just average rainfall data from the past timeset was taken. It is clear that the errors are the higher the more time has to be predicted or vice versa the less information is used for a prediction.
Remark 7.6.4. Of course for an exact prediction in practice one has to update the data sets during the time. The model used here is due to several reductions just in a good
agreement for a short time horizon. Moreover as one can see e.g. in the prediction of the years 2015 and 2016 in Table 7.7 the peaks of the outbreak can not be seen clearly.
Another reason here could be an uncertainty in the climatic input.