would estimate more realistic scenarios, in addition to allowing the testing of mitigation protocols focused on this heterogeneity. In this work, we present a hybrid diffusion model of COVID-19 with the objective of evaluating these complex dynamic scenarios and its difference from conventional homogeneous models based on ODE.
Most epidemiological models currently found in the literature are based on the work of Kermack and McKendrick [4] , known as the Susceptible-Infected-Recovered (SIR) model. The SIR model is appropriate in studies of infectious diseases that confer permanent immunity and consider individuals within the considered region according to their disease status over time. Susceptible individuals, denoted by S, refers to those who do not have the disease at time t, but can be infected through contact with an individual in the infected state, I. The recovered state, R, are those who have recovered from the disease or are dead and can no longer transmit the disease [4, 11] .
The SIR model is described by the following ordinary differential equations:
Where S represents the population susceptible to the disease, while I and R represent the infected and recovered population. and are rates that represent the dynamics of disease spread, the disease transmission rate and the recovery rate, respectively.
Eq. (1) describes the variation in the susceptible population S(t) at time t. Since the infected population I(t) will have contact with susceptible S(t) per unit of time, the number of susceptible people will decrease by S(t)I(t)
individuals and S(t)I(t) susceptible individuals will become infected at time t, eq. (2). At the same time t, as each infected individual will recover at a rate , a total of I(t) individuals will recover, as indicated in Eqs. (2) and
(3).
The simplest epidemics SIR model has a high flexibility for including generalizations and has been usually used as a first step in the study of several diseases. The evolution to more complex models is the natural course in the quest to understand important characteristics of the infection progression over time. These modifications generally include heterogeneities, subdivisions of S, I and R compartments, and recovery rates for population subgroups, or even in the definition of these parameters through the social interaction network of studied individuals [12, 13] .
In addition to the limitations already mentioned in the standard epidemiological model, the equations that describe the dynamics of the disease consider a structure of random pair wise interactions among individuals in the studied population. However, real populations often have specific structures under which these individuals interact. This structure composed of subgroups of individuals and their connections constitute a network where infectious diseases tend to spread [11] .
The information regarding the contagion dynamics between municipalities presented in this work is the result of a hybrid model composed of an epidemiological mathematical model of SIR type, which estimates the disease dynamics for each municipality, and of a flow network among municipalities that estimates the spread of cases among them. Previous studies have shown the importance of the transport flow network in the spread of viruses like dengue [14] .
The hybrid model consists of joining two modeling techniques: a set of ordinary differential equations (ODE) and
a network population diffusion model. At each time step, the model solves multiple ODE for the different municipalities within a state, considering as initial conditions the actual number of cases and recovered for each municipality. After that, the flow network is used to distribute populations of susceptible, infected and recovered among pairs of municipalities and a new population is estimated for the next step. An illustration of the algorithm is presented in Figure 1 . When we include inter-municipal flow interference, equations 1, 2 and 3 can be rewritten respectively by 4, 5 and 6, where compartment variations due to internal contagion, described by EDO, also depend on the variation caused by the flow of individuals.
Were Sᵢ,Iᵢ and Rᵢ represent the containers of the SIR model for each municipality i, βᵢ is the transmission rate adjusted for municipality i, {Kᵢ} is the set of municipalities neighboring municipality i, ΔSⱼ and ΔIⱼ are the flux of susceptible and infected population from municipality j to i. In this way, the sum represents the total number of individuals in each container that is due to the flow network.
In order to illustrate the effects of heterogeneity in multiple scales we used data from Brazil, a country full of heterogeneity, with 8.516.000 km 2 and multiple mitigation policies for different geographical scales. Figure 2 shows the average weekly flow of people in Brazil, corresponding to the full set of road, air, and water transport networks were considered throughout the territory in 2016 [15] . in the country. The value of used for all municipalities was 1/6 [16] . We use the municipalities as the smallest unit in the model, as they are the smallest subdivision in which case notifications were released and in which IBGE's displacement flows are reported. The source codes for the parameter setting can be found at [17] and for the hybrid model at [18] .
Simulations were made by adjusting the cumulative case curves from the occurrence of the first case in each municipality (in Brazil, on February 24, 2020) for different periods that will be described throughout the results. The initial conditions used in the hybrid model were the number of active cases, recovered and susceptible individuals in each municipality on the last day reported in the data set and the respective contact rates i.
J o u r n a l P r e -p r o o f
The results presented in Figure 3a) show the time series of the SIR compartments for Brazil resulting from the ODE and hybrid model. We notice that the hybrid model has a delay in saturation and a significant decrease in the active cases peak, showing that the heterogeneity factor acts as a mitigating mechanism in the spread of the epidemic. In figure 3 b) we see a comparison between the adjusted curve and the real data for the city of Salvador, explicitely a good agreement between the adjusted and the real data. 
We evaluated the effect of the network topology by comparing the model results using the observed network with a random network with the same number of edges and nodes and the same weight distribution. Figure 5 shows that randomizing the network, the infected curve presents a different behavior for each territorial scale. At the city level, the infected curve is affected similarly to Figure 4 , that is, the randomization process only seems to modify the initial configuration of the infected cities that influence the growth rate of the curve in the city of Salvador, thus changing the date of the peak. In the state curve, the results are close to those presented in the simple ODE model, that is, the randomized network that represents an intermediate topology between 
The 
The results presented in figure 6a indicate the occurrence of the flattening effect of the active case curve as the inter-municipal flow is reduced. Decreasing the flow between the municipalities causes a de-synchronization between the most populous cities, which are the first to become infected, and the other municipalities. This can even lead to a curve with two peaks, reducing the maximum number of active infected and delaying the contamination of new municipalities.
In spite of presenting a significant reduction in the number of cases, the generalized reduction in the flow of people among all pairs of municipalities can represent a high social and economic cost. However, we know that there are great differences in population size and connectivity degree of the municipalities, so that they are expected to act with greater or lesser impact in mitigating the problem. This naturally leads to investigating the effect of applying a drastic reduction in the flow of people to the most populous municipalities, which we call a flux blockage.
In this sub-section, we evaluated the effect of the 95% reduction in the inter-municipal flows directed to the 10% of the most populous municipalities in Bahia, maintaining the flow between the other municipalities without restrictions. We evaluated the effect of applying this strategy in the initial and advanced stages of the epidemic, corresponding to 20 days (26/Mar) and 56 days (01/May), respectively, after the first occurrence. In this analysis, we used the same contact rate values i (based on data available until Mar 26, 2020) for both scenarios, since the objective was to assess only the effect of the intervention time and not the evolution of the rates of contact that occurred due to the mitigation policies in each municipality .
In Figure 7a , the case curves with and without flux lockdown intervention scenarios are compared. In the initial stage (Block mar/26), we can notice a reduction of 29% of the peak in cases curve in addition to the occurrence of a second peak. When applying the same strategy at an advanced stage (56 days after, Block may/01), the second peak disappears and there is a reduction of only 9%. This effect occurs since, in the late scenario, many municipalities have already been infected and the envisaged de-synchronization between more and less populous municipalities no-longer exists, as we can see in the graphs of Figure 7b . 
Our results show that the SIR hybrid model leads to important variations in the pandemic behavior predictions when compared with the SIR ODE model. This is due to three important characteristics incorporated in the model: the spatial heterogeneity for the infected populations, the different local contact rates and the network of people flow. In conventional ODE models there is no spatialization of data, that is, the population of all municipalities is treated as a single homogeneous set. In addition, the different cultural, economic, and political elements of the municipalities establish different internal mechanisms for the epidemic growth, represented by their different rates of transmission. This initially promotes a reduction in the cases growth with the spread of the infected population throughout the network. 
One of the measures adopted by official bodies worldwide was to reduce the flow of people circulating during the pandemic. Studies show that the lockdown [19] When we look at the number of infected cities per day, our results are close to those found by Pujari and Shekatkar [9] for India. In this study, applying a hybrid model with time delays incorporated in its flow network, it was found that most Indian cities were infected in the first 100 days of the beginning of the epidemic. As seen in 
The application of mathematical models to describe the dynamics of infectious diseases plays an important role in policy decision-making in the field of Health Sciences. Typical mathematical models, through non-linear systems of ordinary differential equations, aim to predict and provide a guide for disease prevention and control.
However, most of the outputs of these models compromise the policy making process. This is due to the fact that these methods may not properly take into account the population heterogeneity, here reflected on the spatial variability of the susceptible population.
In this work we present an innovative way to analyze data and parameterize the SIR model accommodating spatial 