key: cord-0644888-vvnanf3m authors: Grave, Mal'u; Viguerie, Alex; Barros, Gabriel F.; Reali, Alessandro; Andrade, Roberto F. S.; Coutinho, Alvaro L.G.A. title: Modeling nonlocal behavior in epidemics via a reaction-diffusion system incorporating population movement along a network date: 2022-05-09 journal: nan DOI: nan sha: 3098b7470b5536fde39008990a2cd165cc60bd2f doc_id: 644888 cord_uid: vvnanf3m The outbreak of COVID-19, beginning in 2019 and continuing through the time of writing, has led to renewed interest in the mathematical modeling of infectious disease. Recent works have focused on partial differential equation (PDE) models, particularly reaction-diffusion models, able to describe the progression of an epidemic in both space and time. These studies have shown generally promising results in describing and predicting COVID-19 progression. However, people often travel long distances in short periods of time, leading to nonlocal transmission of the disease. Such contagion dynamics are not well-represented by diffusion alone. In contrast, ordinary differential equation (ODE) models may easily account for this behavior by considering disparate regions as nodes in a network, with the edges defining nonlocal transmission. In this work, we attempt to combine these modeling paradigms via the introduction of a network structure within a reaction-diffusion PDE system. This is achieved through the definition of a population-transfer operator, which couples disjoint and potentially distant geographic regions, facilitating nonlocal population movement between them. We provide analytical results demonstrating that this operator does not disrupt the physical consistency or mathematical well-posedness of the system, and verify these results through numerical experiments. We then use this technique to simulate the COVID-19 epidemic in the Brazilian region of Rio de Janeiro, showcasing its ability to capture important nonlocal behaviors, while maintaining the advantages of a reaction-diffusion model for describing local dynamics. The outbreak of COVID-19, which started in 2019 and is still continuing, has caused unprecedented disruption in terms of both human cost and economic damage. In order to better understand the dynamics of the disease spread, with the hope of ultimately improving policy and public health outcomes, there has been an explosion in the study of mathematical modeling of infectious disease. These models have taken many forms, and a comprehensive survey of the literature is beyond the scope of the current work. However, we note that many different approaches have been used to model the epidemic, including machine-learning and data-driven approaches [33, 29, 17, 9, 7] , models using a classical compartmental approach, together with parameter estimation techniques [11, 26] , delay differential equations [15] , partial differential equations [5, 6, 31, 32] , network-based methods [23, 21, 22, 5] , as well as agent-based [34, 20] and multiscale models [4] . We note that the cited works represent just a small sample of the total literature, and further that the various approaches discussed are not mutually exclusive. For a recent survey, see [21] . The basis for most mathematical models of infectious disease is rooted in the compartmental theory of Kermack and McKendrick [19] (see [8] for a more modern treatment). This approach divides a given population into different states corresponding to disease status, with the underlying equations governing the evolution of the population composition in terms of disease-state. This mechanism forms the backbone behind the classic S-I-R (susceptible-infected-removed) model and its variants: for instance, the S-I-S (susceptible-infected-susceptible) and S-E-I-R (susceptible-exposedinfected-removed), among others. While the initial model presented by Kermack and McKendrick is quite general, and presented in the form of an integral renewal equation, under certain assumptions, namely a mass-action transmission mechanism and exponentially-distributed sojourn times, we obtain the familiar systems of ordinary differential equations (ODE) widely used today. These systems should be recognized, however, as special cases within a far more general framework [19, 8] . We note that many models, including agent-based and evolving network models, even if not expressed as integral or differential equations, still generally utilize some sort of compartmental structure [28, 34] . Such compartment models have historically been applied with great success in the modeling of infectious disease, and offer a clear mechanistic description of epidemic progression. However, in their most widely-utilized ODE formulations, relevant factors may not be taken into account. Among the most important such factors, and the subject of this work, is the introduction of spatial information into the epidemic model. Many approaches have been proposed to describe epidemic evolution in space as well as time. A common approach is to further stratify the compartmental structure by geographic location, such that different compartments also correspond to different geographic areas, often delimitated by political jurisdiction. For COVID-19, this approach has been utilized in e.g. [23, 11, 7, 22] and offers the advantage of computational efficiency and relative ease of implementation. Additionally, the spatial evolution is not inherently local, and nonlocal processes (ie, transmission between individuals in distant regions) can be easily incorporated. However, a continuous description of spatial dynamics with this approach is impossible, and as the spatial resolution increases, the number of necessary compartments can quickly become intractable. In contrast, partial differential equations (PDE) incorporating the familiar compartmental structure have been offered as an alternative. Most commonly, models of this type incorporate evolution in space via a reaction-diffusion equation [24, 18] . Works incorporating PDE can be seen in [31, 30, 17, 5, 6, 12, 14] . Most of these used a reaction-diffusion process, with [5] introducing a multiscale approach utilizing convection for nonlocal processes and diffusion for local processes. In [27] convection is also employed to account for population mobility. PDE models offer the large advantage of a true, spatially continuous description; however, they too have downsides. From the practical point of view, they are generally more difficult to implement and require greater computational resources than ODE models. From the modeling point of view, the ability of a linear diffusion model to accurately describe the spatial evolution of an epidemic in human populations, given the nature of human movement, is questionable. While many of the COVID-19 related applications of these models mitigate this problem somewhat through, for instance, nonlinear diffusion and depensation effects [31, 30, 17, 6, 12, 14] , the model still behaves locally. Although [5] does make an excellent attempt at a nonlocal description, ultimately, despite the presence of multiple scales of spatial evolution and convection, the underlying process is still local, as contagion is still spread along with the intermediate points between origin and destination. It is well-known that bilaplacian [24] or fractional diffusion [10] operators can describe nonlocal dynamics. Some works have indeed applied such techniques to COVID-19 (see [16, 1] ); however, the complexity of these models in terms of computational effort and stability, as well as their still-developing theory, makes their widespread adoption, at least at the current time, difficult. Further, the nonlocal movement in human populations is far from random; indeed, it often follows a predictable pattern [2] . In contrast to these approaches, we propose herein a combination of the reaction-diffusion model first proposed in [31] with network-based methods, akin to those proposed in [23, 21, 22, 5] . Similarly to the modeling paradigm introduced in [5], we employ a PDE model and postulate that diffusion occurs at the local (i.e. within-city) level. However, rather than defining nonlocal dynamics via convection, we instead model such dynamics as sources and sinks within a network. We embed a spatial network structure within the underlying PDE system, and define a network via transport operators between different areas in the spatial domain (i.e., inter-city mobility). These operators then allow portions of the population to, for instance, move from one city to another without diffusing the disease along the way, as is typically the case with vehicular transport. We believe that such a method offers a reasonable, easily implemented, stable approach that maintains the advantage of a PDE-based model for local dynamics while also incorporating important nonlocal effects. We outline the article as follows: we first introduce and provide a brief explanation of the underlying equations. We then discuss the idea behind the population transfer network and its mathematical formulation, as well as important details regarding its numerical implementation. We then provide numerical examples on an idealized problem, as well as a realistic problem in the Brazilian state of Rio de Janeiro, to provide evidence of the model's ability to describe the relevant local and nonlocal dynamics found in human epidemics. The governing equations are based on the a spatio-temporal SEIRD model, presented in [31, 30, 13] . We consider that only the susceptible and exposed compartments are able to move along the network. The system of equations becomes: where s(x, t), e(x, t), i(x, t), r(x, t), and d(x, t) denote the densities of the susceptible, exposed, infected, recovered, and deceased populations, respectively. β i and β e denote the contact rates between symptomatic and susceptible individuals and asymptomatic and susceptible individuals, respectively (units days −1 ), α denotes the incubation period (units days −1 ), γ e corresponds to the asymptomatic recovery rate (units days −1 ), γ i the symptomatic recovery rate (units days −1 ), δ represents the mortality rate (units days −1 ), and ν s , ν e , ν i , ν r are the diffusion parameters of the different population groups as denoted by the sub-scripted letters (units km 2 persons −1 days −1 ). Finally, S s and S e are population transfer operators, defining the movement of populations within the population transfer network as sources and sinks in the reaction-diffusion system. The precise definition of these operators will be discussed in detail in the following section. The system of equations (1)-(5) is supplemented with proper boundary and initial conditions. The movement occurs between paired regions, ie, nodes that share an edge in the underlying network. Movement away from one region results in movement into the other region. In this sense, in a given time, when persons move away from a region, the population transfer operator acts as a sink term; likewise, when persons move into a region, the population transfer operator acts as a source term. In this section we define the population transfer operators S within a network. For simplicity, we consider a single compartment and a network consisting of two nodes; extending the definitions to larger networks and other compartments follows immediately. We consider the equation (1), defined in a domain Ω, with two distinct subsets Ω 1 ⊂ Ω, Ω 2 ⊂ Ω such that Ω 1 ∩ Ω 2 = ∅. We define the directed, weighted, and without self-loops transfer network as N = G(Ω, η), where Ω and η represent, respectively, the set of nodes Ω i and edges η j i between nodes i and j. For the time being, let N consist of only two nodes, i.e. i = 1, 2. We also assume that N is a time-varying graph (TVG) [25] , where the edges η ij depend on time according to the following assumptions (note this dependence on time is assumed if not explicitly denoted): This is equivalent to stating that the probability that an individual moves from Ω 1 to Ω 2 (or vice versa) over T 2 1 time units is p 2 1 . For example, if one assumes that this probability is distributed exponentially with parameter η 2 1 , we find that: which after some straightforward algebra gives: We now define the population transfer operator along the network S s (N ). The dependence of S s (N ) on the network N is assumed and not denoted explicitly in hereafter. Such an operator is not unique, and in this section we outline the properties that an operator of this type must satisfy. We split the operator into two components, a donor operator D s which transfers a population from Ω i to Ω j at a rate η j i and the corresponding receiver operator R s which receives a population from Ω i into Ω j at a rate η j i . S s is then given as the sum of these components: (and its individual components) is a linear operator in [ · ], although it may have additional dependencies on x, any of the compartments, and/or t in general, and need not be linear in these variables. Both conditions are physically motivated; the first ensures that the quantity of persons transferred away from given node Ω i must not exceed the total persons present in Ω i . Additionally, this condition is necessary to guarantee the coercivity (and thereby well-posedness) of the variational problem corresponding to (1)-(5) (see e.g. [3] ). The population conservation condition is necessary to guarantee that the transfer operator only models movement among the different regions, and does not change the overall quantity of persons (or mass) in the system. We now discuss the donor and receiver operators separately. In the movement of populations from Ω 1 to Ω 2 (and vice versa) at the rate η 2 1 (and η 1 2 ), the donor operator D s defines the removal of persons from the origin node. We propose the following, simple definition of D s : where χ Ωi (x) are characteristic functions: Assuming the η j i < 1, it follows immediately from Ω 1 ∩ Ω 2 = ∅ that (I + D s )s > 0, ensuring that the non-negativity condition is satisfied. The receiver operator R s receives the populations taken from the donor region by the donor operator and distributes them in the receiver region, and has the following form: The kernel functions K Ωi (x, t) are defined such that: and Ωi K Ωi (x, t) dΩ i = 1 ∀t. Theorem. The transfer function defined by (9), (11) satisfies the mass-conservation property. Proof. One may verify this with straightforward computation: Transfer operators of the type introduced above are non-unique, and any operator satisfying the necessary properties may be applied in practice. Besides the clearly problem-dependent η j i , the definition of the kernel functions K Ωi (x, t), which corresponds to the distribution of the transferred population within the receiver region, is the most important component in constructing S s . As shown in the preceding sections, as long as the conditions (12)- (13) are satisfied, then one may define an acceptable transfer function. However, these conditions are quite non-restrictive, and many such functions will satisfy them. Defining K Ωi is therefore far from straightforward, and problem-specific considerations, computational constraints, and available data are likely to inform such choices in practice. Below we list two possible, obvious definitions. Uniform distribution. A possible definition is given by: This defines the simplest possible transfer function, which distributes the population uniformly among the receiver domain. It is computationally simple and, therefore, the approach applied in the simulations of this work. However, it is questionable how well this definition captures reality in most instances. Receiver-proportional distribution. Another possible, natural choice is: which distributes the population in proportion to the existing population density in Ω i . This is likely more realistic than the uniform distribution; however, it requires the additional evaluation of an integral, must be updated in time, and introduces a nonlinearity into the problem. Depending on how one chooses to discretize the term (16) in time, it may become necessary to apply some nonlinear iteration scheme at each time step. We note, however, that such schemes are generally necessary when solving (1)-(5) anyway, and this may not necessarily represent a large additional computing cost. In order to check the implementation of the source and sinks algorithm, we use a rectangular domain that is divided into two regions. The dimensions of each region are 1 × 1 for Region 1 and 2 × 1 for Region 2. We start with 1000 people in each region, i.e., Region 1 has 1000 people/km 2 and Region 2 has 500 people/km 2 ( Figure 1 ). The objective is to seasonally send people from Region 1 to Region 2 and people from Region 2 to Region 1. We define a weekly pattern in which people leave Region 1 on Fridays and go back home on Sundays. For instance, assuming day 1 is Monday, people leave Region 1 to Region 2, and on day 7 people leave Region 2 to Region 1. The movement is defined to be 20% of the population in each region and assumed to be distributed evenly throughout the day. Hence: We define K Ω1,2 as in (15) . We use different meshes to study population (mass) conservation. Three meshes are structured with different element mesh sizes (25×50, 50×100, and 100×200). The fourth mesh is unstructured with 18,231 elements of sizes varying between the bigger and smaller elements of the fixed meshes. All discretizations with linear triangles. In Figure 2 we show the coarsest structured and the unstructured meshes. As we are not considering diffusion in this example, the only change in population distribution is due to S s . Therefore, it is possible to determine the analytical solution for this problem. In Figures 3 and 4 we plot the total populations in each time-step for Region 1 and 2, considering the four meshes. Figures 5 and 6 show the relative error between the analytically-computed number of people during each time-step for Region 1 and 2, for each of the four meshes. Finally, Figure 7 shows the total population for each mesh and Figure 8 shows the total population (mass) relative error between the simulations and analytical solution. We briefly note that, as we keep all of the percentages fixed at 20%, we observe population drift over time, with greater numbers of people being moved each time. This occurs due to the fact that, after sending 20% of the population in Region 1, when we send 20% of the population in Region 2 back to Region 1, this population now includes the original persons transferred. Noting that 20% of 1000 is 200, we may avoid this problem by defining the transfer rates instead as: sdΩ2 · day −1 if Sunday, 0/day otherwise. (18) While defining η in such a manner requires integral evaluations, we note that these are already computed in order to define (11) . Depending on the mesh refinement, we start the simulation with a different amount of total population, close to 2000. This happens because of the way that the initial conditions are set. The more refined the mesh, less difference between the real value and the one provided by the simulation. As we increase the refinement, we decrease the errors of conservation of the total population during the simulation. However, when we look to the populations at each region, we see that the unstructured mesh returns bigger errors than the structured ones. In [14] we have simulated the spread of COVID-19 in the state of Rio de Janeiro, Brazil. We have shown difficulties in reproducing the virus spread dynamics in some counties, particularly in Cabo Frio (Region 2) and Campos dos Goytacases (Region 3), for example (see Fig. 9 for their locations). In that simulation, there were no initial infected or exposed populations in these counties (Fig. 10) , making diffusion the only possible way for the virus to reach those areas. Given the large distance between them and the hotspots near the city center of Rio de Janeiro (Region 1, respectively 157/279 km distance from Regions 2 and 3), as well as the large sparsely-populated areas between them (see Fig. 11 ), population-weighted diffusion is not able to properly represent these dynamics. Note that, for all simulations considered in the following, parameters unrelated to the population transfer network are exactly as considered in [14] . Cabo Frio (Region 2) is a coastal county, where there are seasonal visitors during the weekends. Therefore, we define η 2 1 as people going from Rio de Janeiro to Cabo Frio on Fridays and coming back on Sunday, as in the previous simple test problem. The main difference is the percentages of people that comes and goes from each city. Since Rio has 6.7 million inhabitants and Cabo Frio has about 200 thousand, movement of the same amount of persons will correspond to different percentages of the population. We assume that about 30,000 people have this weekly pattern movement, or We apply the same rate to both s and e compartments. We note that, due to the previously-discussed population drifting effects, diffusion among the population, as well as the fact that these η j i are applied to both s and e, they may not correspond to exactly the same number of persons moved at each time period. Given that s is several orders of magnitude larger than other compartments, the simplified rates defined above do not create serious issues in general. Indeed, these differences in number and composition of the population transferred in time may in fact be more realistic, as we would not expect such numbers to remain constant. Nevertheless, one may also define the η j i similarly to (18) in order to ensure that population movements remain consistent in time. Another region that was underestimated in the previous simulations for similar reasons to Cabo Frio is Campos dos Goytacazes (Region 3). We may then postulate movement between Rio de Janeiro and Campos dos Goytacazes, but consider a different weekly pattern, as Campos dos Goytacazes is not a vacation area. Rather, many people in the oil industry work there for one week in the offshore oil platforms, returning to Rio during the next week. We assume that there are about 50,000 people who move in this manner, corresponding to 0.75% the population of Rio de Janeiro, and 10% the population of Campos dos Goytacazes, and adopt the convention of denoting even-numbered and odd-numbered weeks in the obvious way. Then, analogously to (19) : 0.0075 · week −1 on odd-numbered weeks, 0 · week −1 otherwise , η 1 3 = 0.1 · week −1 on even-numbered weeks, for both η j,s i and η j,e i . Lastly, to complete the network, we also consider movement between Campos dos Goytacases and Cabo Frio. We define this movement analogously to (19) ; however, rather than moving approximately 30,000 people between the regions, we move about 10,000, corresponding to about 2% of the population of Campos dos Goytacases and 4.45% of the population of Cabo Frio: In total, we simulate eight different cases, considering different networks as presented in Table 1 and Figure 12 . The mesh of the Rio de Janeiro region contains 11,632 elements and 6,082 nodes, and, although unstructured, the element size is very close to uniform, with a mean area of 1 km 2 and 95% of elements between 0.5 and 1.5 km 2 . Simulations are done considering linear elements. In Figures 13, 14 and 15 we show the comparison of cumulative deaths between the observed data and the simulated cases. For each case, we are considering diffusion. Hence, we do not have an analytic solution as in the previous example. However, we do know that the total overall population in the region should be conserved. In Figure 16 we show the total population relative error with time, which is substantially under 1% throughout the time period. When considering isolated movement between Cabo Frio and Campos dos Goytacazes (Case 3) we see very little contribution to the spread of the virus in these counties because they both have no infected or exposed at the beginning Movement between Case PDE RJ and CF RJ and CG CF and CG 1 Table 1 : Definition of the movement between Rio de Janeiro (RJ), Cabo Frio (CF) and Campos dos Goytacazes (CG) for each case. of the simulation, and no nonlocal propagation of the virus occurs. However, when incorporating movement from Rio de Janeiro, which has many exposed and infected in the early stages of the epidemic (as in Cases 1 and 2), we see a marked improvement in the simulated agreement with reality. In particular, we observe that the nonlocal movement defined by the population transfer network enables the epidemic to spread in these regions. We also note that, as the population in Rio de Janeiro is much larger than in Cabo Frio or Campos dos Goytacazes, the deceased population in Rio shows little change when introducing nonlocal population movement. To compare the effects of the population transfer network, Figure 17 shows the spatial distribution of the deceased population at t=180 days for Case 1 (no population transfer) and Case 8 (full network incorporating population transfer between Rio de Janeiro, Cabo Frio, and Campos dos Goytacazes). In Rio de Janeiro, we see that the purely diffusive model (Case 1) already works well, and the introduction of the population transfer network does not cause any sort of deterioration in accuracy (Case 3). In contrast, in Cabo Frio and Campos dos Goytacazes, the purely diffusive model fails to recreate the observed dynamics. However, upon the introduction of the population transfer network, in Cabo Frio, the results closely match the observed data. In the case of Campos dos Goytacazes, the simulations now match measured data more closely, as the deaths are no longer close to zero; however, they are still a fair deal lower than the observed reality. This may be due to our currently simplistic assumptions about the nature of the population movement; under more realistic assumptions, these results may improve. simple, defining an operator properly representing these dynamics within a reaction-diffusion PDE is nontrivial. We discussed the proper definition of such an operator that both respects sensible physical constraints and maintains the mathematical well-posedness of the underlying system. In particular, the operator must be non-negative and population (mass)-conservative. We then demonstrated conditions that ensure that the transfer operator satisfies these important properties, and proposed some natural candidates for its definition. We then performed two numerical experiments with different goals. The first was designed to confirm that our numerical implementation of the network transfer operator demonstrated the necessary conservation and non-negativity properties as outlined in Section 3. We found that our implementation generally respects conservation, with numerical errors below one percent for all meshes considered, and tending to zero with mesh refinement. Unsurprisingly, we also found that structured meshes generally demonstrate superior conservation properties. Our second numerical experiment was a proof-of-concept designed to demonstrate the potential of this approach on a realistic problem. We considered a three-node network of three distinct regions in the state of Rio de Janeiro in Brazil. Without considering population transfer along the network, the PDE model is able to accurately capture the dynamics in the Rio de Janeiro region. However, it is unable to account for the contagion in the Cabo Frio and Campos dos Goytacazes regions, due to the inherent lack of nonlocal dynamics. We then showed that, under reasonable assumptions, introducing a network structure endowed with a population transfer operator is able to account for these nonlocal dynamics, improving the model's performance in the Campos de Goytacazes and Cabo Frio regions. There are several important directions for future research. We considered only the uniform distribution of the population transfer operator in the present work; however, it is important to know the sensitivity of the model to this specific definition. Given the observed mesh dependence on population conservation, it may also be beneficial to explore adaptive meshing algorithms that take the geography of the network regions into account, as this may improve performance. Perhaps most importantly, we stress that the simulations shown here were also preliminary and designed primarily to show the potential of this methodology to account for nonlocal dynamics, and the considered network was very simple. However, applying the described approach to more complex networks with more nodes is important. Similarly, we defined the edges in a rudimentary way. In future work, these edges must be defined more realistically, incorporating measured data. Such extensions represent an important step toward the development of a realistic model which may ultimately inform public health decision-making. Fractional-order susceptible-infected model: definition and applications to the study of covid-19 main protease Understanding mass transfer directions via data-driven models with application to mobile phone data Well-posedness for a diffusion-reaction compartmental model simulating the spread of covid-19 A multiscale model of virus pandemic: Heterogeneous interactive entities in a globally connected world Hyperbolic compartmental models for epidemic spread on networks with uncertain data: application to the emergence of covid-19 in italy Least-squares finite element method for a meso-scale model of the spread of covid-19 COVID-19 dynamics across the US: A deep learning study of human mobility and social behavior On the formulation of epidemic models (an appraisal of Kermack and McKendrick) Metapopulation network models for understanding, predicting, and managing the coronavirus disease covid-19 The fractional laplacian operator on bounded domains as a special case of the nonlocal diffusion operator Spread and dynamics of the covid-19 epidemic in italy: Effects of emergency containment measures A new convected level-set method for gas bubble dynamics Adaptive mesh refinement and coarsening for diffusion-reaction epidemiological models Assessing the spatio-temporal spread of covid-19 via compartmental models with diffusion in italy, usa, and brazil. Archives of Delay differential equations for the spatially resolved simulation of epidemics with specific application to covid-19 Novel fractional order sidarthe mathematical model of covid-19 pandemic Bayesian-based predictions of covid-19 evolution in texas using multispecies mixture-theoretic continuum models Numerical simulation of a susceptible-exposed-infectious space-continuous model for the spread of rabies in raccoons across a realistic landscape Containing papers of a mathematical and physical character Covasim: an agent-based model of covid-19 dynamics and interventions Is it safe to lift covid-19 travel bans? the newfoundland story Scaling effect in covid-19 spreading: The role of heterogeneity in a hybrid ode-network model with restrictions on the inter-cities flow Mathematical Biology II: Spatial Models and Biomedical Application Components in time-varying graphs Suihter: A new mathematical model for COVID-19. application to the analysis of the second epidemic outbreak in Italy A comprehensive spatial-temporal infection model Progression and transmission of hiv (path 4.0)-a new agent-based evolving network simulation for modeling hiv transmission clusters Coupled and uncoupled dynamic mode decomposition in multi-compartmental systems with applications to epidemiological and additive manufacturing problems Simulating the spread of covid-19 via a spatially-resolved susceptible-exposed-infected-recovereddeceased (seird) model with heterogeneous diffusion Diffusion-reaction compartmental models formulated in a continuum mechanics framework: application to covid-19, mathematical analysis, and numerical study A data-driven epidemic model with social structure for understanding the covid-19 infection on a heavily affected italian province An integrated framework for building trustworthy datadriven epidemiological models: Application to the COVID-19 outbreak in New York City An agent-based computational framework for simulation of global pandemic and social response on planet x The work was partially supported by the Brazilian agencies CAPES (Finance code 001), FAPERJ (ALGAC), and CNPQ (ALGAC and RFSA). M. Grave was partially supported by Fiocruz. A. Reali was partially supported by the Italian Ministry of University and Research (MIUR) through the PRIN project XFAST-SIMS (No. 20173C478N).