Spatial considerations for the allocation of pre-pandemic influenza vaccination in the United States
Joseph T Wu*†, Steven Riley*, Gabriel M Leung*
* Department of Community Medicine and School of Public Health, Li Ka Shing Faculty of Medicine, The University of Hong Kong, Hong Kong SAR, China
† To whom correspondence should be addressed:
Supplementary Information: Text and Figures
Reduction in IAR is a convex function of coverage
We first use the simple homogeneous mixing model to illustrate why the relation between reduction in IAR and coverage is convex as shown in Figure 1B. For simplicity, we assume that the vaccine provides 100% reduction in susceptibility. The attack rate for a given coverage can be obtained by solving the following equation (Diekmann and Heesterbeeck, 2002) :
The second derivative of can be easily obtained as follows:
Since and is the only root of in , is a concave function of . Since the number of infections averted is , where is the population size, it is a convex function of coverage . Computational results suggest that this convex relation also holds in the vaccine response model used in this study (see below): We randomly generate 1,000,000 sets of parameters from all possible parameter values with up to 10 susceptible states () and convexity holds in all these scenarios.
Mathematical details for the example in the Results section: Why is pro-rata the least efficient policy?
Here, we provide detailed mathematical argument as to why pro-rata is the least efficient policy for the example in the introduction. Consider populations of sizes that are isolated from each other. Let be the reduction in IAR achieved by vaccine coverage level . All populations have the same , which is convex and becomes flat only after the critical coverage level has been achieved (Figure 1B). Given doses of vaccines and an allocation vector , where is the proportion of the overall stockpile allocated to population , the total number of infections averted is . Note that we must have . Substituting this relation into, we have . Since , is a stationary point of . Further, since ( is convex when ), we have
Therefore, is a convex function of . Taken together, is minimized at , which is the pro-rata allocation. That is, pro-rata allocation gives the least number of infections averted and is thus the least efficient policy.
Mathematical details for the example in Figure 4: Why is it never optimal to split the stockpile unless at least one population has achieved control?
Consider the example in Figure 4 where Population 1 and 2 are non-interacting. For , let be the number of infections averted in Population i if doses of vaccines are allocated to it, and be the number of doses needed to reach critical coverage for Population i. Note that and (see Figure 1B) where is the same reduction in attack rate for coverage defined above. Given doses of vaccines, let be the total number of infections averted if a proportion of the overall stockpile is allocated to Population 1. The goal is to minimize the total number of infections. When is not large enough to achieve control simultaneously in both populations, which is necessary for the allocation problem to be worthwhile, we only need to consider allocations where and . In this range, and are strictly convex (i.e. ), so . This means that all stationary points of are minima and is maximized at either or. If , either all vaccines are allocated to Population 2 () or Population 2 has achieved control (). Similarly, if , either all vaccines are allocated to Population 1 () or Population 1 has achieved control (). Taken together, these observations imply that it is never optimal to split the stockpile unless at least one population has achieved control. These observations also explain concentration of resources and choppiness (discontinuity in allocation arises when the optimal allocation switches between and).
The Meta-population Model Without Vaccination
Let , , and be the number of susceptible, infectious, and removed individuals in populationwhere . Given a mixing matrix, where is the average proportion of time that a resident of populationspends in population(see next section), the model without vaccination is defined by
for , whereis the mean infectious duration. Since we are only concerned about the final attack rates, omitting the exposed class, which corresponds to the latent period, has no effects on our results.
Inter-population Mixing: Construction of the Mixing Matrix
Thematrix. We constructed a mixing matrix whereis the average proportion of time that a resident of populationspends in populationover one year. To this end, we first construct the following matrices:
· whereand, , is the number of air person-trips from populationto populationover one year.
· whereand,, is the number of individuals who reside in populationand commute to work in populationon a daily basis.
· whereand,, is the total duration (in days) of non-air person-trips from populationto populationover one year with travel distance at least 100 miles.
· whereand ,, is the average duration of an air person-trip from populationto population.
Construction of these matrices is documented below. We estimated the matrix by
Thematrix: Air Travel Volume. We estimate the US domestic air travel volume using the 2005 T-100 Domestic Market and the 2005 Airline Origin and Destination Survey (DB1B) from the Bureau of Transportation Statistics (www.transtats.bts.gov/Tables.asp?DB_ID=125&DB_Name=Airline%20Origin%20and%20Destination%20Survey%20(DB1B)). DB1B is 10% sample of airline tickets from reporting carriers collected by the Office of Airline Information of the Bureau of Transportation Statistics. Data include origin, destination, and other itinerary details of passengers transported for each quarter. On average, about 70% of these itineraries (weighted by the number of passengers in each itinerary) are roundtrips. A unique destination can be identified for more than 99.5% of these roundtrips (the others have multiple stops before returning to the origin). For these trips, letbe the number of passengers traveled from population to population and be the total number of coupons during quarter . Let be the total number of passengers during quarter (as reported in T-100). We estimated the matrixby
.
Thematrix: Workflow Volume. We estimate workflows using the United States Census 2000 County-To-County Worker Flow Files (www.census.gov/population/www/cen2000/commuting.html). We choose from the database only residence-workplace pairs for which the distance between the two counties is less than 200 km; these account for 99% of all residence-workplace pairs and 83% of those for which residence and workplace are in different states. We compute the distances using the latitudes and longitudes provided in the Census 2000 US Gazetteer Files (www.census.gov/geo/www/gazetteer/places2k.html).
Thematrix: Long-Distance Non-Air Travel Volume. We estimate the long-distance non-air travel volume using the 1995 American Travel Survey (www.bts.gov/publications/1995_american_travel_survey/index.html). We set if the number of reportednon-air person-trips is positive but less than 15; these origin-destination pairs account for only 1.2% of non-air person-trips in the database.
Thematrix: Air-Travel Trip Average Duration. We estimate the interstate air-travel trip average durations using the 1995 American Travel Survey (www.bts.gov/publications/1995_american_travel_survey/index.html). Letbe the set of pairs for which the number of reported air person-trips is positive but less than 15; these account for 10.8% of the air person-trips in the database. Letandbe the average duration of air person-trips with stateas the origin and destination, respectively, in the database withexcluded. Letbe the average duration of all air person-trips other than those in. We estimate,, by
ORMAT Vaccine Efficacy Model
After an individual is vaccinated, his susceptibility and infectiousness (if he becomes infectious) are reduced by a factor of and, respectively. We assume that andare independent random variables. We denote the different states of andin ascending order by,and,, respectively. We say that an individual is of class, and, if his relative susceptibility and relative infectiousness areand, respectively. Letbe the probability of an individual becoming classupon vaccination.
The Meta-population Model With Individual Vaccine Response
Let be the number of classsusceptible in population, where, , . Define and similarly. The full meta-population model is
where, , .
Algorithm for Calculating Attack Rates
We develop an algorithm to calculate the attack rates. This algorithm allows us to obtain the attack rates by solving a system of nonlinear equations regardless of the number of immune states introduced in the vaccine response model. In the following, we first illustrate the algorithm for a single population and then generalize to the meta-population model.
Attack Rate in a Single Population. Consider a single population of size . Let and be the initial number of classindividuals and final number of infections in class, respectively. The attack rate can be obtained by solving for’s in
(1)
Given a vaccine coverage of, the initial number of class individuals are
The numbers of infections ’s can be computed using the following algorithm:
- Express in terms of as follows:
. (2)
This relation is obtained using the relation which is true for all and.
- Solve for in the Equation (1)
where the double summation can be expressed as a function of only using Equation (2) as follows:
where .
- Using the solution of obtained in Step 2, the total number of infections (again, using Equation (2)) is
.
From this algorithm, we see that ’s, and therefore the attack rate, depend on via only and .
Attack Rates in the Meta-Population. Let and be the initial number of classindividuals and the final number of infections in class in population, respectively. As in the case for a single population, the attack rate in a meta-population can be obtained by solving for’s in the following system of equations:
The numbers of infections’s can be computed using the algorithm for a single population but with the following modification: In Step 2, solve for’s in the following system of equations:
.
The Next Generation Matrix and Reproductive Number
The next generation matrix is used to compute the reproductive number. For a completely susceptible meta-population, the next generation matrix is amatrixwhereis the expected number of secondary infections in the fully susceptible populationresulting from a randomly selected infectious individual in population(Heesterbeek, 2002). For our meta-population model,
.
The basic reproductive number is equal to the largest eigenvalue of the next generation matrix. With vaccination, the next generation matrix becomes a matrix whereis the expected number of secondary infections among classindividuals in populationfrom a randomly selectedinfectious individual in population, where :
.
Minimizaton of Post-vaccination Reproductive Number
Pro-rata is Near-Optimal. We first use a special case of our meta-population model to illustrate why pro-rata almost minimizes the post-vaccination reproductive number. Suppose there are only two populations with sizes and . The two populations have the same local transmissibility with a mixing matrix
for an arbitrary between 0 and 1. Suppose the vaccine provides 100% reduction in susceptibility and let be the proportion of a vaccine stockpile of size allocated to Population 1. The next generation matrix after vaccination is
where is the coverage in Population . The reproductive number is the maximum eigenvalue of this next generation matrix. Therefore, the reproductive number is minimized by finding a that minimizes this eigenvalue. Straightforward calculations show that the reproductive number is minimized when , which is the pro-rata allocation. Supplementary Figure 5A shows that in the 10-region US model, the reproductive number under pro-rata is larger than the optimal by no more than 15%. Supplementary Figure 5B shows that allocations that minimize the attack rate do not minimize the reproductive number.
Meta-population with the Same Reproductive Number can have Different Attack Rates. To show this, consider a simple 2-population epidemic system where the two populations have the same size:
We randomly generate 10,000 transmission matrix as follows: where ’s are i.i.d. U(0,1) random variables and is a scaling factor such that the resulting matrix gives a reproductive number of 1.8. Supplementary Figure 6 shows that these simple systems have different attack rates even though they have the same reproductive number.
Supplementary Figure 1. as a function of basic reproductive number(y-axis) and coverage (x-axis) in the 4-region and 49-state model. The black contour marks the boundary below which 6%. A. The 4-region model. B. The 49-state model.
Supplementary Figure 2. Reduction in attack rate for the 100% discretionary policy in the 10-region model when there is no inter-regional mixing (i.e. regions are isolated from each other). The black contour marks the boundary below which 6%.
Supplementary Figure 3. Univariate (A-C) and multivariate (D) sensitivity analyses for the 50% discretionary policy in the 10-region model. (A-C) Coverage is fixed at 20%. The black contours in each panel mark the boundary at which and the black dashed vertical line indicates the base case parameter values. A. drops as inter-regional mixing () increases. B. drops as antigenic match () deviates from the base case value. can become negative when antigenic match is near perfect. C. is mostly unaffected when transmissibility is overestimated () but drops to significantly negative values when transmissibility is overestimated (). D. Multivariate sensitivity analysis. When base case assumptions (used to construct the discretionary policies) are not satisfied, the actual is typically lower than the expected in the base case (see caption of Figure 5B on how the multivariate analysis is performed).