Abstract
Objective We aimed to use mathematical models of SARS-COV-2 to assess the potential efficacy of non-pharmaceutical interventions on transmission in the parcel delivery and logistics sector.
Methods We developed a network-based model of workplace contacts based on data and consultations from companies in the parcel delivery and logistics sectors. We used these in stochastic simulations of disease transmission to predict the probability of workplace outbreaks in this settings. Individuals in the model have different viral load trajectories based on SARS-CoV-2 in-host dynamics, which couple to their infectiousness and test positive probability over time, in order to determine the impact of testing and isolation measures.
Results The baseline model (without any interventions) showed different workplace infection rates for staff in different job roles. Based on our assumptions of contact patterns in the parcel delivery work setting we found that when a delivery driver was the index case, on average they infect only 0.18 other employees, while for warehouse and office workers this went up to 0.93 and 2.58 respectively. In the large-item delivery setting this was predicted to be 0.83, 0.94, and 1.61 respectively. Nonetheless, the vast majority of simulations resulted in 0 secondary cases among customers (even without contact-free delivery). Our results showed that a combination of social distancing, office staff working from home, and fixed driver pairings (all interventions carried out by the companies we consulted) reduce the risk of workplace outbreaks by 3-4 times.
Conclusion This work suggests that, without interventions, significant transmission could occur in these workplaces, but that these pose minimal risk to customers. We found that identifying and isolating regular close-contacts of infectious individuals (i.e. house-share, carpools, or delivery pairs) is an efficient measure for stopping workplace outbreaks. Regular testing can make these isolation measures even more effective but also increases the number of staff isolating at one time. It is therefore more efficient to combine these measures with social distancing and contact reduction interventions, as these reduce both transmission and the number of people needing to isolate.
IMPORTANCE During the COVID-19 pandemic the home-delivery sector has been vital to maintaining people’s access to certain goods, and sustaining levels of economic activity for a variety of businesses. However, this important work necessarily involves contact with a large number of customers as well as colleagues. This means that questions have often been raised about whether enough is being done to keep customers and staff safe. Estimating the potential risk to customers and staff is complex, but here we tackle this problem by building a model of workplace and customer contacts, from which we simulate SARS-CoV-2 transmission. By involving industry representatives in the development of this model, we have simulated interventions that have either been applied or considered, and so the findings of this study are extremely relevant to decisions made in that sector moving forward. Furthermore, we can learn generic lessons from this specific case study which apply to many types of shared workplace as well as highlighting implications of the highly stochastic nature of disease transmission in small populations.
I. INTRODUCTION
Demand for home-delivery services has spiked globally during the COVID-19 pandemic, as people have stayed at home to reduce transmission [1]. In the UK, non-essential retail shops were closed for much of 2020 and 2021, increasing the demand for online retail and home delivery. Additionally, stay-at-home orders brought new demand for large items such as furniture and white goods as many people adjusted to spending more time at home [2]. This new and displaced demand has, on the whole, been successfully absorbed and managed by the delivery and logistics sector, due in no small part to the efforts of the key workers in those sectors to keep business moving, while adapting to a changing work environment. Meanwhile, key workers in all sectors were disproportionately exposed to transmission of SARS-CoV-2 [3]. In the delivery sector, drivers and warehouse workers were also at risk, given their exposure to a large number of contacts, the likelihood of asymptomatic transmission in SARS-CoV-2, and the potential economic impact of absence due to the prevalence of flexible or zero-hours contracts in this sector. Furthermore, studies from other countries indicate that delivery drivers there could be at much greater risk [4, 5] than the general population, and so is a sector that requires greater attention.
Mathematical models have been central to understanding transmission of SARS-CoV-2 and in predicting the impact of various interventions. As more data has become available, models have been developed for a number of specific settings, including schools, hospitals, prisons and workplaces [6–9], to take into account the nuances and unique features of each setting. In this paper we present a model of delivery sector that has been uses to assess the impact of various measures that some companies have taken, as well as measures that were under consideration. One unique feature of these settings is the high number of brief contacts that delivery drivers have with members of the public, who themselves may otherwise have very limited contacts. Another feature in the delivery of heavy or large items is the safety requirement for employees to handle and deliver goods in pairs, often requiring prolonged close contact and entry into customers’ properties. Finally, there is still the poorly understood route of fomite transmission that has the potential to be important in this setting, due to the large volume of packages being handled. The model we present considers all of these aspects, and where data is unavailable or uncertain (e.g. for risk of fomite transmission), we consider a wide range of possible scenarios.
We have developed an agent-based network model with stochastic transmission. Therefore, each realisation of the simulation represents a possible chain of transmission within a workplace, and so conclusions can only be drawn from the aggregated results of many simulations. There are commonalities with several models in the literature, including the network models for COVID-19 transmission in workplaces [10]. The stochastic infection and isolation model is similar to other agent-based and branching process models [8, 11]. The model was developed based on a combination of epidemiological data and data and qualitative information gained from consultations with companies in the logistics sector in the UK.
As vaccines are rolled out globally, we expect the most severe impacts of COVID-19 on public health to be curtailed. However, containing the spread of new variants is likely to require good surveillance testing. There has been considerable debate around the usefulness of Lateral Flow Device (LFD) tests that can be self-administered and give rapid results [12–14]. Primarily, this centres around the lower sensitivity of LFD against Polymerase Chain Reaction testing (PCR), particularly at low viral loads [15], and the potential impact of false positives. However, recent data suggests that LFD specificity may be at least 99.9% [16], suggesting that false positives will have a negligible impact. Furthermore, culturable SARS-CoV-2 virus is only found, at most, in the first 8-10 days following symptom onset [17–19], when viral load is higher. This suggests that lower sensitivity tests may still be useful at detecting people when they are most infectious. However, the way tests are performed (e.g. self-administered vs. trained tester) can have an impact on sensitivity [20], plus the method of rollout (e.g. supervised vs. unsupervised testing) can affect the adherence to the testing policy. The model we present accounts for these various factors.
The aim of this paper is to estimate the efficacy of different workplace interventions with a model particularly tailored to the home-delivery sector. We considered several interventions and scenarios based on formal consultations with company representative from this sector. A secondary aim is to estimate the potential impact of presenteeism (working while sick) with COVID-19 symptoms. Flexible or ‘gig’ contracts are common in the home delivery sector, as well as the use of self-employed couriers, all of which are factors associated with increased presenteeism [21], so this is an important factor to consider.
II. METHODS
A. Data Collection and Company Consultations
We carried out recorded consultations via teleconference with representatives from six companies between July and August 2020 (Round 1), and May and June 2021 (Round 2), three of these companies were interviewed in both rounds. Each semi-structured interview lasted 60-90 mins and was based around a set of open-ended questions regarding how the pandemic had impacted on the operations of the business and what measures had been put in place to protect staff and customers. As part of these consultations we asked questions regarding the number of staff working at typical sites and the frequency of contacts between employees and the public. Additionally, two companies provided data on staff numbers and deliveries, which are detailed in Appendix A 1. A summary report was sent to each company for comments and corrections to verify that we had interpreted their answers accurately, and the data correctly. Further details on the consultations are published in [2].
We also used data from an online contact survey aimed at delivery drivers in the UK [article in preparation], which received 170 responses (104 of which were from the workers involved in the delivery of small packages and/or large items). This survey was elective so was not statistically representative. The results of this survey are to be published elsewhere but a few results are utilised in this paper. Namely, only 5.3% reported working while having symptoms of COVID-19 or with a member of their household having a suspected or confirmed case of COVID-19. Conversely, 17.2% reported having isolated with symptoms of COVID-19 or due to a member of their household having a suspected or confirmed case of COVID-19. This suggests approximately 1 in 4 failing to isolate for one of these reasons. For this reason we consider two pisol values (0.5 and 0.9) as ‘low’ and ‘high’ isolation rates, noting the likely caveat of reporting biases. Staff reported large numbers of daily contacts (mean 15.0) at their place of work, which, tallying with the results of consultations, we interpreted as a result of repeated interaction within a work cohort (with only rare random interactions on top). Hence our assumed cohort size for drivers of ≈ 13.
B. Workplace Network Model
In this section we present an overview of the model details, with further details supplied in appendix B. The model we use is a stochastic agent-based network model of disease transmission. The parameters and symbols used in the following section are all described in table I.
The model considers contacts between all employees working in a home delivery depot (i.e. engaged in business-to-consumer delivery or B2C) that has a warehouse and onsite offices. The workplace is populated by 3 groups: drivers, who deliver packages from the warehouse to customers; pickers, who transport and load packages within the warehouse; and office/admin staff, who work in the same building but in shared offices. There exists a pool of ND drivers, NL pickers and NO office staff available for work each day. Workforce turnover is ignored, as it is assumed negligible over the time scales considered, however it may play a role over long time periods.
1. Employee work schedules
The model network consists of all within-workplace contacts between employees, as well as contacts between employees due to house-sharing or carpooling, in order to simulate workplace outbreaks in detail. Thus we assume, unless they share a household, employees only make contact with other employees if they are both at work on that day. We use an idealised model for the work schedule, whereby the number of employees in work depends on the day of the week, this pattern was calculated from data from two UK logistics firms (see Appendix A for details). For all pickers and drivers, we randomly assign consignments (i.e. deliveries/packages) for loading and delivery (as detailed in Appendix A). We assume that each consignment is first handled by pickers, then subsequently by drivers, and finally by the customer. Drivers are the only group of employees that have direct contact with members of the public while on shift. For simulation efficiency, repeat interactions with customers are not considered (as contacts via this route have a very low probability of infection, so double counting of infections is very unlikely), but these contacts are simulated and infection ingress/egress through this route is included in the model.
We also consider the case where drivers and pickers work in pairs (i.e. large goods delivery), we round the number of staff required in these roles to the nearest even number, and then assign pairings randomly each day. One intervention simulated is fixed pairings; in this case, these are assigned a priori and we pick the pairs working on a given day at random from those available. A pair is unavailable if either worker in that pair is isolating, therefore this intervention is always used alongside “pair isolation”, where one member of the fixed pair isolates for the same period as their partner (whether or not they develop symptoms).
2. Workplace contacts and infections
Infections are modelled to occur via four routes; face-to-face (F2F) contact with infectious individuals, indirect contact via sharing a space with infectious individuals, household sharing, and fomite transmission via goods handling.
The model generates direct F2F contacts between employees through three different mechanisms, summarised in table II. Table II also lists the parameters for the different contact routes simulated. Contacts made via these routes are assumed to be dominated by face-to-face transmission (see appendix A 4 for simulation details).
Indirect aerosol-mediated transmission is taken to occur on a one-to-all basis. Given the well-ventilated nature of warehouses, we assume that this kind of transmission only occurs in offices, or in lunch/break rooms. Finally, fomite transmission via package handling is simulated as a decaying random process, such that the probability of onward transmission depends on the time between package handling events by infectious and susceptible individuals. The details of how these routes are simulated is given in appendix A 4.
The transmission routes between different groups are illustrated in fig 1.
C. Individual characteristics: viral load, infectivity, and test positivity
Viral load trajectories are generated from the individual level data in [26]. The algorithm to generate individual viral load and infectiousness profiles is described in further detail in appendix B and in [27] is available at [28]. It is summarised in appendix A 3, with values given in table III.
D. Simulation algorithm
We simulate two types of scenario:
Point-source outbreak: A single index case is chosen and we assume that there are no other introductions during the simulation. All other employees are susceptible at the simulation start (i.e. zero prevalence). The simulation terminates when there are no infectious cases remaining. This type of scenario is modelled in section III A and in appendix sections C 0 a–C 0 c.
Continuous-source outbreak: No index cases are chosen initially and introductions occur randomly (Poisson process) based on the community incidence and prevalence in March-June 2020. The simulation runs for a fixed time window, and the number of customer contacts and packages delivered follow the pattern of demand experienced during that period of time (see appendix A). This scenario is modelled in section III A 1.
In the point-source outbreak scenarios, in order to define a ‘successful’ outbreak we arbitrarily set a threshold of a final attack rate of 5%. Note that we choose this as it is a low-threshold, like the epidemiological definition of an outbreak as a single linked secondary case. However, by defining it as a percentage of workplace size this makes the results from the two different settings more comparable. Therefore, if R − 1 > 0.05(ND + NL + NO), where R is the number of recovered individuals at the end of the simulation, then we record this simulation as a successful outbreak. The fraction of simulations where a successful out-break occurs is then used as an estimate of the probability of an index case resulting in an outbreak.
For continuous-source outbreaks, there is random ingress of new cases, so instead we compare the number of workplace infections (ignoring introductions) as well as the number of isolation days to measure impacts on productivity. Introductions can occur in these simulations through two routes:
Community ingress: Each susceptible individual in the workplace has probability I(t) of being infected outside of work, where I(t) is the community incidence at time t.
Customer ingress: For each delivery a driver makes, there is probability P (t) that the customer is currently infectious, where P (t) is the community prevalence at time t. When a susceptible driver interacts with an infectious customer, there is probability pcust = 1 − exp(−ciβF 2F τdoorstep) of an infection.
This is a very simple model of case ingress and does not account for household structure, the geographical/individual variability in the wider population, or repeat deliveries to customers.
The testing strategy we model here is non-directed mass testing, i.e. all employees are tested regularly every τp days. A random day in the period [1,τp] is drawn as the first test day, and all subsequent test days follow sequentially τp days after the previous. Following a positive test, an individual cannot be tested again for τpause days after their positive test. Other testing strategies may be beneficial, particularly if looking to reduce the burden on employees or because of affordability, and we address some of these in the discussion.
The simulation follows an SIR-type structure, such that individuals who have previously been infected cannot be re-infected. This is a reasonable assumption over the timescales of up to 3 months that we consider here. An example visualisation of a single realisation of the simulation is shown in figure 2. The source code for the simulations can be found at [29].
III. RESULTS
Appendix C contains results pertaining to the baseline transmission dynamics in the two work settings that we focus on. This includes an investigation of changing parameters that affect the mixing rates in the settings, including changing work cohort size and close-contact work in pairs. This shows that, in the parcel workplace, we predict that office size and occupancy is a more important potential factor in workplace outbreaks than transmission between drivers at the workplace, even though office workers are in the minority (section C 0 a). In the large-item delivery workplace, close-contact working pairs (primarily delivery pairs, who share a vehicle for much of the day) were predicted to be very important, and keeping these pairs fixed had a significant impact on reducing workplace spread (section C 0 b). Finally, we also present the effect of presenteeism, which in this model we define as workers with symptomatic COVID-19 attending work, which we find can have a notable effect on transmission, particularly when coupled with other measures to isolate close-contacts of symptomatic individuals (section C 0 c).
Furthermore, appendix D contains a sensitivity analysis of some of the pertinent baseline model parameters, to show how these affect predicted transmission rates.
In this following section, we focus on the impacts of testing and the combination of different workplace interventions to analyse their potential effectiveness.
A. Impact of mass testing on point-source outbreaks
Given the long incubation period of COVID-19 (compared to flu) and the significant proportion of asymptomatic cases, regular mass testing has been proposed and deployed in various settings to screen asymptomatic and pre-symptomatic cases. Figures 3 and 4 show the reduction in outbreak probability resulting from testing at different frequencies with different test types in the parcel and large-item delivery work settings respectively. Overall, the results show that even with a 2-day turnaround, PCR has a bigger effect on transmission that LFD tests, due to the increased sensitivity, however the effect of LFD tests is still significant considering their relative low cost. Note this estimate for the sensitivity of LFD tests is lower than some previous models [20, 30, 31] as it is based on data from self-swabs in a real-world setting rather than lab data.
In each figure, two cases, representing idealised behaviours, are shown. In the first case testing is voluntary meaning 90% of people do 60% of the required tests on average, while the other 10% do no tests (3(a) and 4(a)). This therefore reduces the potential benefits of testing. In the second case testing is enforced(3(b) and 4(b)) meaning that all workers test and report their results. This is the theoretical maximum effect that we could expect testing to have.
Comparing figures 3 and 4 shows that testing has a greater proportional impact in the large-item delivery delivery setting. Furthermore, in this case testing can still have a noticeable effect even when performed as infrequently as 14 days. With total compliance to testing, the probability of outbreaks in the large-item workplace model is reduced by approximately 50-60%, and around 40-50% in the parcel model by twice weekly LFD tests (which have been deployed in other sectors). The reason for this difference is the same as in the previous section on presenteeism, namely the combination with the fixed-pair isolation intervention in the large-item delivery setting.
To conclude, we have found that regular testing, particularly in combination with close-contact isolation, can have a very significant effect on workplace transmission. Any testing intervention needs to be weighed against potential costs, at low community prevalence the vast majority of tests are likely to be negative, and those that are positive are more likely to be false positives and so the intervention may not represent good value for money. Alternatively, at high community prevalence, testing and close-contact isolation could result in many isolations, some of which are only precautionary, which can have a huge impact on business. The latter case is not well described by the point-source outbreak considered in this section, as introductions into the workplace are more likely to occur in quick succession. Therefore in the following section we look at the impacts of a range of interventions in the case of a continuous-source outbreak.
1. Impact of interventions in a real-world context
In this section we model each workplace in the context of realistic community SARS-CoV-2 incidence rates. We used incidence rates inferred from deaths and hospitalisations in the UK during the period 1st March 2020 until 31st May 2020 (see appendix A 2). We then applied ran simulations with different interventions in place, for each scenario we added an extra intervention to the ones applied before, the interventions are:
Symptom isolation only: People who develop symptoms self-isolate with probability pisol = 0.5.
Improved isolation: To mimic the impact of pandemic messaging, isolation probability is increased to pisol = 0.9.
Distancing: All F2F interactions, except those involved in pair work, have interaction distance x = 2m.
Cohort Size Reduction: In the parcel delivery setting, the number of cohorts for all job types is doubled.
House share isolation: All employees who share a household isolate when one self-isolates.
Fixed-pairings: In the large-items delivery setting, driver and picker pairs are fixed and both self-isolate if one self-isolates.
Office WFH: Office staff do not enter work, they only make contact other employees if they share a household.
Testing: Twice weekly lateral flow testing is introduced for all employees.
Enforced testing: Testing becomes mandatory so no tests are missed.
Car share isolation: If a person travels to work with someone who self-isolates, they self-isolate.
Cohort isolation: If one member of the cohort isolates, all people in the cohort isolate.
In the model, introductions due to customer interactions only had a small but noticeable effect meaning that drivers were slightly more exposed than other employees (mean 0.11 introductions per driver for both work settings, averaged over all scenarios vs. 0.09 for other staff respectively). Nonetheless, over the period, around 10% of the workforce is infected purely due to the imposed prevalence and incidence.
Figure 5 shows the cumulative impact of interventions on secondary cases and isolations in the parcel delivery workplace. The interventions are applied in approximately the sequence that was reported by companies we consulted. The intervention “Distancing” increases all “cohort” and “random” contacts to 2m interactions and has a large effect. Reducing cohort size and office staff working from home (“Office WFH”) have a big impact on reducing transmission since this model predicts that outbreaks are most likely to start in this group. Interventions beyond “enforced testing” are predicted to increase isolation levels without much greater impact on transmission, particularly “cohort isolation” which likely causes a great deal of disruption despite these groups being unlikely to be infected. Note this becomes a much more viable option though if cohorts are much smaller, which is one major benefit of reducing cohort size if possible. Comparing the two graphs in 5 we see that is a slightly more efficient to have contact reduction measures in place before adding isolation-based measures, as these reduce the number of workers who will need to isolate. When isolation measures are implemented alone, we see an increase in the predicted number of isolations even though the relative reduction in transmission is similar.
The impact of interventions for the large-item delivery workplace is very similar (see figure 6). We see that the “fixed pairings” intervention (which includes pair isolation) has a marked effect on transmission. The extra benefit gained from testing is clearly visible too, but again isolation measures beyond this appear have little further effect.
To conclude, figures 5 and 6 demonstrate some of the trade-offs for different intervention measures in terms of their impact on transmission and their impact on the number of isolating employees. Certain interventions act to reduce both (social distancing, Office staff WFH) but potentially have other costs for business/feasibility issues that need to be considered. When there are employees that still need to be in close-contact (e.g. driver and picker pairs in this model) the combination of fixed pairings, pair isolation, and regular testing is highly effective for reducing transmission. However, in cases where there is high rates of household-sharing (or, more generally, any contacts between employees outside of work during isolation) then this can continue to drive transmission between employees and is difficult to distinguish from workplace transmission.
IV. DISCUSSION
In this paper we have developed a stochastic model of SARS-CoV-2 spread in small/medium size workplaces. The contact patterns simulated were designed to represent warehouses/depots in the home-delivery sector, particularly those focusing on B2C delivery. To our knowledge this is the first model to consider SARS-CoV-2 transmission in this sector specifically. While the parameterisation of these models has significant uncertainty, we have been able to test the relative impact of various interventions that companies in this sector deployed to reduce SARS-CoV-2 transmission over a range of scenarios and parameter regimes.
The results predict that workplace transmission in this sector is modest, due to drivers working alone most of the day, and warehouse staff working in large open environments. Without any interventions there is predicted to be a small risk to customers for an individual delivery, but in workplaces undergoing an outbreak, home-installation of items can pose a significant risk to customers without interventions. The companies we consulted discontinued home-installation during in the spring of 2020, but later re-introduced it with social distancing measures. Parcel delivery companies have switched to “contactless” delivery, meaning that signatures are no longer required, essentially eliminating the only route of transmission to customers. Overall, this suggests that this sector played a key role in reducing community transmission of SARS-CoV-2, as it allowed people to stay at home during periods of high-prevalence. Quantifying this impact is more difficult though as the counterfactual situation (i.e. how people would have behaved if this sector failed to keep up with increased demand) is unknown.
Safeguarding the key workers in this sector was a broader challenge and companies reported implementing multiple measures based on government guidelines and their own judgement. A key result of this paper is the finding that the potential impact of self-isolation and testing policies is dependent on how much transmission between employees occurs in work vs. at home. Thus, identifying high-risk contacts (due to e.g. shared accommodation or work tasks requiring prolonged close-contact) is very important and forms the basis of contact-tracing interventions. Workplaces have an extra advantage over contact-tracers in that they have control and knowledge over some of the contacts that employees are required to make in their line of work. Therefore, high-risk contacts can be limited by using fixed pairs for large-item delivery, reducing then number of people sharing an office, and reducing the occupancy in shared spaces. This then allows efficient isolation policies to be implemented based on knowledge of the limited number of high-risk contacts people have made (e.g. workers are given paid isolation leave if they share accommodation or work in a delivery pair with an employee who has tested positive or reported COVID-like symptoms). Figures 12 and 4 show that this combination can be extremely effective in small workplaces.
Rates of presenteeism, which in this case we define as those who do not self-isolate when they develop symptoms, has been shown to be much less likely if fully-paid sick leave is offered [32, 33]. Therefore, in order to be effective, such isolation policies could incur considerable costs to a business as well as reducing productivity. Similarly, company-backed testing and isolation interventions will incur further costs and mean that asymptomatic cases are detected, potentially resulting in even more isolations. Therefore, companies may be apprehensive about deploying such strategies. Here, we showed that these strategies are much more efficient when close-contacts of index cases can be identified and isolated too, particularly when such contacts are necessary for the job. This can mean that chains of transmission are quickly shut down, and outbreaks are much less likely to occur. Furthermore, combining these measures with social distancing, WFH, and similar interventions that reduce transmission (e.g. masking) can reduce the number of isolations since workplace outbreaks become less likely.
The model developed has several limitations which are important for the interpretation of the results presented. First, the contact model has been developed based on a mix of quantitative (survey data, staff numbers and demand levels) and qualitative data (consultations). Novel insight was gained by speaking directly to representatives from the sector, but their position was not objective and so there may have been implicit biases in the descriptions of the nature of workplace contacts and some potential routes of infection contacts could have been missed. Furthermore, simplifying assumptions, such as all contact durations being identical for the same mode of contact, mean that this model is idealised compared to reality. Third, the transmission rate and the modifiers used for different types of contact are uncertain, and is based on a combination of peer-reviewed [23, 24, 34] and non-peer-reviewed literature [22]. Improvements to this transmission model from the adaption of more mechanistic modelling approaches that predict explicitly the infectious dose associated with different modes of contact [24, 34–41], as well as updating with data on new variants and vaccines, will mean that this model could be applied to numerous future workplace scenarios to test the impact of different non-pharmaceutical interventions.
There are also some complicating factors we choose to ignore in this model, that may be important to consider when interpreting these results. First, we do not model severe illness, which can impact results by increasing the time away from work of individuals with COVID-19. Second, we do not model the complex relationship between interventions and behaviour. It is possible that as more interventions are introduced, adherence with other interventions wanes so the expected impact of combined interventions may not be as high as predicted. This behavioural change is difficult to predict, and so would need to be monitored by companies to gauge whether interventions are working as expected.
One major benefit of the model presented here is that it incorporates the dynamics and variability in individual viral load, and simulates its impact on test sensitivity and infectivity. This, means that the correlation between test-positivity and infectivity is incorporated, meaning that impacts of these interventions can be more accurately estimated. Thus we were able to estimate not only the effect on average transmission rates, but also the frequency of rare superspreading events. This has highlighted the importance of stacking interventions that reduce transmission through different mechanisms. The source code is open access [29] and the underlying network transmission model is malleable enough to be applied to any small closed populations.
In conclusion, this paper has shown that the multiple interventions put in place by the logistics and home delivery sector during the early stages of the pandemic are likely to have reduced the risk of community and workplace transmission by safeguarding customers and staff. The availability of lateral flow tests is another valuable layer of protection that could be added, and that this would most effective when combined with isolation measures that target the most high-risk contacts.
Data Availability
All relevant data for reproducing the results are given in the manuscript. Open access source code for the models used is also cited in the manuscript.
https://github.com/CarlWhitfield/Viral_load_testing_COV19_model
https://github.com/CarlWhitfield/Workplace_delivery_transmission
Appendix A: Sector-specific data collected and derivation of parameters
1. Staff and delivery numbers
To obtain realistic figures for the workplace size and number of consignments, we received data from a parcel delivery company and a logistics company that deliver large items, both in the UK.
From the parcel company, we received company-wide figures for number of consignments (i.e. delivery drops) and number of drivers on-shift from 01/01/2020 to 01/09/2020. First, we fitted a linear generalised additive model (GAM) to the number of consignments, in order to extract a smoothed curve P (p)(t) representing the demand over this period (see figure 7(a)). Next, we used negative binomial regression to extract the weekday dependence for both the number of drivers and the number of parcels. The model means (µ for drivers and λ for deliveries) were parameterised as
such that w(t) returns the day of the week (Monday = 0, Tuesday = 1, etc.) and t = 0 is a Monday. The weekday dependence for each quantity is shown in 7(b). Finally, we were also provided with weekly data for each site (again for number of consignments and number of drivers). We again fitted a negative binomial regression to this data over the same period, with categorical parameters for each site to account for the differences in demand. The mean of this model was parameterised as
where u is the time in weeks (e.g w/c 01/01/2020 is u = 0), W (p)(u) is the weekly demand (i.e. the total demand P (p)(t) summed over the days in week u) and s is the site number (indexed from 0 for all Ns sites).
Similar data was provided by a large-items logistics company. In this case, we were provided with daily data for the number of consignments and number of vans for each site in the company from 01/04/19 to 01/09/20. First, we fitted a linear GAM to the total number of consignments over all sites to extract a smoothed demand curve P (l)(t) (shown in 7(a)). Next, we fitted a negative binomial model to both account for weekday effect and site, with mean
We assume that, in the absence of COVID isolations, staff numbers on average follow the same daily pattern as driver numbers shown in figure 7(b). The selection of employees in work for a given job j on a given day is a three-step process:
All workers in COVID-related isolation are removed.
Of those available, each is assigned a random absence with probability pabs and removed from the availability pool.
Of those still available, each has probability exp of being in work on day t, where is the negative binomial parameter corresponding to the day with highest occupancy. In this notation i = l if in a large-item delivery workplace and i = p is the parcel delivery workplace.
The fitted weekday coefficients are used to model the day-to-day variation in the number of deliveries, so that for i ∈ {l, d} We use different values for the background demand based on the type of simulation.
Scenario simulation: For the parcel delivery workplace, the weekly negative binomial model can be used to approximate the number of deliveries over this period such that where is the index of the median site in the parcel delivery dataset. In the large-item delivery workplace then is the index of the median site in the large-item delivery dataset.
Outbreak simulation: We assume fixed pre-pandemic levels of demand (approximately equal to the scenario simulation case at March 1st 2020). In this case we use . in the parcel and large-item workplaces respectively. In this way, represents the number of deliveries that the parcel workplace completes in a week. For outbreak simulations we fix this at
In the case of outbreak simulations, we fix where NP = 3000 deliveries per week for the parcel company and NP = 210 deliveries per week for the logistics company.
To make the number of staff in work each day consistent with the data for a median workplace, we defined the total number of drivers as where pP D is the average number of parcels delivered per driver per day under normal conditions (85 for parcel delivery and 15 for large-item). The number of warehouse and office staff were then defined relative to this.
In the in situ scenarios, we use the negative binomial fitted demand such that where is the demand coefficient for the median site and i ∈ p, l indicates whether it is the parcel or logistics company.
2. Community incidence and prevalence
The estimates community incidence data used for the continuous-source outbreak simulations are shown in figure 8. They are based on data on hospitalisations and deaths from 01/03/2020 – 31/05/2020 inclusive. Community prevalence was simply estimated as the cumulative sum of the incidence from the previous 10 days inclusive.
3. Viral load, infectivity, and test-positive probability
The viral load and infectivity trajectories are generated at random based on the mechanistic model fits in [26]. In order to use this data, we fitted the viral load trajectories with a piecewise exponential model where Vp is the peak viral load (in copies/ml), r is the exponential growth rate, d is the exponential decay rate, and τp is the time since infection to the peak viral load. These parameters are drawn from a multivariate distribution that is inferred from the data in [26]. The parameters of this distribution are given in table III.
There is probability psymp of each individual developing symptoms that would be identified as symptomatic COVID-19. For those that do, symptom onset time is drawn from a gamma distribution (see table III for parameter values), which is truncated so that symptom onset occurs within ± 2 days of peak viral load [42, 43].
In order to model infectivity, we use “infectious virus shed” as fitted in [26] as a proxy. The infectious virus shed is a Hill function of viral load We fitted a bivariate lognormal distribution to the random parameters of this model (Jp, h), given in table III. Note that, in [26], the magnitude of Jp is given in arbitrary units. In order to make this consistent with the transmission rate used here, we scale the distribution of Jp such that the population mean infectivity over the average duration of infection is 1.
The probability of testing positive is also assumed to be deterministically linked to viral load. To model PCR testing, we use a logistic function with a hard-cutoff to account for the cycle-threshold (Ct) cutoff To model LFD testing, we use data from regular LFD and PCR testing in the social care sector in the UK. Since this data is based on positive PCR results, we assume PLFD(Vk) = theta(Vk)PPCR(Vk). The parameters of these testing models are given in III.
Finally, individuals then have probability pisol of being an “isolator”, meaning that if they develop symptoms or test positive they will isolate. For “non-isolators”, we assume these individuals do not participate in voluntary testing, but do participate in mandatory testing. In voluntary testing, those that do participate have probability pmiss of missing a test.
4. Transmission rates
A contact of type i between an infectious individual k and a susceptible person kl is assumed to have transmission rate where βF 2F is the average transmission rate for face-to-face contact while talking at 1m separation, ci is the transmission modifier for the type of contact, taken from [22], τi is the contact duration, and xi is the contact distance in metres. Skt and Jk are the susceptibility and infectivity of individuals kl and k respectively on day t.
The baseline transmission rate βF 2F is estimated using multiple sources. First, in the [22], the probability of transmission while talking at 1m separation for one hour is given as 0.06 for the original SARS-CoV-2 virus and 0.15 for the delta variant. This suggests a βF 2F of 0.07−1 or 0.15−1 respectively. Test and trace data [23] suggests a workplace/school secondary attack rate of ∼ 3% between close contacts, defined as either < 1m for more than a minute, 1-2m for > 15minutes or sharing a vehicle/skin-to-skin contact. If we take the average contact to be between 1-2m for 15-60min this gives a plausible range of βF 2F =0.03– 0.24h−1. Finally, in [24] the probability of infection while talking at 1m separation for 15 mins is ∼ 0.03, suggesting a βF 2F of 0.13h−1. Therefore, we take βF 2F = 0.15h−1 in this model as a “best guess” from this range of estimates.
For household contacts, we use the transmission rate (per day) of where βHH = 0.07. This was set empirically using 106 realisations of infectivity trajectories so that the mean probability of household infections was 40% (approximately twice the observed household attack rate from the original SARS-CoV-2 strain [4, 22, 48, 49]). Comparing to (A16), we see that this means the household transmission rate is approximately equivalent to 30 mins of F2F conversation indoors at a distance of 1m per day.
For package mediated fomites, we use a two-step process such that a fomite infection between workers occurs with probability where t1 < t2 is the time the package is handled by an infectious individual k and t2 is the time the package is handled by susceptible individual kl. The parameter λ is the half-life of the viable virus on the packaging, so the first exponential term represents this decay process. The second term is the same as used for all other contacts except that βF OM is effectively the probability of transmission if an infected person handles a package that is then immediately handled by another employee.
Package handling is assumed to occur as follows
A picker (or picker-pair) k handles all of their lk packages at times uniformly distributed in a time window of τL hours.
Each infectious packages n is then assigned a random integer rn between 1 and P, which then determines which driver (or driver-pair) they get assigned to, such that is assigned to driver k.
Drivers are assumed to handle all of their dk packages twice. First they handle all at time τL, and second at times distributed between τL and τL + τD (where τD the time is takes to complete all of the deliveries).
Finally, for shared-spaces (not including F2F contact), we use the dose-response relationship originally derived for SARS-CoV [50] and applied to SARS-CoV-2 [51]. Assuming that the number of aerosolised infectious plaque-forming units (PFUs) exhaled per hour is n0, the concentration in the surrounding air (assuming fast mixing) will be where V is the volume of the room and f is the air exchange rate (air changes per hour, ACH), and λ = ln(2)/t1/2 is the decay rate of the aerosolised virions (such that t1/2 is their half-life).
Therefore, the average concentration of infectious PFUs in the room while occupied by a single infectious person for duration τ is To convert this to a transmission rate, we assume that an average individual has an inhalation rate (volume of air inhaled per hour). Their average dose, if sharing the space for the full period τ, is then .
The exponential dose-response model used in [50] has the same form as our infection model, so the transmission rate between individuals sharing the same space for period τ is given as such that the probability of infection via this route is pss = 1 − exp(−βssτ). The infectious dose IDe is highly uncertain, as is the dose exhaled n0. However, the ratio n0/IDe is approximated in [34] as ∼16–18 per hour at peak infectivity. Given that, in our model, peak infectivity is ∼4 times mean infectivity, we use n0/IDe = 4.0 here.
In order to more directly compare this to the F2F interactions, we derive an “effective interaction distance” that is implied by this transmission rate. I.e. how far away would the two contacts have to be in a (no-talking) F2F interaction to experience an equivalent transmission rate in our model. This is given by xeff = 1 + log(βF 2F /5βss), since βF 2F /5 is the F2F transmission rate at 1m distance with the “no-talking” modifier. The parameters used to estimate xss are given in table IV.
a. Modes of face-to-face contact
We include several routes through which face-to-face contacts occur in the workplace, which are summarised in table II. Below we give some further details about how these are simulated:
Cohort contacts: Daily F2F interactions between employees working the on the same shift and in the same area. For details of the rationale between how these cohorts are assigned, see Appendix A. It is assumed that, each employee contacts all others within their cohort for the same duration, such that the total face-to-face contact time for any given individual in a cohort is fixed. I.e. if a cohort consists of M people, each has a F2F contact with the other M − 1 people in their cohort for τcoh/(M − 1).
Random interactions: On a given day, each person with job j at work will randomly contact another employee with job jl with probability
The factor ρD is to account for the reducing mixing in drivers due to the time they spent not at the workplace during a shift. Note, pDD = ρD and not ρ2 because we are assuming that drivers are on site at approximately the same time as each other.
Explicit workplace pairings (includes large-item handling, pair delivery, and pair delivery routes): These are necessary for execution of picking and delivery of large items. When these pairings exist, all drivers and pickers are paired at the start of each shift, and this contact is guaranteed to occur during the shift. If the “fixed pairings” intervention is applied, then staff in pairs always have the same partners, otherwise the pairs are assigned randomly at the start of each day.
Appendix B: Simulation algorithm
At initialisation, the following model features are generated
Agents are assigned job roles and susceptible status.
The cohort contact network is generated as follows:
Within each job role, each agent is assigned a cohort at random.
For each driver cohort, one member of warehouse staff (a ‘picker’) is assigned to that cohort (representing a supervisory role).
Within each cohort, all agents are connected by edges with identical weight, which are active only when both agents are in work.
The house-share network is generated as follows:
NH = l(ND + NL + NO)/(1 + H)l households are generated.
Each agent is assigned randomly to a household.
Within each household, all agents are connected by edges with identical weight that are active every day.
The car-share network is generated as follows:
NC = lNH/(1 + C)l car-shares are generated.
Each household is assigned to a car-share, with all agents in that household assigned to the same car-share.
Within each car-share, all agents are connected by edges with identical weight that are active only when both agents are in work.
If simulating an outbreak scenario, an index case is infected on a random start day from Monday to Sunday. Upon infection, an agent is assigned the following:
Viral load and infectivity trajectories.
Symptom onset time.
Adherence to self-isolation (true with probability pisol).
If testing: a test positive probability trajectory (section A 3).
The main simulation loop is executed for each day of the simulation, and proceeds as follows:
Update infectious state of all individuals moving any to ‘Recovered’ status who have reached the end of their infectious period.
Perform testing (if a testing day). For all positive tests generate an isolation time from the current day as lτd + u01l where u01 ∼ U (0, 1) is a number uniformly distributed beteween 0 and 1, and l.l indicates rounding to the nearest integer.
Update isolation status for any who are due to isolate on this day.
Randomly generate any new introductions due to community incidence.
Select employees in work, all others are excluded from any workplace contacts for the rest of the day.
Generate all successful workplace infection events.
Generate customer infections and introductions.
Remove any infections where the infection target is not in the ’susceptible’ state. For any viable targets that are subject to more than one successful infection event, select the recorded infection event at random.
Record all infection events, and for every individual infected change their status to ’infected’ and their infection time to the current day.
Increment the day and return to step 1 unless the maximum number of days has been simulated. If an outbreak simulation and all individuals are in the ’recovered’ state, terminate the simulation.
Appendix C: Baseline modelling results
At baseline (assuming all default parameters in table I, but high symptomatic isolation rate of pisol = 0.9), our model predicts that the average number of cases the index case will infect is 0.181 ± 0.005, 0.934 ± 0.016, and 2.58 ± 0.03 when the index-case is a driver, picker, or office worker respectively. In the large-items delivery work setting, it is predicted to be 0.833 ± 0.012, 0.935 ± 0.015, and 1.608 ± 0.019 when the index-case is a driver, picker, or office worker respectively. Note that the error bounds here are 1.96 standard errors of the mean, and do not account for parameter uncertainty (only stochasticity of the simulations).
Given the numbers of workers in each job role (see I) this means that if the index case is selected at random we expect an R number of ≲ 1 in the parcel work-setting, and ≈ 1 in the large-item setting. Therefore, without any interventions our model predicts that these setting are very likely to see small outbreaks. Given the uncertainty in the underlying parameters, these baseline characteristics need to be considered as a “best guess” case. Bearing this in mind, there remains a great deal to be learned from quantifying the relative impact of various non-pharmaceutical interventions.
The average number of secondary cases in each job role and via each contact route is summarised for each workplace in figures 13 and 14.
a. Parcel delivery: Role of cohort size and mixing
Figure 9 shows the effect of increasing the number of work cohorts (and thereby reducing cohort size) in the parcel delivery work setting, as well as changing the rate at which staff switch between cohorts. Cohort-size only makes the large difference in the case where the index case is a member of office staff. This is because separating into different offices reduces aerosol transmission. In the case of driver/picker cohorts, the underlying model assumptions mean that increasing the cohort size does not change the total amount of contact time (it is simply spread over more contacts). Therefore cohort size only makes a difference to outbreak probability when there is significant enough transmission via this route that saturation effects can occur. In our model this is the case for picker staff, where the effect is small, but not drivers.
Similarly, figure 9 also shows that increasing the rate at which employees move between cohorts (fc) has very little effect on outbreak probability. Drivers are predicted to be much less likely to cause an outbreak due to their reduced time spent in the workplace.
Thus we predict that, for these settings, cohort size effects are most important when considering shared indoor spaces such as offices. In well ventilated spaces, where workers are generally spread out and only make contact intermittently, the dominant factor to consider is the total F2F contact time that each employee has with colleagues (not the number of distinct contacts). Note this ceases to be true for longer interactions where transmission risk is higher, as will be demonstrated in the following section.
b. Large-item delivery: Role of pair work and longer customer exposure
We now consider the case of workplaces that deliver large-items (e.g. furniture, white goods, etc.), such that drivers and pickers work in pairs to perform their tasks. This means they are in close-contact with one another for long periods of the day, increasing the risk of transmission. In figure 10 we consider 4 intervention scenarios:
No intervention: pairings are picked randomly at the start of each shift.
Fixed pairings: pairs are fixed and isolation of one partner automatically triggers isolation of the other.
Windows open: transmission rate between drivers sharing a cabin assumed to be the same as outside transmission rates, rather than inside.
Combination of scenarios 2 + 3.
Figure 10(a) shows the sizeable impact that the fixed pairings intervention has on the baseline dynamics. In particular, when the index case is a driver, the baseline dynamics are dominated by pair transmission (as shown in supplementary figure 14), and more-so for driver pairs due to the time spent in a shared cabin, which we assume to be higher risk than close-contact transmission in better-ventilated settings, such as the warehouse or on the doorstep. The randomly switching pairs mean that all drivers are connected by this transmission route, resulting in a larger number of secondary cases. Fixed pairings (scenario 2) close off these chains, meaning drivers can infect at most one person via this route, and that person then cannot infect anybody else via this route.
On the other hand, when the index case is a warehouse worker, or when we assume ”open windows” for drivers, figure 10(a) shows that changing from random to fixed pairings does not have a large effect on transmission. This indicates that, daily contact within picker-pairs (or driver-pairs with open windows) are not predicted to be so high-risk to see a saturation effect. Note that this result is sensitive to the underlying parameterisation, and so has significant uncertainty associated with it.
Figure 10(b) compares the risk to customers related to a point-source outbreak with a driver index case in the large-item delivery workplace vs. the parcel workplace. It is clear that, for each individual customer there is increased risk (due to the longer contact time, chance of home entry needed). Nonetheless, the total number of customers likely to be infected is comparable in the two cases, owing to the large number of deliveries carried out by each driver in the parcel delivery setting. In both cases, even without interventions, this probability is predicted to be small (0.18 ± 0.005 and 0.35 ± 0.008 customers infected per simulation in the parcel and large-item settings respectively). However, in the event of a successful workplace outbreak (i.e. more than 5% of the workforce infected), this increases to 0.90 ± 0.10 and 1.01 ± 0.02 customer infections for the parcel and large-item delivery workplaces respectively.
We conclude that pair-work can pose a high-risk of transmission in these settings, in which case fixed pairings in combination with pair isolation is an effective measure to break transmission chains and reduce the probability of a workplace outbreak. Generally, we predict that home delivery is unlikely to seed a significant number of new infections in the community, except in cases where there is a workplace outbreak. Thus, extra measures to protect customers (contactless delivery, mask wearing indoors etc.) will have likely made a meaningful impact when prevalence in the workplace was high.
c. Effects of presenteeism
In this model we define presenteeism as working with symptoms of COVID-19 (as we do not model other sickness absences explicitly). On average, 50% of the employees in our simulations develop symptoms relevant for isolation, and if they isolate then they do not attend work for the following 10 days (but household transmission between employees is assumed unaffected). Note that we do not model the effects of illness severity (e.g. prolonged symptoms after recovery, hospital stays, or even death).
Figure 11 shows the effect of increasing self-isolation rates among symptomatic individuals in both work settings, for various values of H, the house-sharing factor. Outbreak probability reduces linearly with isolation adherence, and the proportional effect is larger when house-sharing is rare.
Figure 12 shows the impact of pisol in the large-item delivery setting. The impacts are similar in this workplace to the parcel setting. Comparing the impact with and without fixed-pair isolation policies shows that this interventions results in a slightly larger relative reduction in outbreak probability as presenteeism is reduced. This shows that the combination of this policy with efforts to reduce presenteeism has a larger combined effect than either measure independently, but only in the case of driver index cases. However, we do not see any impact when the index case is a picker or an office staff member, reflecting that the transmission dynamics in these groups are somewhat indepedent of the drivers in our model.
To conclude, adherence to isolation measures are predicted to make a substantial difference to transmission dynamics, particularly in the case where close-contact pairs are both isolated when one develops symptoms. We saw a less-strong effect for household isolation in our model even when a large fraction of workers share accommodation, this is because we assumed that this measure does not prevent household transmission, which was the dominant transmission mode in this scenario. This stresses the importance of messaging around measures to reduce within-household transmission as well isolation measures.
Appendix D: Sensitivity to assumptions and parameter choices
1. Baseline transmission rates by job role and contact type
The average number of secondary cases from a points source outbreak (using the baseline parameters used in section C) is broken down by job role of the index case, the secondary case and the contact route for each work setting in figures 13 and 14.
2. Transmission parameters
The baseline transmission rate for F2F contacts is derived from multiple sources, but these are largely observational studies and involve a number of assumptions. It is therefore useful to observe how some of the key results are affected if the transmission rate is higher (e.g. due to more transmissible variants, or due to error in the transmission rate used). Figure 15 shows how the histogram of secondary case numbers changes as the transmission rate βF 2F is increased (note we assume that shared-space transmission also increases proportionally). Assuming a random index case, we see a transition from an exponential-type distribution to a multi-modal distribution, as the workplace R-number exceeds 1 and most introductions result in several secondary cases. In these cases, NPIs can have a much more significant impact, as they can push the workplace R-value back to below 1.
Second, the fomite transmission parameter is an unknown in the model, so it is important to see for what values it must take for fomite transmission to be significant. In figure 16 we see that increasing the fomite transmission by a factor of 10 makes a relatively minor difference to the outcome of point-source outbreaks originating with a picker (it has a much smaller effect for other index cases). However, the change is substantial when increasing by another factor of 10, suggesting that, if this were the case, it would be a dominant mode of transmission and a significant contributor to overall workplace transmission. Note that, even though the effect is much less significant in the large-items workplace, the βF OM value might actually be larger in this environment, since these items require more force to move (which can affect surface deposition rates) and may require more overall close-contact with the item.
3. Number of employees
The influence of overall workplace size is largely controlled by how we assume the contact networks scale in warehouses with fewer/more staff. In order to test this we use simulate workplace sizes corresponding to the upper and lower quartiles of the datasets in section A 1. We scale the workplace size by factor 0.5 ≤ r ≤ 2 in the parcel setting and 0.5 ≤ r ≤ 3 in the pairs setting. We assume that the parameters scale as follow
while the definition of pc remains the same (and so decreases with workplace size).
Figure 17 shows the results of point-source outbreak simulations in the two settings across the range of feasible workplace size scalings. We see that the number of secondary cases resulting from the outbreak is increases with workplace size, although this is mostly explained by the increase in office cohort size (the sharp reduction occurs when the number of teams is increased). This is because, for the parameters used here, most outbreaks become extinct. Therefore this means that the outbreak probability, as defined here, actually drops as workplace size increases, but this is an artifact of choosing 5% as the threshold for an outbreak.
It is worth noting though that, with greater mixing rates than simulated here, larger workplaces could pose a larger risk of outbreak because of the presence of a larger susceptible population.
4. Mixing rates
We vary the amount of mixing in the workplace by changing the parameter pc. As pc increases, not only do more contacts occur but also contacts between different job roles become more frequent. The parameter ρD remains fixed at 0.1, so drivers still have less contact than the other job roles. Figure 18 shows how the number of secondary cases in a point-source outbreak increases as pc is increased for each workplace.