Priority and age specific vaccination algorithm for the pandemic diseases: a comprehensive parametric prediction model

Background There have been several destructive pandemic diseases in the human history. Since these pandemic diseases spread through human-to-human infection, a number of non-pharmacological policies has been enforced until an effective vaccine has been developed. In addition, even though a vaccine has been developed, due to the challenges in the production and distribution of the vaccine, the authorities have to optimize the vaccination policies based on the priorities. Considering all these facts, a comprehensive but simple parametric model enriched with the pharmacological and non-pharmacological policies has been proposed in this study to analyse and predict the future pandemic casualties. Method This paper develops a priority and age specific vaccination policy and modifies the non-pharmacological policies including the curfews, lockdowns, and restrictions. These policies are incorporated with the susceptible, suspicious, infected, hospitalized, intensive care, intubated, recovered, and death sub-models. The resulting model is parameterizable by the available data where a recursive least squares algorithm with the inequality constraints optimizes the unknown parameters. The inequality constraints ensure that the structural requirements are satisfied and the parameter weights are distributed proportionally. Results The results exhibit a distinctive third peak in the casualties occurring in 40 days and confirm that the intensive care, intubated, and death casualties converge to zero faster than the susceptible, suspicious, and infected casualties with the priority and age specific vaccination policy. The model also estimates that removing the curfews on the weekends and holidays cause more casualties than lifting the restrictions on the people with the chronic diseases and age over 65. Conclusion Sophisticated parametric models equipped with the pharmacological and non-pharmacological policies can predict the future pandemic casualties for various cases.

constraints on fighting the pandemic diseases. In this case, the state authorities seek policies that optimize the respective priorities such as reducing the deaths, easing the curfews, lifting the restrictions, and opening of the schools. This paper proposes a comprehensive parametric model with the priority and age specific vaccination policy which can be used for the prediction and analysis of the future casualties under the constructed policy. In this paper, the healthcare staff constitutes the highest priority group, and the elderly people are located into the risk groups based on their ages. Accurate models can be a useful tool to understand the dynamics of the pandemic diseases and to identify the role of the internal (mutation) and external (pharmacological and non-pharmacological policies) impacts on the pandemic casualties [3]. Modelling of the pandemic diseases can be achieved by using the non-parametric (statistical and machine learning), and the parametric (mathematical) approaches. Statistical approaches usually reveal a mean and a standard deviation which can be used for characterizing the pandemic properties such as the incubation period and the infectious rate of the pandemic diseases. Overton et al. produced data for the incubation period of the COVID-19 by assuming that it has the Gama distribution and fitted the data with the maximum likelihood estimator [4]. This research stated that the majority of the infected people develop symptoms in 14 days. Hong and Li estimated a time-dependent reproduction number of disease with the Poisson model having a removal rate to account for the random uncertainties in the reported casualties [5]. It is concluded that China, Italy, Sweden, and the United States of America (USA) have high COVID-19 reproduction numbers since they were unable to control the spread of the virus. Oehmke et al. determined the speed, acceleration, jerk, and 7-day-lag in the COVID-19 transmission for the USA and determined the parameters with the Arellano-Bond statistical estimator [6]. It is expressed that there were significant differences in the spread of the virus among the states of the USA due to lack of a national non-pharmacological policy.
In terms of the machine learning based modelling approaches, Pinter et al. proposed an adaptive networkbased fuzzy inference systems (ANFIS) and a multilayered perceptron-imperialist competitive algorithm (MLP-ICA) to estimate the infected individuals and the mortality rate [7]. Tuli et al. considered the generalized Inverse Weibull distribution combined with the cloud computing to predict the growth of the epidemic and to design control strategies for the COVID-19 spread [8]. Aydin and Yurdakul evaluated the policies of the 142 countries against fighting the COVID-19 by using the k-means clustering, decision trees, and random forest algorithms [9]. The research revealed that the economic welfare, smoking rates, and the diabetes rates are not directly related to the effectiveness level of the countries. Rustem et al. utilized the linear regression (LR), the least absolute shrinkage and selection operator (LASSO), the support vector machines (SVM), and the exponential smoothing (ES) to estimate the threating factors of the COVID-19 [10]. The results confirmed that the ES outperforms the others while the SVM performs poorly. Bird et al. evaluated the country-level pandemic risks and classified the preparedness of the countries in terms of the transmission, the mortality, and the inability to test by applying the stack of gradient boosting, the decision trees, the stack of SVM, and the extra trees [11]. It is concluded that the geopolitics and the demographic attributes shape the risks caused by the COVID-19.
Parametric modelling approaches suit their purpose and are also parameterizable by the available data, since they have a certain model structure representing the mathematical relationships as simple as possible [12]. Goel and Sharma proposed a mobility-based susceptible, infected, recovered (SIR) model covering the population distribution and the lockdowns [13]. It is observed that the infected casualties are delayed and decreased in the presence of the lockdowns. Piovella provided a simplified analytical solution of the susceptible, exposed, infected, recovered (SEIR) model to predict the casualty peaks and asymptotic cases without iteratively solving the ordinary differential equations [14]. Even though the numeric and analytical solutions are close, there exist biases around the peak values. Piccolomini and Zama proposed a forced susceptible, exposed, infected, recovered, dead (fSEIRD) model with two different piecewise time-dependent infection rates [15]. It is stated that the model fits the data and makes reliable predictions for Italy. However, even though the SIR, SEIR, fSEIRD models are simple and require few parameters, they do not consider the pharmacological and non-pharmacological policies which play important roles on the dynamics of the pandemic diseases. In addition, they do not include the hospitalized, intensive care, and intubated pandemic casualties. Lee et al. modelled an optimal age specific vaccination policy against the H1N1 pandemic influenza in Mexico [16]. The model suggested that the optimal vaccination can be achieved by allocating more vaccines for the young adults age between 20 and 39. Recently, we developed a suspicious, infected, recovered (SpID) model with the second order difference equations rather than the first order ordinary differential equations as in the SIR, the SEIR, and the SEIRD models [17]. The results confirmed that the SpID model can represent the higher order properties such as the peak in the COVID-19 casualties. In our further research, we proposed a SpID-N model with the non-pharmacological policies (N) including the curfews, restrictions, and lockdowns [2]. The results highlighted the role of each non-pharmacological policy on the COVID-19 casualties. In addition, recently we performed a research to analyses the linear and non-linear dynamics of the COVID-19 by only considering the pharmacological policies [18]. This paper developed three model structures from linear to strongly non-linear and optimized their parameters with the mathematical optimization and machine learning approaches. As an alternative to the model-based control of the pandemic casualties, an artificial intelligence approach, which is implicitly a model free approach, was constructed to generate the multi-dimensional non-pharmacological policies [19]. This artificial intelligence algorithm allowed to weight each non-pharmacological policy together with each pandemic casualty under a certain vaccination policy. It firstly aimed to stabilize the pandemic casualties and then minimize them in time. Zhao et al. recently built an age-specific transmission model to quantify the transmissibility in different age groups [20]. Matrajt et al. developed an optimal vaccine allocation algorithm aiming at reducing the deaths, infections, and hospitalizations [21]. This paper proposes a susceptible (S c ), suspicious (S p ), infected (I n ), hospitalized (H), intensive care (I t ), intubated (I b) , recovered (R), death (D) with the priority and age specific vaccination (V) and non-pharmacological (N) policies (S c S p I n HI t I b RD-VN model). The key contributions of the paper can be summarized briefly as • A comprehensive S c S p I n HI t I b RD-VN model has been constructed by referring the known relationships among the COVID-19 casualties illustrated in Fig. 1. • Priority and age specific vaccination policy has been formulated and incorporated into the S c S p I n HI t I b RD-VN model together with the nonpharmacological policies. • Constrained recursive least squares (RLS) optimizer has been modified to learn the unknown parameters of the S c S p I n HI t I b RD-VN model by satisfying the structural and proportional contribution requirements of the design. • An extensive analysis has been performed to assess the role of the priority and age specific vaccination policy and the non-pharmacological policies. It is important to note that even though this paper mostly refers the COVID-19 pandemic, it can be implemented to all the pandemic diseases having the architecture shown in Fig. 1, which is constructed based on the epidemiological facts. In the rest of the paper, the proposed model structures, the proposed S c S p I n HI t I b RD-VN model, the constrained RLS for the multi-dimensional models, and the analysis of the model have been provided.

The proposed model architecture
Individuals in the susceptible S c k group are vulnerable to the pandemic diseases where the suspicious S p k ones leave the group (Fig. 1, number 1) and the non-infected I nn k ones re-join the susceptible group (number 5). The vaccinated (V k )(number 2), the recovered (R k ) (number 22), and the death (D k )(number 21) become non-susceptible S nc k (number 3) and leave the susceptible S c k group where the remained ones constitute the current susceptible S c k+1 (number 4) group. The individuals in the suspicious S p k group, who are tested and/or quarantined, either move to the infected I n k (number 6) group or the non-infected I nn k (number 7) group where some individuals in the infected I n k group can return the suspicious S p k group again. Also, since the infected I n k individuals spread the virus until they are isolated, they act like as an excitation signal I na k−1 (number 8) on the suspicious casualties. Individuals in the infected I n k group can be in the hospitalized (H k )(number 9) group or in the non-hospitalized H n k (number 10) group where the non-hospitalized (H k ) individuals join the recovered (R k )(number 15) group after a quarantine period. The individuals in the hospitalized (H k ) group can union with the intensive care I t k (number 11), the intubated I b k (number 12), the death (D k )(number 13), or the recovered (R k )(number 14) groups. The individuals in the intensive care I t k group can move to the intubated I b k (number 16), the death (D k )(number 17), or the hospitalized (H k )(number 18) groups. Similarly, the individuals in the intubated I b k group can join either the intensive care I t k (number 19), or the death (D k )(number 20) groups. The non-pharmacological policies (u k )(number 23) and priority and age specific vaccination policy V * k (number 24) act like an external inhibitor on all the casualties at varying rates.

The ScSpInHItIbRD-VN model
This section initially formulates the parametric sub-models, and then the vaccination and the non-pharmacological policies of the S c S p I n HI t I b RD-VN model.

The S c S p I n HI t I b RD-VN sub-models
This sub-section constructs the parametric models of each sub-model illustrated in Fig. 1.

The susceptible S c k sub-model
Considering the connections coming in and leaving out the susceptible S c k group in Fig. 1, one can formulate the S c k sub-model with a difference equation. We can initially write the difference equation of the non-susceptible S nc k group shown in Fig. 2 as where S nc k represents the non-susceptible individuals who have gained immunity and also the individuals who lost their lives, R k represents the recovered individuals, D k represents the dead individuals, V k represents the vaccinated individuals, a 14 ,a 15 ,c 1 are the unknown parameters.
The representation of the susceptible S c k+1 group in Fig. 2 is where S c k represents the individuals who may be infected and have a lack of immunity, S p k represents the suspicious individuals, I nn k represents the non-infected individuals, a 11 ,a 12 ,a 13 are the unknown parameters.Substituting Eq.
(1) in Eq. (2) yields All the parameters in Eq. (3) are unknown and will be learned from the available data with the RLS algorithm subject to the inequality constraints in the next section.
The next sub-section provides the modelling steps of the suspicious S p k sub-model.

The suspicious S p k sub-model
Some of the susceptible S c k individuals become suspicious S p k as they exhibit symptoms or contact an infected individual, or return from the regions where the pandemic disease is a threat. These individuals are either tested or quarantined for a time duration. In this paper, we define the suspicious S p k individuals as the number of the people tested daily. Therefore, the model can predict the number The susceptible S c k sub-model of the required tests in the future. We can represent the S p k sub-model shown in Fig. 3 as where I na k represents the individuals who can become suspicious again and excitation effect of the infected individuals on the suspicious casualties (related to filiation time), u k is the non-pharmacological policy, V S p k is the vaccination policy, a 21 , a 22 , a 23 , a 24 , a 25 b 2 , c 2 are the parameters.
The next sub-section presents the modelling steps of the infected I n k sub-model.

The infected I n k sub-model
Some of the suspicious S p k individuals becomes infected I n k where they either become hospitalized H k or non-hospitalized H n k , who are quarantined for a period of time, as illustrated in Fig. 1. We can formulate its model by considering the corresponding connections in Fig. 4 as where V I n k is the vaccination policy, a 31 , a 32 , a 33 , a 34 , a 35 b 3 , c 3 are the parameters.
The next sub-section introduces the hospitalized H k sub-model.

The hospitalized H k sub-model
Some of the infected I n k individuals requiring standard treatments join the hospitalized H k group. The hospitalized H k individuals can join the intensive care I t k , the intubated I b k , the recovered R k , or the death D k groups as shown in Fig. 5. We can formulate the hospitalized model as where V H k is the vaccination policy, a 41 , a 42 , a 43 , a 44 , a 45 , a 46 , b 4 , c 4 are the parameters, The next sub-section presents the formulation of the intensive care I t k sub-model.

The intensive care I t k sub-model
Some of the hospitalized H k individuals move to the intensive care I t k group where some of them move back to the hospitalized H k group as shown in Fig. 6. Similarly, some of the intensive care I t k patients become intubated I b k where some of them re-join the intensive care I t k group, and the rest join the death D k group. We can construct the intensive care I t k model as where V I t k is the vaccination policy, a 51 , a 52 , a 53 , b 5 , c 5 are the parameters.
The next sub-section provides the intubated I b k sub-model. Fig. 3 The suspicious S p k sub-model Fig. 4 The infected I n k sub-model Fig. 5 The hospitalized H k sub-model The intubated I b k sub-model Some of the hospitalized H k individuals and the intensive care I t k patients become intubated I b k as shown in Fig. 7. A number of the intubated I b k patients move back to the intensive care I t k unit while the rest join the death D k group. We can construct the intubated model as where V I b k is the vaccination policy, a 61 , a 62 , a 63 , a 64 b 6 , c 6 are the parameters.
The next sub-section formulates the recovered R k sub-model.

The recovered R k sub-model
A number of the hospitalized H k and the non-hospitalized H n k individuals join the recovered R k group who become non-susceptible S nc k as illustrated in Fig. 8. We can formulate the recovered R k sub-model as where V R k is the vaccination policy, a 71 , a 72 , a 73 , a 74 b 7 , c 7 are the parameters.
The next sub-section expresses the death sub-model.

The death D k sub-model
Some of the hospitalized H k , the intensive care I t k , and the intubated I b k individuals join the death D k group and become non-susceptible S nc k as illustrated in Fig. 9. We can form the death D k model as where V D k is the vaccination policy, a 81 , a 82 , a 83 , a 84 , a 85 b 8 , c 8 are the parameters.
The next sub-section formulates the vaccination policy V * k and reviews the non-pharmacological u k policies.

The vaccination V * k and non-pharmacological u k policies
This section firstly introduces the priority and age specific vaccination policy V * k and reviews the non-pharmacological policies u k that we have developed recently for the first time in the literature [2]. Fig. 6 The intensive care I t k sub-model Fig. 7 The intubated I b k sub-model Fig. 8 The recovered R k sub-model Fig. 9 The death D k sub-model The priority and age specific vaccination policies V * k The * in the priority and age specific vaccination policy V * k represents the S c ,S p ,I n ,H,I t ,R , and D in the sub-models given by Eqs. from (4) to (10). The priority and age specific vaccination policy basis V b k is defined in terms of the number of the daily vaccinated people in each group as where H s k is the healthcare staff, A 80+ Similarly, we can construct the priority and age specific vaccination policy V * k for the other sub-models by following the same steps introduced in this section. The next sub-section provides the revised non-pharmacological policies u k .

The non-pharmacological policies u k
The authorities impose various curfews and restrictions to confine the spread of the virus. The most common ones are the curfews on the people age over 65, age under 20, and people with the chronic diseases which have been parametrized in [2] (since there is no available data) as (11) where u s k is the response of the curfew (in closed form solution), n s is the number of the people under the curfew, k is the number of the days and k i is the start day of the curfew, α is the discount factor of the response, where α k ≈ 0 for α = 0.71 and k = 14 (quarantine duration), σ s k is the random non-parametric uncertainty in the response.
The other common precaution is the curfews on the weekends and holidays, which has a transient ascent part as where u wh i,k is the response of the curfews on the weekends and holidays, n wh is the number of the people under the curfews on the weekends and holidays, σ wh i,k is the random uncertainty in the response.
Its transient descent part is modelled as The overall response u wh k is In terms of the closure of the schools and universities, it is not a curfew as it only hinders mass gatherings of the students; hence, they can come together in smaller groups. Therefore, the response has a transient ascent part as in Eq. (15) and transient descent part as in Eq. (16). These parts are essentially for removing the negative impacts of the schools being open. Then an uncertain saturated part u sat represents the small gatherings after the closure of the schools. After the transient ascent and descent parts, the saturated part can be represented as where n su is the number of the students, σ su k is the random uncertainty in the response, k n = k su i + k su n where k su i is the start day and k su n is the duration of the closure.

Comparison of the prediction models
One can summarize the main advantages of the constructed S c S p I n HI t I b RD-VN model over the well-known models such as the SIR, SEIR models in terms of the solution and analysis as • It has difference equations rather than the differential equations. Therefore, it can be solved iteratively without requiring an ordinary differential equation solver. • It has coupled and linear dynamics instead of the slightly coupled nonlinear dynamics. Thus, the mathematical analysis of the parametric model is straightforward. • Its unknown parameters are learned from the reported data by using the well-known multi-dimensional optimization approaches rather than the single dimensional statistical approaches.
The next section forms the RLS approach with the inequality constraints to learn the unknown parameters of the S c S p I n HI t I b RD-VN model.

The constrained RLS algorithm
In this paper, the constrained optimization is considered for two reasons: The first one is that the sub-models have certain parameter structures together with the corresponding parameter signs and the second reason is to reflect the contributions of the data having huge magnitude differences (for example, while the susceptible S c k group covers millions of the individuals, the hospitalized H k group covers only thousands of them). In this section, initially we will divide the optimization problem in terms of the estimated sub-model casualties (outputs) and the real casualties. Then, the RLS algorithm with the inequality constraints are modified to learn the unknown parameters.

The estimated sub-models
We can represent the estimated sub-models ŷ * k in terms of the known basis vector b * k and the unknown parameter vector w * k , where the * is denoted for the S c ,S p ,I n ,H,I t ,R , and D in the sub-models given by Eqs. from (3) to (10) as For example, the basis b S c k of the estimated susceptible ŷ S c k sub-model is formed with respect to the left hand side of Eq. (3) as (19) And the corresponding unknown parameter vector w S c k of the estimated susceptible ŷ S c k sub-model with respect to the right hand side of Eq. (3) is The other estimated sub-models, their bases and parameter vectors are formed by following the same procedures as in Eqs. (19), (20), and (21), respectively. The next sub-section introduces the modified RLS algorithm with the inequality constraints to learn the unknown parameter vectors w * k .

Learning the Unknown Parameters with the Constrained RLS
The reported casualties are the outputs of the S c S p I n HI t I b RD-VN sub-models and we call them as the real outputs y * k . For example, the real output of the susceptible sub-model is the left hand side of Eq. (3), which is S c k+1 . The objective function is constructed with the 2-norm of the instant estimation error defined as where α is the inequality constraints which are the lower bound of the parameters. We can construct the Lagrange multipliers used for solving the optimization problems as Getting partial derivative of L w * k , with respect to the w * k yields Getting partial derivative of L w * k , with respect to the gives Re-organizing Eq. (24) as w * k is on the left and the rest are on the right, and then substituting it in Eq. (25) yields The Lagrange multiplier from Eq. (26) is obtained as (21) w S p k = a 11 a 12 a 13 a 14 a 15 c 1 Then by reinserting Eq. (27) into Eq. (24), the unknown parameter vector w * k can be attained. The next section extensively analyses the S c S p I n HI t I b RD-VN model.

Results
This section initially presents the parameters of the proposed S c S p I n HI t I b RD-VN model and then analyses the training and prediction results. Table 1 provides the parameters of the model. The next sub-section compares the real and estimated COVID-19 casualties with the constrained RLS algorithm. Figure 10 shows the real (reported) and estimated casualties for Turkey.

Real and estimated casualties
As can be seen from Fig. 10, the estimated casualties with the S c S p I n HI t I b RD-VN model closely follow the real casualties. The model can track the steep peaks and also the daily variations in the casualties even though the constructed parameter spaces are limited (in machine learning approaches, we randomly manipulate the parameter spaces until we have close estimations). The casualties in Fig. 10 have two distinctive peaks and estimated future casualties in Fig. 13 shows the third peak, which will occur in 40 days. It is clear that the susceptible S c k casualties have noticeable reduction with the initiation of the vaccination process. It seems that this vaccination process has affected the other casualties since they sharply decrease as well. The decrease in the casualties has also been supported with the non-pharmacological policies u k . Figure 11 shows the mean errors and the corresponding standard deviations in the estimates.
As illustrated by Fig. 11, even though all the mean errors are small, the standard deviations are quite large. This is due to existence of the steep peaks shown in Fig. 10. These peaks occurred in December when there were not any active pharmacological and non-pharmacological policies. This implies that the character of the casualties is largely shaped based on the external impacts such as the pharmacological and non-pharmacological policies. Our recent work highlighted that the pharmacological and non-pharmacological policies have damping impact on the casualties whereas they also have natural frequency determined by the internal and coupling dynamics. The next sub-section presents the priority and age specific vaccination policy results. Figure 12 shows the priority and age specification vaccination policy for the hospitalized V H k = w H T k V b k and death V D k = w D T k V b k sub-models.   Table 2 provides the background parameters in Fig. 12.

Priority and age specific vaccination policy
Since there is no reported data for the healthcare staff H s ('Staff ' in Table 1), the largest values are assigned for them as they are in the highest risk group. Therefore, the vaccination of the healthcare staff has the largest hospitalized w H k and death w D k policy values ( have the smallest population among the all age groups, they have the smallest regions on the horizontal axis representing the number of days k. Figure 13 shows the future casualties estimated by the S c S p I n HI t I b D-VN model.

The estimated future casualties
The future estimates are obtained under the assumption that 100.000 peoople are vaccinated daily. It is clear that the number of the susceptible S c k individuals reduces almost linearly since the vaccinated V k individuals leave the group. Since the susceptible S c k and the suspicious S p k groups are strongly coupled, reduction in the susceptible S c k group is reflected onto the suspicious S p k group as well. Figure 13 also clearly shows the third peak in the COVID-19 casualties. In addition, it is noticaeble that even though the suspicious S p k , the infected I n k , and the hospitalized H k converge to zero in 500 days, the intensive care I t k , the intubated I b k , and the death D k converge  to zero around 120 days. This fast convergence is due to priority and age specific vaccination policy which focus on vaccination of the people in the high risk groups. Figure 14 shows the number of the daily vaccinations and the corresponding average future casualties. As can be clearly seen from Fig. 14, all the COVID-19 casualties reduce depending on the number of the daily vaccinations. The largest reductions occur in the number of the intensive care I t k , the intubated I b k , and the death D k when the number of the vaccination rises from 50.000 to 100.000. The further noticeable reduction occurs in the number of the suspicious S p k and the infected I n k when the daily vaccination number rises from 100.000 to 200.000. Figure 15 shows the role of the non-pharmacological policies on the casualties.

Analysis of the non-pharmacological policies
As can be seen from Fig. 15 when all the non-pharmacological policies are in place, all the casualties are small and they increase when the curfews are lifted (blue bar). Removing the restrictions on the people age over 65 and people with the chronic diseases has limited effects as they are in the priority group and most of them have been already vaccinated (orange bar). With respect to the partial opening of the schools, since the majority of the students are not attending the schools, their impact is bounded as well (yellow bar). However, curfews on the weekends, holidays, and nights cover the whole population; henceforth, their roles on the casualties are distinctive (purple bar).

Conclusions
This paper developed a comprehensive parametric S p S c I n HI t I b RD-VN model to analyse and estimate the role of the priority and age specific vaccination policy and the non-pharmacological policies. The model has a structure constructed by using the key insights about the pandemic diseases. To satisfy the model structural requirements and avoid dominant effects of the large susceptible and suspicious data, a constrained RLS algorithm has been formed. The results clearly show the importance of the priority and age specific vaccination policy on all the casualties. The future hospitalized, intensive care, intubated, and death casualties converge to zero before the other casualties since they have larger importance. However, the future susceptible, suspicious, infected, and recovered casualties are large due to the people in the lower risk groups are not vaccinated yet. The paper also addresses the relationships among the various daily vaccinations, non-pharmacological polices, and the corresponding COVID-19 casualties. The results confirm that the curfews on the weekends and holidays has an overwhelming role on reducing the casualties.

Limitations of the work
Effects of the non-pharmacological policies on each age and chronic diseases group are not weighted. Moreover, vaccine effectiveness for each age group has not been added the model. Besides, climate and environmental effects are not considered.

Future works
The non-pharmacological policies of the S c S p I n HI t I b RD-VN model should be also modified to consider the priority and age specific impacts on each casualty. In addition, effectiveness of the different brand of vaccines on each age group should be considered in the model. Moreover, the model can be expanded by considering the unknown uncertainties. Finally, a toolbox should be constructed and provided for free to help the researchers when applying the proposed model.