Ten Years Graph
First Five Years Graph
Second Five Years Graph
Test R code for a single SIR epidemic
Test results for a single SIR epidemic
Test results for a single SIR epidemic (with MLE2 method)
Single SIR model for 1989-1990 epidemic
Double SIR Epidemic Model
Weston's R code for the two SIR model
Weston's R code for the two SIR model (with MLE2 method)
Triple SIR Epidemic Model
Double SIR model for 1989-1990 epidemic with "BFGS" method (different alpha values)
Double SIR model for 1989-1990 epidemic with "SANN" method (different alpha values)
Double SIR model for 1989-1990 epidemic with "Nelder-Mead" method (different alpha values)
Double SIR model for 1989-1990 epidemic with "BFGS" method (same alpha values)
AICc to compare the "BFGS" method model with the same and different alpha values
Double SIR model for 1989-1990 epidemic with "BFGS" method (with alpha values equal to 1/7) (Junling: Note that alpha=1/7 means the infectious period is 7 weeks, not 7 days, which means you are not likely to get reasonable fit, as illustrated by your results)
From the MLE2 Method
Double SIR model for 1989-1990 epidemic with "BFGS" method (different alpha values, mle2)
Three SIR model for 1989-1990 epidemic with "BFGS" method (different alpha values, mle2)
Double SIR model for 1992-1993 epidemic with "BFGS" method (different alpha values)
Double SIR model for 1992-1993 epidemic with "BFGS" method (same alpha values)
First eight of the models (1990-1991 same except different data)
1989-1990 Epidemic 2 SIR model baseline Different Alpha BFGS
1989-1990 Epidemic 2 SIR model baseline Different Alpha NM
1989-1990 Epidemic 2 SIR model baseline Same Alpha BFGS
1989-1990 Epidemic 2 SIR model baseline Same Alpha NM
1989-1990 Epidemic 3 SIR model baseline Different Alpha BFGS
1989-1990 Epidemic 3 SIR model baseline Different Alpha NM
This code restricts all of the parameters to be positive, as well as having N be larger than the total, otherwise it will return a small likelihood (1e10).
I tried using 1/r as the parameter and the mle2 optimizer was not very happy with it, many times giving the response
Error in .local(fitted, …) :
Hessian is ill-behaved or missing,can't find an initial estimate of std. error(consider specifying std.err in profile call)
Error in .local(object, parm, level, …) :
Problem with profiling: Error in .local(fitted, ...) : Hessian is ill-behaved or missing,can't find an initial estimate of std. error(consider specifying std.err in profile call)
Calls: SIR_fit → confint → confint → .local Execution halted
1989-1990 Epidemic 2 SIR model baseline Different Alpha BFGS Jul 9 (worked)
1989-1990 Epidemic 2 SIR model baseline Different Alpha NM Jul 9 (did not work, “Execution halted”)
1989-1990 Epidemic 2 SIR model baseline Same Alpha BFGS Jul 9 (did not work, “Execution halted”)
1989-1990 Epidemic 2 SIR model baseline Same Alpha NM Jul 9 (did not work, “Execution halted”)
1989-1990 Epidemic 3 SIR model baseline Different Alpha BFGS Jul 9 (did not work, “Execution halted”)
1989-1990 Epidemic 3 SIR model baseline Different Alpha NM Jul 9 (did not work, “Execution halted”)
1989-1990 Epidemic 3 SIR model baseline Same Alpha BFGS Jul 9 (did not work, “Execution halted”)
1989-1990 Epidemic 3 SIR model baseline Same Alpha NM Jul 9 (did not work, “Execution halted”)
I think all of these other models did not work because I changed the starting values to reflect the good values that seemed to be working with the
“1989-1990 Epidemic 2 SIR model baseline Different Alpha BFGS Jul 9”.
Additionally, I started the parameter r at for all the others 760, whereas I started r at 2 for the model that worked. I think it's probably better to let the optimizer “find” the value of r, than make it so high initially.
SIR equations written in C code
Test code in R to make sure the SIR equations in C returned the same values as if they were written in R.
Test code for SIR Equations in C
The code is correct as evidenced by the results of the test code.
Results of the Test Code for SIR Equations in C
Method for Even Faster C Code (ode solver also in C using gsl)
SIR Equations and ODE Solver in C
Test code for SIR Equations and ODE Solver in C
Unfortunately, the above test code does not work. The same error is repeated
Error in dyn.load(“SIR_solver.so”) :
unable to load shared library '/home1e/wcr/SIR_solver.so': /home1e/wcr/SIR_solver.so: undefined symbol: gsl_odeiv_step_rkf45
Execution halted
It appears that other people have had this problem too with linking the gls library from C code to R code.
https://stat.ethz.ch/pipermail/r-devel/2006-December/043927.html
http://www.mail-archive.com/r-help@r-project.org/msg76240.html
The PDF below for linking is good but complicated and in C++ rather than C. http://www.google.ca/#hl=en&source=hp&q=%24+PKG_LIBS%3D%22-lgsl+-lblas%22&btnG=Google+Search&aq=f&aqi=&aql=&oq=%24+PKG_LIBS%3D%22-lgsl+-lblas%22&gs_rfai=&fp=77a98ce8ed84c158
I tried
$ PKG_LIBS = “-L/user/local/lib/libgsl -lgsl” R CMD SHLIB SIR_solver.c
on the command line and nothing happened. I know that gsl is installed on the “euler”.
Thanks for your help! Wes
Result with Poisson distribution after r > 10,000:
Result from 1989-1990 epidemic 2 SIR model Baseline Different Alpha BFGS Jul 9
Code for the model with only the Negative Binomial option:
Code for 1989-1990 Epidemic 2 SIR model Baseline Different Alpha BFGS Neg Bin
Result with only Negative Binomial option:
Result 1989-1990 Epidemic 2 SIR model Baseline Different Alpha BFGS Neg Bin
July 16th
Code for the model using C code with Poisson distribution after r > 10,000 (same code as Jul 9th model with new faster C code):
Code for 1989-1990 Epidemic 2 SIR model Baseline Different Alpha BFGS (poi/neg C code)
Result using C code with Poisson distribution after r > 10,000:
Result 1989-1990 epidemic 2 SIR model Baseline Different Alpha BFGS Jul (poi/neg C code)
July 19th (very quick to finish 3 and 1/2 hours versus 3-4 days with July 16th code)
Code for the model using lsoda, Jacobian, and C code with Poisson distribution after r > 1,000,000:
Code for 1989-90 Epidemic 2 SIR model Baseline Different Alpha BFGS (poi/neg fast C code)
Result using very fast C code with Poisson distribution after r > 1,000,000:
Result 1989-1990 epidemic 2 SIR model Baseline Different Alpha BFGS Jul (poi/neg fast C code)
Code for the model using lsoda, Jacobian, and C code with only Neg. Bin. distribution:
Code for 1989-90 Epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result using very fast C code with only Neg. Bin. distribution:
Result 1989-1990 epidemic 2 SIR model Baseline Different Alpha BFGS Jul (neg fast C code)
AICc value 175.1937
Best so far
July 21st
Code for 1989-90 Epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1989-1990 epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 169.8512
July 28th
NA value(s) in the CI
Code for 1989-90 Epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1989-1990 epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 182.0587
Best so far
July 21st
Code for 1990-91 Epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1990-91 epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 190.7815
July 26th
Code for 1990-91 Epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1990-91 epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 210.1084
Best so far
July 21st
Code for 1991-92 Epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1991-92 epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 189.7440
July 28th
NA value(s) in CI
Code for 1991-92 Epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1991-92 epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 197.7814
July 28th
NA value(s) in CI
Code for 1991-92 epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1991-92 epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 196.0226
Best so far
July 22nd
NA value(s) in the CI
Code for 1992-93 Epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1992-93 epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 202.8748
July 22nd
Code for 1992-93 Epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1992-93 epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 207.8083
July 22nd
Code for 1995-96 Epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1995-96 epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 275.9575
July 28th
NA value(s) in the CI
Code for 1995-96 epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1995-96 epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 275.9059
Best so far
July 25th
Code for 1995-96 Epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1995-96 epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 270.8821
July 22nd
Code for 1996-97 Epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1996-97 epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 352.8229
Best so far
July 25th
Code for 1996-97 Epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1996-97 epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 342.7335
July 22nd
Code for 1997-98 Epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1997-98 epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 327.6377
Best
July 28th
NA value(s) in the CI
Code for 1997-98 Epidemic 2 SIR model Baseline Same Alpha SANN/BFGS (neg fast C code)
Result 1997-98 Epidemic 2 SIR model Baseline Same Alpha SANN/BFGS (neg fast C code)
AICc value 325.0415
July 23rd
Code for 1997-98 Epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1997-98 epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 326.4857
July 25th
Code for 1997-98 Epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1997-98 epidemic 3 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 330.13
July 23rd
Code for 1998-99 Epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1998-99 epidemic 2 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 308.5382
Best so far
July 30th
NA value(s) in the CI
Code for 1998-99 Epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
Result 1998-99 epidemic 2 SIR model Baseline Same Alpha BFGS (neg fast C code)
AICc value 307.3100
July 23rd
Code for 1998-99 Epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
Result 1998-99 epidemic 3 SIR model Baseline Different Alpha BFGS (neg fast C code)
AICc value 323.2764