Iterated filtering

Iterated filtering algorithms are a tool for maximum likelihood inference on partially observed dynamical systems. Stochastic perturbations to the unknown parameters are used to explore the parameter space. Applying sequential Monte Carlo (the particle filter) to this extended model results in the selection of the parameter values that are more consistent with the data. Appropriately constructed procedures, iterating with successively diminished perturbations, converge to the maximum likelihood estimate. ^{[1] [2] [3]} Iterated filtering methods have so far been used most extensively to study infectious disease transmission dynamics. Case studies include cholera, ^{[4] [5]} Ebola virus, ^[6] influenza, ^{[7] [8] [9] [10]} malaria, ^{[11] [12] [13]} HIV, ^[14] pertussis, ^{[15] [16]} poliovirus ^[17] and measles. ^{[5] [18]} Other areas which have been proposed to be suitable for these methods include ecological dynamics ^{[19] [20]} and finance. ^{[21] [22]}

The perturbations to the parameter space play several different roles. Firstly, they smooth out the likelihood surface, enabling the algorithm to overcome small-scale features of the likelihood during early stages of the global search. Secondly, Monte Carlo variation allows the search to escape from local minima. Thirdly, the iterated filtering update uses the perturbed parameter values to construct an approximation to the derivative of the log likelihood even though this quantity is not typically available in closed form. Fourthly, the parameter perturbations help to overcome numerical difficulties that can arise during sequential Monte Carlo.

Overview

The data are a time series ${\displaystyle y_{1},\dots ,y_{N}}$ collected at times ${\displaystyle t_{1}<t_{2}<\dots <t_{N}}$. The dynamic system is modeled by a Markov process ${\displaystyle X(t)}$ which is generated by a function ${\displaystyle f(x,s,t,\theta ,W)}$ in the sense that

${\displaystyle X(t_{n}^{})=f(X(t_{n-1}^{}),t_{n-1}^{},t_{n}^{},\theta ,W)}$

where ${\displaystyle \theta }$ is a vector of unknown parameters and ${\displaystyle W}$ is some random quantity that is drawn independently each time ${\displaystyle f(.)}$ is evaluated. An initial condition ${\displaystyle X(t_{0})}$ at some time ${\displaystyle t_{0}<t_{1}}$ is specified by an initialization function, ${\displaystyle X(t_{0})=h(\theta )}$. A measurement density ${\displaystyle g(y_{n}|X_{n},t_{n},\theta )}$ completes the specification of a partially observed Markov process. We present a basic iterated filtering algorithm (IF1) ^{[1] [2]} followed by an iterated filtering algorithm implementing an iterated, perturbed Bayes map (IF2). ^{[3] [23]}

Procedure: Iterated filtering (IF1)

Input: A partially observed Markov model specified as above; Monte Carlo sample size ${\displaystyle J}$; number of iterations ${\displaystyle M}$; cooling parameters ${\displaystyle 0<a<1}$ and ${\displaystyle b}$; covariance matrix ${\displaystyle \Phi }$; initial parameter vector ${\displaystyle \theta ^{(1)}}$

for ${\displaystyle m_{}^{}=1}$ to ${\displaystyle M_{}^{}}$

draw ${\displaystyle \Theta _{F}(t_{0}^{},j)\sim \mathrm {Normal} (\theta ^{(m)},ba^{m-1}\Phi )}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle X_{F}(t_{0}^{},j)=h{\big (}\Theta _{F}(t_{0}^{},j){\big )}}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle {\bar {\theta }}(t_{0}^{})=\theta ^{(m)}}$

for ${\displaystyle n_{}^{}=1}$ to ${\displaystyle N_{}^{}}$

draw ${\displaystyle \Theta _{P}(t_{n}^{},j)\sim \mathrm {Normal} (\Theta _{F}(t_{n-1}^{},j),a^{m-1}\Phi )}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle X_{P}(t_{n}^{},j)=f(X_{F}(t_{n-1}^{},j),t_{n-1}^{},t_{n},\Theta _{P}(t_{n},j),W)}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle w(n,j)=g(y_{n}|X_{P}(t_{n}^{},j),t_{n}^{},\Theta _{P}(t_{n},j))}$ for ${\displaystyle j=1,\dots ,J}$

draw ${\displaystyle k_{1}^{},\dots ,k_{J}^{}}$ such that ${\displaystyle P(k_{j}^{}=i)=w(n,i){\big /}{\sum }_{\ell }w(n,\ell )}$

set ${\displaystyle X_{F}(t_{n}^{},j)=X_{P}(t_{n}^{},k_{j}^{})}$ and ${\displaystyle \Theta _{F}(t_{n}^{},j)=\Theta _{P}(t_{n}^{},k_{j}^{})}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle {\bar {\theta }}_{i}^{}(t_{n}^{})}$ to the sample mean of ${\displaystyle \{\Theta _{F,i}^{}(t_{n}^{},j),j=1,\dots ,J\}}$, where the vector ${\displaystyle \Theta _{F}^{}}$ has components ${\displaystyle \{\Theta _{F,i}^{}\}}$

set ${\displaystyle V_{i}^{}(t_{n}^{})}$ to the sample variance of ${\displaystyle \{\Theta _{P,i}^{}(t_{n}^{},j),j=1,\dots ,J\}}$

set ${\displaystyle \theta _{i}^{(m+1)}=\theta _{i}^{(m)}+V_{i}(t_{1})\sum _{n=1}^{N}V_{i}^{-1}(t_{n})({\bar {\theta }}_{i}(t_{n})-{\bar {\theta }}_{i}(t_{n-1}))}$

Output: Maximum likelihood estimate ${\displaystyle {\hat {\theta }}=\theta ^{(M+1)}}$

Variations

1. For IF1, parameters which enter the model only in the specification of the initial condition, ${\displaystyle X(t_{0})}$, warrant some special algorithmic attention since information about them in the data may be concentrated in a small part of the time series. ^[1]
2. Theoretically, any distribution with the requisite mean and variance could be used in place of the normal distribution. It is standard to use the normal distribution and to reparameterise to remove constraints on the possible values of the parameters.
3. Modifications to the IF1 algorithm have been proposed to give superior asymptotic performance. ^{[24] [25]}

Procedure: Iterated filtering (IF2)

Input: A partially observed Markov model specified as above; Monte Carlo sample size ${\displaystyle J}$; number of iterations ${\displaystyle M}$; cooling parameter ${\displaystyle 0<a<1}$; covariance matrix ${\displaystyle \Phi }$; initial parameter vectors ${\displaystyle \{\Theta _{j},j=1,\dots ,J\}}$

for ${\displaystyle m_{}^{}=1}$ to ${\displaystyle M_{}^{}}$

set ${\displaystyle \Theta _{F}(t_{0}^{},j)\sim \mathrm {Normal} (\Theta _{j},a^{m-1}\Phi )}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle X_{F}(t_{0}^{},j)=h{\big (}\Theta _{F}(t_{0}^{},j){\big )}}$ for ${\displaystyle j=1,\dots ,J}$

for ${\displaystyle n_{}^{}=1}$ to ${\displaystyle N_{}^{}}$

draw ${\displaystyle \Theta _{P}(t_{n}^{},j)\sim \mathrm {Normal} (\Theta _{F}(t_{n-1}^{},k_{j}^{}),a^{m-1}\Phi )}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle X_{P}(t_{n}^{},j)=f(X_{F}(t_{n-1}^{},j),t_{n-1}^{},t_{n},\Theta _{P}(t_{n},j),W)}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle w(n,j)=g(y_{n}|X_{P}(t_{n}^{},j),t_{n}^{},\Theta _{P}(t_{n},j))}$ for ${\displaystyle j=1,\dots ,J}$

draw ${\displaystyle k_{1}^{},\dots ,k_{J}^{}}$ such that ${\displaystyle P(k_{j}^{}=i)=w(n,i){\big /}{\sum }_{\ell }w(n,\ell )}$

set ${\displaystyle X_{F}(t_{n}^{},j)=X_{P}(t_{n}^{},k_{j}^{})}$ and ${\displaystyle \Theta _{F}(t_{n}^{},j)=\Theta _{P}(t_{n}^{},k_{j}^{})}$ for ${\displaystyle j=1,\dots ,J}$

set ${\displaystyle \Theta _{j}=\Theta _{F}(t_{N}^{},j)}$ for ${\displaystyle j=1,\dots ,J}$

Output: Parameter vectors approximating the maximum likelihood estimate, ${\displaystyle \{\Theta _{j},j=1,\dots ,J\}}$

Software

"pomp: statistical inference for partially-observed Markov processes"  : R package.

References

1. Ionides, E. L.; Breto, C.; King, A. A. (2006). "Inference for nonlinear dynamical systems". Proceedings of the National Academy of Sciences of the USA. 103 (49): 18438–18443. Bibcode:2006PNAS..10318438I. doi: 10.1073/pnas.0603181103 . PMC   3020138 . PMID   17121996.
2. Ionides, E. L.; Bhadra, A.; Atchade, Y.; King, A. A. (2011). "Iterated filtering". Annals of Statistics. 39 (3): 1776–1802. arXiv: 0902.0347 . doi:10.1214/11-AOS886. S2CID   6527480.
3. Ionides, E. L.; Nguyen, D.; Atchadé, Y.; Stoev, S.; King, A. A. (2015). "Inference for dynamic and latent variable models via iterated, perturbed Bayes maps". Proceedings of the National Academy of Sciences of the USA. 112 (3): 719–724. Bibcode:2015PNAS..112..719I. doi: 10.1073/pnas.1410597112 . PMC   4311819 . PMID   25568084.
4. King, A. A.; Ionides, E. L.; Pascual, M.; Bouma, M. J. (2008). "Inapparent infections and cholera dynamics" (PDF). Nature. 454 (7206): 877–880. Bibcode:2008Natur.454..877K. doi:10.1038/nature07084. hdl: 2027.42/62519 . PMID   18704085. S2CID   4408759.
5. Breto, C.; He, D.; Ionides, E. L.; King, A. A. (2009). "Time series analysis via mechanistic models". Annals of Applied Statistics. 3: 319–348. arXiv: 0802.0021 . doi:10.1214/08-AOAS201. S2CID   8400632.
6. King AA, Domenech de Celles M, Magpantay FM, Rohani P (2015). "Avoidable errors in the modelling of outbreaks of emerging pathogens, with special reference to Ebola". Proceedings of the Royal Society B. 282 (1806): 20150347. doi:10.1098/rspb.2015.0347. PMC   4426634 . PMID   25833863.
7. He, D.; J. Dushoff; T. Day; J. Ma; D. Earn (2011). "Mechanistic modelling of the three waves of the 1918 influenza pandemic". Theoretical Ecology. 4 (2): 1–6. doi:10.1007/s12080-011-0123-3. S2CID   2010776.
8. Camacho, A.; S. Ballesteros; A. L. Graham; R. Carrat; O. Ratmann; B. Cazelles (2011). "Explaining rapid reinfections in multiple-wave influenza outbreaks: Tristan da Cunha 1971 epidemic as a case study". Proceedings of the Royal Society B. 278 (1725): 3635–3643. doi:10.1098/rspb.2011.0300. PMC   3203494 . PMID   21525058.
9. Earn, D.; He, D.; Loeb, M. B.; Fonseca, K.; Lee, B. E.; Dushoff, J. (2012). "Effects of School Closure on Incidence of Pandemic Influenza in Alberta, Canada". Annals of Internal Medicine. 156 (3): 173–181. doi: 10.7326/0003-4819-156-3-201202070-00005 . PMID   22312137.
10. Shrestha, S.; Foxman, B.; Weinberger, D. M.; Steiner, C.; Viboud, C.; Rohani, P. (2013). "Identifying the interaction between influenza and pneumococcal pneumonia using incidence data". Science Translational Medicine. 5 (191): 191ra84. doi:10.1126/scitranslmed.3005982. PMC   4178309 . PMID   23803706.
11. Laneri, K.; A. Bhadra; E. L. Ionides; M. Bouma; R. C. Dhiman; R. S. Yadav; M. Pascual (2010). "Forcing versus feedback: Epidemic malaria and monsoon rains in NW India". PLOS Computational Biology. 6 (9): e1000898. Bibcode:2010PLSCB...6E0898L. doi: 10.1371/journal.pcbi.1000898 . PMC   2932675 . PMID   20824122.
12. Bhadra, A.; E. L. Ionides; K. Laneri; M. Bouma; R. C. Dhiman; M. Pascual (2011). "Malaria in Northwest India: Data analysis via partially observed stochastic differential equation models driven by Lévy noise". Journal of the American Statistical Association. 106 (494): 440–451. doi:10.1198/jasa.2011.ap10323. S2CID   53560432.
13. Roy, M.; Bouma, M. J.; Ionides, E. L.; Dhiman, R. C.; Pascual, M. (2013). "The potential elimination of Plasmodium vivax malaria by relapse treatment: Insights from a transmission model and surveillance data from NW India". PLOS Neglected Tropical Diseases. 7 (1): e1979. doi: 10.1371/journal.pntd.0001979 . PMC   3542148 . PMID   23326611.
14. Zhou, J.; Han, L.; Liu, S. (2013). "Nonlinear mixed-effects state space models with applications to HIV dynamics". Statistics and Probability Letters. 83 (5): 1448–1456. doi:10.1016/j.spl.2013.01.032.
15. Lavine, J.; Rohani, P. (2012). "Resolving pertussis immunity and vaccine effectiveness using incidence time series". Expert Review of Vaccines. 11 (11): 1319–1329. doi:10.1586/ERV.12.109. PMC   3595187 . PMID   23249232.
16. Blackwood, J. C.; Cummings, D. A. T.; Broutin, H.; Iamsirithaworn, S.; Rohani, P. (2013). "Deciphering the impacts of vaccination and immunity on pertussis epidemiology in Thailand". Proceedings of the National Academy of Sciences of the USA. 110 (23): 9595–9600. Bibcode:2013PNAS..110.9595B. doi: 10.1073/pnas.1220908110 . PMC   3677483 . PMID   23690587.
17. Blake, I. M.; Martin, R.; Goel, A.; Khetsuriani, N.; Everts, J.; Wolff, C.; Wassilak, S.; Aylward, R. B.; Grassly, N. C. (2014). "The role of older children and adults in wild poliovirus transmission". Proceedings of the National Academy of Sciences of the USA. 111 (29): 10604–10609. Bibcode:2014PNAS..11110604B. doi: 10.1073/pnas.1323688111 . PMC   4115498 . PMID   25002465.
18. He, D.; Ionides, E. L.; King, A. A. (2010). "Plug-and-play inference for disease dynamics: measles in large and small towns as a case study". Journal of the Royal Society Interface. 7 (43): 271–283. doi:10.1098/rsif.2009.0151. PMC   2842609 . PMID   19535416.
19. Ionides, E. L.. (2011). "Discussion on "Feature Matching in Time Series Modeling" by Y. Xia and H. Tong". Statistical Science. 26: 49–52. arXiv: 1201.1376 . doi:10.1214/11-STS345C. S2CID   88511724.
20. Blackwood, J. C.; Streicker, D. G.; Altizer, S.; Rohani, P. (2013). "Resolving the roles of immunity, pathogenesis, and immigration for rabies persistence in vampire bat". Proceedings of the National Academy of Sciences of the USA. 110 (51): 20837––20842. Bibcode:2013PNAS..11020837B. doi: 10.1073/pnas.1308817110 . PMC   3870737 . PMID   24297874.
21. Bhadra, A. (2010). "Discussion of "Particle Markov chain Monte Carlo methods" by C. Andrieu, A. Doucet and R. Holenstein". Journal of the Royal Statistical Society, Series B. 72 (3): 314–315. doi: 10.1111/j.1467-9868.2009.00736.x .
22. Breto, C. (2014). "On idiosyncratic stochasticity of financial leverage effects". Statistics and Probability Letters. 91: 20–26. arXiv: 1312.5496 . doi:10.1016/j.spl.2014.04.003. S2CID   122694545.
23. Lindstrom, E.; Ionides, E. L.; Frydendall, J.; Madsen, H. (2012). "Efficient Iterated Filtering". System Identification. 45 (16): 1785–1790. doi: 10.3182/20120711-3-BE-2027.00300 .
24. Lindstrom, E. (2013). "Tuned iterated filtering". Statistics and Probability Letters. 83 (9): 2077–2080. doi:10.1016/j.spl.2013.05.019.
25. Doucet, A.; Jacob, P. E.; Rubenthaler, S. (2013). "Derivative-Free Estimation of the Score Vector and Observed Information Matrix with Application to State-Space Models". arXiv: 1304.5768 [stat.ME].