Sweden's infection rate may have peaked at less than 20% - good for Sweden but a deeply troubling number for the workhorse "SIR" epidemic model. It implies that the ratio of infection rate to recovery rate is 5:4 and much lower than estimates of 2:1 inferred from the spread of the disease.Drawing more from the financial literature than epidemiology, I provide a simple model for early stage growth of an epidemic that is more likely to be consistent with the time series and also more consistent with observed levels of peak infection.
A key benefit of leaning on the financial theory is that convexity adjustments are built into the pricing of bonds, and these may prove critical to understanding why many forecasts of COVID-19 were overly dire.
1. Negative Interest Rates
Imagine a world where the instantaneous rate of interest on your savings account is negative (some countries don't have to imagine slightly negative interest rates, although we are going to stretch the thought experiment). Every day the balance in your bank account goes down by a considerable amount. The negative interest rates are large and also very volatile. On some days you lose five percent of your principle. On other days you lose ten percent. Distressed, you consider hiding money under the mattress but legislation is passed which threatens to confiscate cash.
With some regret you decide to purchase a zero coupon bond figuring that you will be able to guarantee receiving the principle at maturity - and that might be better than nothing. The bond will return you $100 at the date of your choosing. To your horror you discover that bonds are hideously expensive. It will cost you $2000 to buy a bond from the U.S. Treasury with a one month maturity. Since this locks in a loss of $1900, you figure you might as well contemplate riskier options as well. You ask for a plot of bond prices as a function of maturity in order to survey your options. You are not pleased by the fact that the bank needs to use logarithmic paper to represent the price you might pay.
The Vasicek Model
If you take a careful look at the bond prices I have plotted you will notice that the expected interest rate today is not the same as the anticipated interest rate of a week from now, or two weeks from now. That's not terribly profound in interest rate modeling. People buying and selling bonds in our fictitious world have been thinking about the longevity of the negative interest rate crisis we imagine, and they don't think it will last forever.
Notice that the rate at which the price of a bond rises on the logarithmic plot starts out with a steep slope but further along the curve it is flatter. The tangent lines reflect this. Traders price bonds as if they expect the interest rate to start out at a high level and mean revert to a lower interest rate - though there is a heathy literature on the subtleties of that last assertion.
One way to account for term structure is through the use of the Vasicek interest rate model. In this model the rate of interest follows a random walk that is drawn towards a mean reversion level.
The walk starts off at an interest rate which could be set equal to the slope of the first tangent you see in the plot. The eventual average rate is a parameter also, which gives you a chance to decide what the mean growth will be later on. There is a speed parameter kappa that in our case might be on the order of a few percent, if we measure time in days. And there is a variance parameter that determines the extent of the randomness of the walk.
Now nobody is claiming that the Vasicek model is the best interest rate model, but it will give you a better chance to fit data than this one
which is analogous to the constant parameter SIR model for disease, at least as far as early stage growth is concerned. Pro tip: flat term structures are not the recipe for the next great fixed income hedge fund.
Variation in Infection Rates
I did not take those bond price plots you see from a Bloomberg screen. As I think you may have anticipated, they are COVID-19 growth plots (from Wikipedia).
Let's connect the dots between the way bond prices grow with maturity and the way epidemics grow with time. First, in the Vasicek stochastic differential equation you can add and subtract the same constant from the drift.
In the SIR epidemic model, early stage growth is driven by the difference between infection rate (beta) and recovery rate (gamma). You can see from the above that if gamma is constant and infection rate follows a "Vasicek" process (we should say Ornstein-Uhlenbeck) then so does the difference between beta and gamma. We will use r(t) to denote this differential. Both growth r(t) and infection rate are Ornstein-Uhlenbeck processes.
There are good reasons for a choice of mean reverting, initially decaying, infection rate. The initial fast drop in infection rate can model the impact of repeat contacts - you can only infect someone else once. This was noted in my prior article and the accompanying working paper where a continuous agents model is solved by a compartmental model. When you run the delayed differential equations for that system, you get infection rates shown in the wavy lines below. Our Vasicek model can approximate this.
Empirical surveys have also suggested slowing growth, prompting Chowell et al 's to suggest the exact same form:
which is the deterministic version of a Vasicek infection rate. A very recent paper by Ives and Bozzuto includes these R0 estimates for New York (paper)
Mean reverting infection rate may help model what you see. However if we want to understand the relationship of this pattern to herd immunity we also need variation.
Multiple Copies at Once
We might want to model infection rates stochastically because undoubtedly they are. But there is another more useful interpretation in which all of the paths are run at once. Each neighborhood gets a different Brownian motion. Under these assumptions, the aggregate statistics across neighborhoods are achieved by the same averaging that would be performed in a Monte Carlo experiment. Both kinds of averaging are identical mathematically.
As you would expect from someone who has written numerous articles about COVID-19 modeling with the common thread of heterogeneity, your author views the ability to generate convexity adjustments (corrections for the "flaw of averages" if you will) as critical to matching both early stage growth and herd immunity. The modeling setup I present gives you several opportunities to generate convexity and you can spread out the convexity adjustments as you see fit.
If you need convincing that variation is real, here's another plot from Ives and Buzzoto:
Clearly variation should be modeled. Brave modelers might even consider variation amongst cohorts who are not spatially removed and clearly interact - though there might be better frameworks for that.
A New Use for an Old Bond Price Model
To make further progress we need to understand the predictions of this model. The most obvious question we might ask is what the R value is. Because we left the recovery rate constant there is a simple relationship between R and growth. But here growth refers to aggregate growth across all the neighborhoods so this must be computed.
Unfortunately one cannot simply use the average instantaneous growth across neighborhoods. Instead one has to add up all the people who get infected and then back into an overall growth rate. This means taking an average of the exponential of the integral of a mean reverting Ornstein-Uhlenbeck process, where the average is across all possible paths taken by our random walk.
The Vasicek bond price formula for a bond is a calculation of an average value of ... the exponential integral of an Ornstein-Uhlenbeck process. What a convenient coincidence, except that the sign is reversed but that's easily fixed. We will denote the interest rate by x to avoid confusion with r (which will take the opposite sign) and likewise x0 is the initial value it takes. The bond price is
and D(t) is another function of time that does not involve x0.
This is an example of an affine bond price because the logarithm of the bond price is linear in the initial value x0. To recycle Vasicek's bond price formula to calculate mean growth rates for our epidemic model we reflect r(t) below the x-axis and treated as an interest rate x(t) - one that happens to be negative. The point is that
which is a trivial but important identification. Notice that this only works if everything is reflected however. We need to plug into the Vasicek bond price formula a negative initial interest rate and a negative mean reversion level.
The Vasicek interest rate model was always one of my personal favorites even though it had a seeming flaw: it allows the interest rate to go negative. We see that this bug is now a great feature. It gives us the mean growth across the entire population in our epidemic model and it looks like this:
where r0 we recall is the difference between the initial value of infection rate and recovery rate, and the other terms do not involve r0. The function tau grows at the same rate at time initially. Here B just the Vasicek term
except with x=-r. The negatives cancel out so rewriting as a slightly different B(r) is cleaner.
To keep the analogy to finance rolling along let's introduce the notion of yields. First some useful parameters
The first parameter is a ratio of the amount by which initial growth rate exceeds the eventual growth rate, which related to the extent to which initial infection rate exceeds eventual infection rate - though we have to be careful about the denominator. The second parameter is the ratio of the eventual standard deviation of growth rate to eventual mean growth rate (similar comment applies to the asymptotic standard deviation of infection rate).
Using this clean parametrization we can think of the growth equation in a different way. The ratio of mean growth to eventual instantaneous growth can be written
The basis functions look like this:
There is a flat term in blue which is all you get from the stock standard epidemiological model. There is an attenuation term suggested by Chowell (paper) and my own work. And there is a curvy line that captures the effect of stochastic infection or, in our interpretation, variation in infection rates across different neighborhoods.
To emphasize this is a non-orthogonal basis for epidemic yields. But you can put these yields together in mostly any combination. There is an exception because notice that you subtract the green curve (and scale it by the square of the standard deviation ratio). I've scaled down the green line so it is easy to see, but it will usually attract a small coefficient.
Resolving the Swedish puzzle
Let's return to the Swedish puzzle. A peak in infection seems to have occurred, consistent with a low value of R. But a low value of R cannot account for how fast it spread in the first place.
Sweden's controversial strategy does not equate to doing nothing at all (and we can be sure that behavioral changes are a big part of the disappearing growth) but the relatively lax approach probably tightens the relationship between early stage R values and implied R values at infection peak, at least compared to other countries. Sweden gives us more of a genuine puzzle than other countries, for now.
Taking a more conservative R estimate from a somewhat later stage doesn't make our task much easier. This paper by Tom Britton from the Department of Mathematics at Stockholm University was authored April 15th. It estimated early stage doubling times in the Stockholm region of 3.5 days prior to interventions, but the author's contemporaneous estimate mid April, a full month after measures were put in place, corresponded to a considerably longer doubling time of around 9 days. So let us take that later estimate and see if there is still something to be resolved.
On April 15th the author suggested a forward looking R value of 1.58 be used (considerably lower than his estimate R=2.5 in the earlier stages). His parameters for the Stockholm region were:
Again gamma is the recovery rate, beta is the infection rate. The reproduction number R is the ratio of beta to gamma. The typical generation time is the inverse of gamma, equal to 7.5 days. Peak infection with these parameters means 1/R of the population will not be infected, which is to say that 37% of the population will be infected. That is still looking too high by a factor of two. I refer the reader to Britton's paper for an explanation of how growth rates translated into R values, but his approach was very standard.
The story of Sweden is yet to play out fully and perhaps there are more twists and turns to come. But the important point is that there may be no puzzle at all if the Vasicek model or something like it is employed. Let's say that even after accounting for first moment effects as with Chowell's suggestion we still have a missing ratio. For arguments sake we need an extra 70% early stage growth to be plausible. Then the missing factor might be right there in front of us made up of convexity adjustments that work in tandem:
Both the terms on the right hand side are 1.0 in the SIR constant parameter model, so no wonder there is a puzzle. In words:
This might resolve the paradox.
If not there is also a third effect which you can also include if you wish, though it is harder to get a grip on analytically (at least for me). Due to variation in infection rates, individual peaks occur at different times. So what looks like a low peak might, in the absence of variation, have been a higher one. This is a superposition effect, if you will. One halfway decent approach to quantification of it might use the first passage time of a lognormal process, or the analytic solution to the SIR model.
Still got a puzzle? Take a closer look at the generation time. This is usually assumed to be exponentially distributed but if not, we need another correction.
For Next Time: Tractable Meta-Population Models
It is also possible to take mixtures of models where each model is a copy of the model with independent neighborhoods as above. We introduce a mixing variable omega that can vary the parameters.
There are two interpretations of this mixing. We could imagine that omega is drawn from some abstract probability space. It is a die cast before the start of each simulation. But it is probably more useful to imagine that the simulations take place simultaneously, and omega is interpreted as indexing a larger geography - and capturing additional spatial variation. Then we have a collection of populations, each of which is a collection of populations. So for instance omega might be a city whereas a path of Wt is a neighborhood.
There is a subtle interpretation of this last mixture model that makes it extremely attractive. This is the subject of this article, which is Part II of this post.
On the topic of Vasicek and epidemiological convexity adjustments:
My current best bet for improving COVID-19 prediction is taking shape at www.microprediction.org and you are welcome to contribute algorithms there.
AcknowledgementsThanks Graham Giller for pointing me to Ives and Buzzuto. Graham has a number of interesting posts on COVID-19 modeling at @stattrader on Medium.