Mathematical modelling plays a significant role in understanding the complexities of infectious diseases, which is important for disease control and prediction. Severe acute respiratory syndrome (SARS), as a severe infectious disease, was rapidly spread in more than 28 countries, infected 8098 people and caused 774 deaths worldwide in late 2002 and early 2003. SARS is an acute respiratory illness caused by infection by the SARS virus and the key signs and symptoms are fever, respiratory compromise, chilis, muscle aches and headache. The etiological agent of SARS is a coronavirus, which could be spread predominantly by droplet and by direct or indirect contacts. Beijing, as the capital city of China with population 17,800,000 was one of the major cities suffer from the disease with 2434 people infected and 147 deaths. Several models were created to study the outbreak of SARS including compartmental models, small world network models, system dynamic models.
In this project, we want to build mathematic model to investigate the transmission of SARS in Beijing, China and determined the basic reproductive number for this model. The data will be analyzed with comaprtmental model.
Chinese government update the SARS data every day online from April 19th to August 16th, 2003. After August 16th, these is no new infection case found in China and it has been reported that the SARS virus was fully eliminated on September 2nd worldwide. The data used in this project was pulled from the everyday report of Ministry of Health of the People’s Republic of China from China Education and Research Network website www.edu.cn . The dataset contains the number of people infected and dead from SARS in Beijing and other parts of China.
At first, we will look at the data about infection, recovery and death cases from April 19th to August 16th in Beijing. The cases were recorded daily with no missing data.
From the above plots, there very few new SARS infection cases found in Beijing after 60 days. Most people were recovered from the infection before 60 days too. Data after 60 days are mainly about very few people slowly recovered from the infection in hospital. So, data after 60 days are not considered in this model as it may not be important for the model.
Our model has the earmarks of susceptible-infectious-recovered(SIR) model. As the signs of SARS infection is clear and easy to identify, there is no death directly coming from susceptible group.
Model analysis was conducted using the POMP framework. POMP is useful for non-linear processes, which is useful for epidemiological system. POMP gives more flexibility in model creation, which is useful for complex process such as the spread of a disease. The model codes used in this project are heavily borrowed from Note 12 as a similar model structure was used. The model was coded with R and processed under FLUX. The original R code and results and creates from FLUX were included, while R-markdown file was created with directly linking the results and figures created by R-code.
At first, several parameters were chosen to compare simulation results with the data. The simulated result is close to data, which will be used in the future. \[ \beta = 2 \] \[ \mu_I= 0.04 \] \[ \rho = 0.3\] \[ \mu_R1 = 0.05\] \[mu_R2 = 0.1\] \[mu_R3 = 0.1\] The infection vs date plot was shown below:
## [1] "B" "C" "K"
Our initial step was to perform particle filtering on these maximum likelihood estimates in order to find a suitable set of initial parameters. As we do not have a good filter available due the access of references, a arbitory filter was created based the flue exmaple in the course. As it’s not a real data obtained previous studies, the local search using the iterated filtering algorithm specified above does not work well.
With the maximum loklik parameters, the simulation is not close to the real data.
Next a global search using the same iterated filtering algorithm was done. This time, the initial parameters were randomly and uniformly selected from a specified range by mif2 algorithm.
The parameters range is shown below:
\[ \beta = c(0.001, 0.1) \] \[ \mu_I= c(0.5, 10) \] \[ \rho = (0.5, 1)\]
The parameters with maximum liklihood is shown as below:
\[logLik = -708.9\]
\[LogLik se = 0.58\]
\[ \beta = 0.007 \] \[ \mu_I= 65.3 \] \[ \rho = 0.935\] \[\mu_R1 = 0.044\]
\[\mu_R2 = 0.063\]
\[\mu_R3 = 0.355\]
From the simulation with the maximum liklihood parameters from the global search, the result is better than the local search, which means the arbitory parameters are far from the real parameters. The filter parameters should be updated with the global search results. And the parameter range should be changed in future study.
The following are the mif2 maximization procedure plots. From the plots, we can see that \(\beta\) and \(\mu_I\) fail to converge during the run. It might be becuase the parameter range is not good or the filter is not good. We did see \(\rho\) converges to 0.2 which suggest that a low reporting rate. After a few iterations, the mif2 likelihood displays a slow downward trend. This is because mif2 adds extra variability to the model in the form of random perturbations to the parameters. to allow one model to converge to the other.
In this project, we proposed mathematical model for SARS infection and recovery in Beijing with pomp process. Model parameters were determined from local and global search. But as we do not have previous data as filters, so next step we need to update the filters with our results from global search and rerun the models. Aother future direction is that we might change the model structures. We can add high risk group to the model, as they have close contaction with infected patients. In the goverment reports, these high risk group of people will be separated for about 2 weeks and about 20-30% of these were found to be infected. The new structure model was shown below: