Ignoring transmission dynamics leads to underestimation of the impact of a novel intervention against mosquito-borne disease

Wolbachia is an intracellular bacterium that many hope could have a major impact on dengue and other mosquito-borne diseases that are notoriously difficult to control. The balance of future investments in Wolbachia versus other public health needs will be informed to a great extent by efficacy estimates from large-scale trials, which can be affected by multiple sources of bias. We used mathematical models to quantify the possible magnitude of these biases, finding that efficacy would have been severely underestimated in a recent trial in Indonesia if the spatial scale of clusters had been smaller than it was. We also identified a previously unrecognized source of bias owing to the coupled nature of transmission dynamics across clusters. This too led to a consistent underestimate of the protection afforded by Wolbachia. Taken together, our findings suggest that this intervention may be even more promising than currently recognized.

across geographical contexts, as has been recently done for the RTS,S/AS01 vaccine 11,12 and the endectocide ivermectin 13 for malaria. If failing to account for such transmission dynamics contributes to an underestimated biological effect of Wolbachia on DENV, we risk incorrectly assessing its broader impact.
In this study, we used a simple mathematical model of DENV transmission to gain insight into the possible magnitudes of these three sources of bias. Our approach involved translating model inputs of the basic reproduction number (R 0 ), population susceptibility (S), the spatial scale of human movement (b), and the proportional reduction in R 0 afforded by Wolbachia-infected mosquitoes (ε) into outputs of the infection attack rate (IAR) in control and treatment arms of a trial, in accordance with a structured, susceptible-infectious-recovered (SIR) model 14 . We assumed a checkerboard pattern of control and treatment arms of 1 km 2 to approximate the design used in the AWED trial, which covered the entire city of Yogyakarta with neighboring areas assigned to one arm or another in an (approximately) alternating pattern (Fig.   1A) 7 . We used the outputs of IAR in treatment and control arms (IAR t and IAR c , respectively) to obtain an estimate of the odds ratio of infection, and thereby an estimate of the efficacy of the intervention, Eff = 1 -OR. This is comparable to the measure of efficacy used in the AWED trial, which also used 1 -OR, but in a test-negative design with virological-confirmed symptomatic dengue as the primary end point (see Online Methods for a justification of the comparability of these approaches). We constructed six different model versions for estimating efficacy, each of which includes different combinations of the three biases, all of them, or none of them.
Henceforth, we refer to the efficacy observed in the AWED trial as "observed efficacy," and the efficacy estimated by a given model and ε as "estimated efficacy." Finally, we quantify each bias as the difference in the efficacy estimated by a model including that bias and a model which does not include that bias (see Online Methods for more details of our methods).
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021.  Table S1 unless otherwise stated.
The relationship between the efficacy estimated by the model with all three forms of bias (the estimated efficacy) and the reduction in R 0 (ε) was dependent on the amount of time people spent in their allocated arm (Fig. 1C)-the less time individuals spent in their allocated arm, the higher the reduction in R 0 that was needed to recreate the observed efficacy from the AWED  1D).
When we fixed ε to the value that reproduces the observed efficacy in the AWED trial for b = 60 m, and increased human movement between arms by increasing b, the estimated efficacy by the model accounting for all three forms of bias decreased (Fig. 1E). For example, increasing the average distance in one direction between transmission pairs (b) from 60 m to 120 m caused a relative reduction of 25.0% in estimated efficacy, highlighting the sensitivity of efficacy to the spatial scale of human movement. This effect occurs for two reasons: firstly, as people spend less time in their allocated arm, the time that people spend under the intervention becomes more similar between arms; and secondly, as fewer people are infected overall, the local force of infection is reduced over time in both arms. Relatedly, estimated efficacy depended on the dimensions of the trial clusters, which we set to 1 km 2 by default (Fig. 1F).
When we reduced the cluster dimensions to 500 m x 500 m, estimated efficacy dropped from 77.1% to 57.8%, representing a 25.0% relative reduction. This effect occurs because as the cluster dimensions reduce, people spend less time in their home cluster and, hence, the time spent in each trial arm approaches parity (i.e., 50%). Increasing cluster dimensions above 1 km 2 . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint had less of an effect on estimated efficacy. For example, increasing the cluster dimensions to 2 km x 2 km resulted in an estimated efficacy of 87.9%, a relative increase of 14.0%.
Our approach enabled us to directly and separately model each of the three potential sources of bias: (1) mosquito movement, (2) human movement, and (3) transmission coupling.
When we assumed that ε was equal to 70.3%, allowing for mosquito movement but not human movement produced an estimated efficacy of 100%, because there was no transmission in the intervention arm in that case ( Fig. 2A, Fig. S1). If we allowed for both mosquito movement and human movement, we observed a lower estimated efficacy of 89.3%. Although there was no transmission in the intervention arm in this case, individuals residing in the intervention arm could be infected in the control arm. Additionally, those assigned to the control arm experienced lower overall risk due to their time spent in the intervention arm. When we accounted for transmission coupling between trial arms alongside human and mosquito movement, thereby allowing for transmission in the intervention arm, risk was the most similar across the trial arms of all scenarios, leading to the lowest estimated efficacy of 77.1% when ε was 70.3%.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021.

Fig. 2: Sources of bias in efficacy estimates.
In both panels, yellow refers to mosquito movement, red to human movement, and blue to transmission coupling. A: The relationship between the reduction in R 0 (ε) and the estimated efficacy for the six possible models. The black line here is the relationship for a model with no human movement or mosquito movement. Where a line has more than one color, it represents the model which includes each of the types of bias represented by those colors. The difference between this line and each of the colored lines represents the bias introduced by not accounting for the features present in the model described by that colored line. B: the contribution of each source of bias to the total bias.
We quantified total bias as Eff (hmt) -Eff (0) , where Eff (hmt) is the estimated efficacy under the model with all sources of bias and Eff (0) is the estimated efficacy under the model without human or mosquito movement. We then computed the difference in the bias produced by pairs of models to decompose overall bias into each of its three sources (Fig. 2B, Fig. S2-3, see Online methods for details). At the baseline ε of 70.3%, 2.8% of the total bias was attributable to . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint  2B) and absolute (Fig. S2) contributions of transmission coupling to bias were greatest. When ε was below a value of around 35%, the effective reproduction number at the start of the trial exceeded 1 in both arms. This value of ε varied slightly based on the model used ( Fig. S2-3), but was always approximately 1 -1/(S 0 R 0 ).
When ε was below this critical value, increasing it in the context of coupled transmission reduced incidence in the control arm and caused smaller reductions in incidence in the intervention arm than if transmission had been uncoupled (Fig. S1, e.g. panels D vs. F). This implies that the bias introduced by transmission coupling increases as ε increases up to ~35% under our model's parameterization (Fig. 2B). Increasing ε past this point does not lead to further reductions in incidence in the intervention arm in an uncoupled model, as incidence is already minimized. However, for the model with coupled transmission, incidence can continue to decrease in both arms as ε increases. Hence, the relative difference between arms for a coupled model compared to an uncoupled model increases at higher values of ε. In other words, the bias caused by transmission coupling decreases while the bias caused by human movement continues to increase (Fig. S2).
Our results highlight three sources of bias (human movement, mosquito movement, and transmission coupling) that arise in large, cluster-randomized, controlled trials of interventions against mosquito-borne diseases, and have implications for how to mitigate these biases.
Biases arising due to human movement and mosquito movement are typically able to be addressed through careful statistical analysis of trial data or in the design of the trial 8 . For instance, in the analysis of the AWED trial, Utarini et al. accounted for these two forms of bias by combining self-reported recent travel and local Wolbachia prevalence into an individual-level Wolbachia exposure index 7 . Comparing groups with the highest and lowest Wolbachia exposure did not lead to higher efficacy estimates than their primary analysis. Another recently proposed . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint approach to addressing contamination involves describing the effectiveness of the intervention at the boundary between clusters using a sigmoid function 15,16 . Our results suggest that failure to take steps such as this to account for human and mosquito movement would typically lead to underestimated efficacy, while failure to account for transmission coupling would lead to an even greater underestimate, particularly at intermediate reductions in R 0 .
Bias arising from human and mosquito movement could also be mitigated at the stage of planning the trial. The classical design to achieve this is the 'fried-egg' design, in which a treated buffer-zone is placed between intervention and control clusters 17 . A more recently proposed approach involves excluding a subset of clusters from the trial completely, thereby increasing the distance between clusters and leading to disconnected clusters at less risk of contamination 18 . While both of these approaches do mitigate the risk of contamination directly, they also necessitate a larger trial area, and may be logistically infeasible in a trial taking place in a single city, such as the AWED trial. Another approach could include reducing the number of clusters, but keeping total area fixed, leading individuals to spend more time in their assigned arm and reducing mosquito movement by reducing the boundary between clusters. Our results show that the efficacy estimated from cluster-randomized, controlled trials of interventions against mosquito-borne diseases is highly sensitive to cluster size (Fig. 1F). Had the dimensions of the clusters in the AWED trial been much smaller, then the estimated efficacy may have been substantially lower. However, having fewer large clusters would likely introduce new biases by making the arms less comparable, which may not be an acceptable trade-off.
While bias due to human and mosquito movement can be mitigated through trial design and statistical methods, our results highlight a third source of bias, transmission coupling, that requires other tools to fully address. Accounting for this bias first requires data on the spatial distribution of the intervention and on human movement, similar to that used in the supplementary analysis of the AWED trial. However, it also requires interfacing these data with a dynamical transmission model to account for the fact that, in the presence of movement . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; between arms, incidence in each arm depends on prevalence in both arms 19 . It may also be possible to design a trial in ways that reduce the bias due to transmission coupling, for instance by allocating a greater proportion of the trial area to the intervention arm, with small intervention clusters situated among larger control clusters so that transmission suppression in the intervention arm has less of a population-level effect. The ratio of area allotted to treatment and control would depend on many factors, including the expected strength of the intervention, the local force of infection, and logistical constraints such as the size and length of the trial. Utilizing a dynamical model synthesizing these factors in the design of a trial could aid in understanding how different designs might affect bias due to transmission coupling 19 . More work is needed to understand what types of spatial clustering patterns, among other features of trial design, would minimize this form of bias.
Although our modeling approach allowed us to account for different potential sources of bias and to attribute the total bias to each of those sources, it has at least three limitations. First, our model was deterministic, yet stochasticity could be important for a highly efficacious intervention with potential to reduce transmission to very low levels 20 . This simplification implies that our estimates are likely conservative, as these effects could increase the bias due to transmission coupling if a highly effective intervention increases the probability of transmission fadeouts. Second, our simple model does not reflect all of the complexities of DENV transmission. For example, our model does not account for seasonality, so dengue outbreaks in our model end by susceptible depletion rather than by seasonal reductions 21 . We also did not account for spatial heterogeneties in transmission or interactions between serotypes. Accurately quantifying the contribution of these effects to bias would require a more detailed model, but the qualitative results would likely be similar. Finally, we did not calibrate our model to trial data.
However, our aim here was not to precisely quantify bias in the AWED trial, but rather to highlight some potential sources of bias in trials of that nature and to understand how these biases are influenced by transmission dynamics and human mobility. Moreover, although our . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; estimate of R 0 may not be reflective of the transmission context during the AWED trial, our qualitative results were robust to different values of R 0 .
In conclusion, without accounting for human movement, mosquito movement, and transmission coupling, the efficacy of Wolbachia-infect mosquitoes as an intervention to control dengue is likely to be underestimated. As the estimate of efficacy in the AWED trial was already very high (77.1% [95% CI: 65.3% -84.9%]) 7 and, as we show, likely underestimated, Wolbachia-infected mosquitoes have potential to be a game-changing tool in the fight against dengue. Even as vaccines against dengue become available, a variety of vector control approaches are likely to remain key tools in the fight against dengue 2,22 . Although we focused our study on a trial of Wolbachia-infected mosquitoes, our findings are applicable to any efficacy trial of an intervention that has the potential to contaminate the control arm, such as gene drive mosquitoes or ivermectin as interventions against malaria 23,24 . As trials of these interventions continue, it will be important to learn what lessons we can from transmission dynamic modeling when designing and interpreting future trials, to ensure that we understand the true promise of these interventions.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint

Attack rate models
Let ε represent the effectiveness of the intervention, defined as the reduction in the pre-intervention basic reproduction number, R 0 , when the intervention is applied at full coverage in a treatment cluster. Hence, (S2) 0, = 0 Given the epidemic behavior of an SIR model, our interest is in quantifying the infection attack rate (IAR), , within each cluster during a trial. We do this for three different types of bias: human movement between arms, mosquito movement between arms, and transmission coupling between arms. This results in six different models, as transmission coupling can only occur in the presence of human movement, which are described below.

No bias
In the absence of contamination from human movement or mosquito movement between arms, an implicit solution for t (0) and c (0) can be obtained from ) as described by Miller 14 . In eqs. (S3-S4), S′ is the initial susceptibility of the population, which is assumed to be equal in treatment and control clusters.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint

Bias from human movement
Let ϱ ij represent the proportion of the total time at risk that a resident of cluster i spends in cluster j. We calculate the effect of movement on the infection attack rates as weighted averages of the IARs computed under the no-bias scenario, because tt + tc = 1 and cc + ct = 1.

Bias from mosquito movement
We represent the coverage of the intervention-i.e., the proportion of Wolbachia-infected mosquitoes-in the two arms with C t and C c . In the case of mosquito movement, there may be non-zero coverage of intervention in the control arm (i.e., C c >0), and less than 100% coverage in the treatment arm (i.e., C t < 1). To calculate the attack rate under mosquito movement, we adjust the effectiveness in each arm, yielding This leads to the following implicit solutions for t (m) and c  . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. Here we are assuming that movement of mosquitoes between trial arms does not directly contribute to transmission via movement of infected mosquitoes. This discrepancy can be reconciled by the fact that the spread of dengue virus occurs within a single mosquito generation, whereas the spread of Wolbachia occurs over the course of multiple generations.

Bias from human movement and transmission coupling
Thus far, we have assumed that transmission in each arm is only a function of prevalence in that arm, and not in the other. To relax this assumption, we couple transmission between the two arms according to a dynamic two-patch transmission model. We then follow Miller 14  ) . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

Bias from human movement, mosquito movement, and transmission coupling
Finally, we include all three forms of bias by altering equations S7-S8 to account for differential coverage between arms:

Efficacy calculation
The ratio of the IARs in the treatment and control clusters is an infection risk ratio. However, the AWED trial based their efficacy calculations upon an odds ratio 7 , with symptomatic, virologically-confirmed dengue as the end point. That is, efficacy in the trial was computed as 1-p t n c /p c n t , where p i and n i represent enrolled test-positives and test-negatives, respectively, in trial arm i. To generate a comparable quantity, we computed the efficacy according to model x as , (S9) i is the infection attack rate in trial arm i∈{c,t} for model x∈{0,h,m,hm,ht,hmt}. Here, we are assuming that the ratio of infections to enrolled test-positives does not differ between arms (i.e., p i =k p i for i∈{c,t}) and similarly the ratio of those uninfected to enrolled test-negatives also does not differ between arms (i.e., n i =k n (1i ) for i∈{c,t}). If either of these assumptions . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review) preprint
The copyright holder for this this version posted November 21, 2021. ; were violated, for instance if the intervention affected either the proportion of dengue infections that were symptomatic, then our estimate of efficacy would be less comparable to the estimate used in the AWED trial.

Bias calculation
We calculated the bias due to a particular source as the difference in the efficacy between a model with that source of bias and a model without that source of bias. As biases appear in multiple models, this led to three ways to embed the models, and three corresponding ways to quantify each bias. The three embeddings are: A) no bias→mosquito movement→human movement + mosquito movement→full model; B) no bias→human movement→human movement + mosquito movement→full model; and C) no bias→human movement→human movement + transmission coupling→full model. The difference between efficacy estimates for adjacent models in an embedding will lead to an expression for the bias which differs between the two models. Hence the three possible ways to calculate each of the three sources of bias yield the following equations: . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review) preprint
The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint (S12B) We then calculate the average total bias caused by each source of bias as , (S13) ( ) = ∈{ , , } ∑ ( ) 3 where i∈{h,m,t}. Note that bias (t) A =bias (t) B and bias (h) B =bias (h) C , but it is necessary to include each as a separate term so that each of the three model embeddings is included equally.

Apportionment of Time at Risk
We considered a checkerboard arrangement for the treatment and control clusters in a trial across a two-dimensional landscape (Fig. 1A, Fig. S4). Under this scenario, we assume that the population density per unit area is constant and that transmission potential, as captured by R 0 , is homogeneous across the landscape prior to initiation of the trial.
At the core of this derivation is the assumption that the location where an individual j Under the checkerboard arrangement, we considered alternating squares of width δ corresponding to treatment and control clusters within a contiguous urban area (Fig. 1A).
Although any such area would have borders in reality, we ignored any possible edges effects and assumed that the extent of interactions between squares of type t and c in the interior of the checkerboard provide a suitable characterization of overall interaction between individuals residing in t and c, as summarized by ρ tt and ρ cc . Because the area and arrangement of t and c squares are identical, ρ tt = ρ cc and ρ tc = ρ ct (Fig. 1A). (S19) We also need to calculate the time at risk in a non-adjacent interval of width 3 whose edge is spaced distance 2 away from the nearest edge of the interval where the individual resides, which has width 1 . Applying similar reasoning as in eq. (S17), we obtain . (S20)

Calculation of Initial Susceptibility
To obtain an estimate of initial susceptibility, we followed ten Bosch et al. 26 and calculated the proportion of the population exposed to n serotypes, ∀n∈{0,1,2,3,4}, as a function of age. We used an estimate of the force of infection of 0.0318/yr for dengue for Yogyakarta,Indonesia 27 and assumed that the force of infection was constant across space and serotype. Following ten Bosch et al. 26 , we defined e i (a) as the proportion of individuals of age a that have been exposed to i serotypes and r i (a) as the proportion of individuals of age a experiencing temporary heterologous immunity following exposure to i serotypes. The dynamics of how individuals progress through these classes as they age follows In eqs. (S22-S23), Λ is the force of infection, which we set to 0.0318/yr, and is the rate at σ which individuals lose heterologous immunity, which we set to 0.5/yr 26 .
We computed the proportion of the population in Yogyakarta, Indonesia that is of age a using estimates from the United Nations World Population Prospects database 28

Proportion of Wolbachia-Infected Mosquitoes
To determine the proportion of Wolbachia-infected mosquitoes in the treatment and control clusters, we extracted data on the establishment of Wolbachia in each cluster of the AWED trial 7 . We computed and as the average proportion of Wolbachia-infected mosquitoes across all clusters and time points for the treatment and control clusters, respectively. Under this calculation, was equal to 0.950 and was equal to 0.170.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint

Spatial Scale of Human Movement
Our calculations of the apportionment of time at risk depend upon a value of , the spatial scale of human movement, a quantity that is challenging to parameterize. To inform our selection of , we calculated the relationship between (i.e., the reduction in R) and estimated efficacy ε calculated using the model that combined the effects of mosquito movement, human movement, and transmission coupling. We compared the observed efficacies obtained for different values of b to the observed efficacy from the AWED trial. We chose b = 60 m, because the mean and 95% confidence interval limits of observed efficacy from the AWED trial were the only ones that implied realistic values of . This value of corresponded to and ε ≤ 1 = 60 ρ = ρ = 0. 887 for the checkerboard arrangement. ρ = ρ = 0. 113

Sensitivity Analysis
We performed sensitivity analyses to determine how the choice of parameter values affected the conclusions that we reached concerning the contribution of the three effects (i.e., mosquito movement, human movement, and transmission coupling) to the bias in efficacy. To do so, we either varied each parameter one at a time while holding all parameters at their default value, or varied all parameters simultaneously. In both cases, we calculated the total bias attributable to each source of bias as described previously. We report the parameter ranges and default values in Table S1, and values within these ranges were sampled using a Sobol design 29 . We then quantified the global sensitivity of the measured bias to each parameter by computing Sobol indices using the sensobol package in R 30 . We used the Saltelli estimator for first-order indices and the Glen estimator for the total-order indices, and a sample size of 10 5 parameter sets.

Supplementary Text -Sensitivity analysis
We undertook a global sensitivity analysis to understand how the level of bias is affected by our assumed model parameters (Figs S4-S6). For the model with only mosquito movement, the prevalence of Wolbachia-infected mosquitoes in the treatment arm and the value of ε had the biggest effect on the extent of bias, explaining 14% and 38% of the variance in bias, respectively (Fig. S5). For the model with only human movement, 70% of the variance in bias was explained by the proportion of time that individuals spent in their allocated arm (Fig. S6).
For the full model, no single parameter could explain more than 2% of the variance in the bias on its own. Instead, interactions between these parameters explained the variance ( . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint

Fig. S1: Infection attack rates for each of the six models, delineated by control and intervention arms.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint  . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint

Fig. S5: One-at-a-time sensitivity of bias to different parameter values for the model with only mosquito movement (red lines), and scatter plots showing output of global sensitivity analysis (black points).
Baseline parameter values are shown with the red dashed lines. Sobol indices are shown above the plots; the first order index indicates the proportion of the variance in the bias is attributable to variance in that parameter alone, and the total order index indicates the proportion of the bias is attributable to variance in that parameter including all interaction terms in which that parameter appears.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint

Figure S6: One-at-a-time sensitivity of bias to different parameter values for the model with only human movement (red lines), and scatter plots showing output of global sensitivity analysis (black points).
Baseline parameter values are shown with the red dashed lines. Sobol indices are shown above the plots; the first order index indicates the proportion of the variance in the bias is attributable to variance in that parameter alone, and the total order index indicates the proportion of the bias is attributable to variance in that parameter including all interaction terms in which that parameter appears.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint

Figure S7: One-at-a-time sensitivity of bias to different parameter values for the full model (red lines), and scatter plots showing output of global sensitivity analysis (black points).
Baseline parameter values are shown with the red dashed lines. Sobol indices are shown above the plots; the first order index indicates the proportion of the variance in the bias is attributable to variance in that parameter alone, and the total order index indicates the proportion of the bias is attributable to variance in that parameter including all interaction terms in which that parameter appears.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) preprint The copyright holder for this this version posted November 21, 2021. ; https://doi.org/10.1101/2021.11.19.21266602 doi: medRxiv preprint