Applied Mathematics Letters 111 (2021) 106617 Contents lists available at ScienceDirect Applied Mathematics Letters www.elsevier.com/locate/aml Simulating the spread of COVID-19 via a spatially-resolved susceptible–exposed–infected–recovered–deceased (SEIRD) model with heterogeneous diffusion Alex Viguerie a, Guillermo Lorenzo e, Ferdinando Auricchio a, Davide Baroli f, Thomas J.R. Hughes e, Alessia Patton a, Alessandro Reali a, Thomas E. Yankeelov d,e, Alessandro Veneziani b,c,∗ a Dipartimento di Ingegneria Civile ed Architettura, Università di Pavia, Via Ferrata 3, Pavia, PV 27100, Italy b Department of Mathematics, Emory University, 400 Dowman Drive, Atlanta, GA 30322, USA c Department of Computer Science, Emory University, 400 Dowman Drive, Atlanta, GA 30322, USA d Departments of Biomedical Engineering, Diagnostic Medicine, and Oncology, Livestrong Cancer Institutes, The University of Texas at Austin, 107 W. Dean Keeton St., Austin, TX 78712, USA e Oden Institute for Computational Engineering and Sciences, The University of Texas at Austin, 201 E. 24th Street, Austin, TX 78712-1229, USA f Aachen Institute for Advanced Study in Computational Engineering Science (AICES), RWTH Aachen University, Schinkelstraße 2, 52062 Aachen, Germany article info Article history: Received 2 May 2020 Received in revised form 1 July 2020 Accepted 1 July 2020 Available online 15 July 2020 Keywords: Partial differential equations Mathematical epidemiology Compartmental models COVID-19 Mathematical modeling Mathematical biology abstract We present an early version of a Susceptible–Exposed–Infected–Recovered– Deceased (SEIRD) mathematical model based on partial differential equations coupled with a heterogeneous diffusion model. The model describes the spatio- temporal spread of the COVID-19 pandemic, and aims to capture dynamics also based on human habits and geographical features. To test the model, we compare the outputs generated by a finite-element solver with measured data over the Italian region of Lombardy, which has been heavily impacted by this crisis between February and April 2020. Our results show a strong qualitative agreement between the simulated forecast of the spatio-temporal COVID-19 spread in Lombardy and epidemiological data collected at the municipality level. Additional simulations exploring alternative scenarios for the relaxation of lockdown restrictions suggest that reopening strategies should account for local population densities and the specific dynamics of the contagion. Thus, we argue thatyes’;font-family:’Arial Unicode MS’;mso-fareast-font-family:MicrosoftYaHei; font-size:11.0000pt;” > © 2020 Elsevier Ltd. All rights reserved. ∗ Corresponding author at: Department of Mathematics, Emory University, 400 Dowman Drive, Atlanta, GA 30322, USA. E-mail address: [email protected] (A. Veneziani). https://doi.org/10.1016/j.aml.2020.106617 0893-9659/© 2020 Elsevier Ltd. All rights reserved. 2 A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 1. Introduction The outbreak of COVID-19 in 2020 has caused widespread disruption throughout the world, leading to substantial damage in terms of both human lives and economic cost. To arrest the spread of the disease, governments have enacted unprecedented measures, including quarantines, curfews, lockdowns, and suspension of travel. The wide-reaching ramifications of such measures, deemed by many experts as necessary, are driven in part by a lack of clear information about the spatio-temporal spread of COVID-19. Indeed, the absence of reliable data regarding disease transmission has necessarily led to cautious responses. These recent events, which have required important decisions based on forecasts, have demonstrated more than ever the need for reliable tools intended to model the spread of COVID-19 and other infectious diseases [1]. A particularly urgent need is the geo-localization of outbreaks, as this may allow a more effective allocation of medical resources. Several notable models of this outbreak have been presented; indeed, at the time of this writing there are over 1000 COVID-19 articles on MedRXiv, many of which address the modeling of disease spread. Some models aim at offering specific evaluations of policy responses based on the implementation of different social distancing measures with combined compartmental and empirical approaches [2–4]. Rather than adopting a deterministic, mechanism-based model, Zhang et al. employed a statistical approach to analyze the spatio- temporal dynamics of COVID-19 [5]. In Gatto et al. a combined statistical and compartmental approach was employed, in which spatial dependence is addressed by dividing the region of interest (Italy) into local communities connected by a network structure [3]. Here, we propose an alternative approach, using a partial-differential-equation (PDE) model designed to capture the continuous spatio-temporal dynamics of COVID-19. We leverage a compartmental SEIRD (susceptible, exposed, infected, recovered, deceased) model that incorporates the spatial spread of the disease with inhomogeneous diffusion terms [6–9]. The rationale is that the diffusion operator, properly tuned to account for local natural or social inhomogeneities (e.g., mountains, rivers, highways) may describe the local movement of the different populations in a deterministic way, as the limit of a Brownian motion [10]. This is critical to accurately account for information relevant to the outbreak dynamics, such as local population densities, which vary in space and time. While a mathematical description of non-local dynamics is still possible in terms of fractional differential operators [11], we postpone this approach to a follow-up of the present work. Hence, our modeling approach is more appropriate for the local dynamics on mesoscales, such as regions within Italy. Thus, to evaluate the model efficacy, we run a simulation study of the COVID-19 outbreak in the Italian region of Lombardy, which has been severely impacted by the COVID-19 crisis between February and April of 2020 and for which the necessary data was available. Also, the high density of Lombardy’s population and transportation network is specifically suitable to our modeling approach. Our simulations show a remarkable qualitative agreement with the reported epidemiological data. We further explore various reopening scenarios, obtaining contrasting results that highlight the importance of considering local population densities and contagion dynamics. The paper outline is as follows. In Section 2, we describe the SEIRD model. Then, Section 3 addresses the numerical implementation of our model and Section 4 presents the results of the simulation study in Lombardy. We conclude in Section 5 by examining the shortcomings observed from our simulations and discussing the additional work required to improve model accuracy and practical relevance. 2. Model Let Ω ⊂ R2 be a simply connected domain of interest and [0, T] a generic time interval. We denote the densities of the susceptible, exposed, infected, recovered and deceased populations as s(x,t), e(x,t), i(x,t), A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 3 r(x,t), and d(x,t) respectively. Also, let n(x,t) denote the sum of the living population; i.e., n(x,t) = s(x, t) + e(x, t) + i(x, t) + r(x, t). Then, our model is comprised of the following system of coupled PDEs over Ω × [0, T ]: ∂ts = αn − (1 − A/n) βisi − (1 − A/n) βese − μs + ∇ · (n νs∇s) (1) ∂te = (1 − A/n) βisi + (1 − A/n) βese − σe − φee − μe + ∇ · (n νe∇e) (2) ∂ti = σe − φd i − φri − μi + ∇ · (n νi∇i) (3) ∂tr = φri + φee − μr + ∇ · (n νr∇r) (4) ∂td = φd i, (5) where α is the birth rate, σ is the inverse of the incubation period, φe is the asymptomatic recovery rate, φr is the infected recovery rate, φd is the infected mortality rate, βe is the asymptomatic contact rate, βi is the symptomatic contact rate, μ is the general (non-COVID-19) mortality rate, and νs, νe, νi, and νr are diffusion parameters respectively corresponding to the different population groups. Each of these parameters may depend on time, space, or the model compartments. We also consider the Allee effect (depensation), characterized by the parameter A. In this particular setting, the Allee effect serves to model the tendency of outbreaks to cluster towards large population centers. Specific parameter selection as well as initial and boundary conditions are discussed in Section 3. Fig. 1 shows the dynamics of contagion between the compartments in our model. We remark that our model accounts for asymptomatic transmission, which is considered a pivotal driver of the COVID-19 pandemic [4,12,13]. Eqs. (1)–(2) show that exposed asymptomatic patients may transmit COVID-19 to susceptible individuals at contact rate βe. This aligns with recent studies suggesting that patients may transmit COVID-19 almost immediately after exposure [4,12,13]. Additionally, Eqs. (2) and (4) involve a fraction φe of exposed patients who do not develop symptoms and move directly into the recovered population. We also assume that recovered patients are immune, as we do not include any backflow from Eq. (4) to Eq. (1). (We note that this is a current source of debate, but consistent with the existing literature for the time scale of months considered here [14]). The spatial movement over a large population is described by an inhomogeneous random walk, which in the limit tends to a second order differential operator [10]. The diffusivity coefficient is proportional to the population and can be locally adjusted to incorporate geographical or human-related inhomogeneities [7]. 3. Numerical implementation We use a finite-element spatial discretization of the Italian region of Lombardy, consisting of an unstructured mesh containing 30,407 triangles (a mesh convergence analysis was first performed; data not shown). We use the backward-Euler method for time integration and solve each time step fully implicitly with a Picard iteration for stability. The resulting linear systems are solved by the GMRES algorithm using a Jacobi preconditioner. Initial conditions for the subpopulations s, e, i, r, and d in the model are defined by means of Gaussian circular functions centered at the latitude and longitudinal coordinates of each municipality with 10,000 or more inhabitants, and weighted by the municipality’s population size and geographic area. These initial conditions correspond to the data provided by Lab24 [15] from the date 27 February 2020, featuring a severe outbreak in the province of Lodi, and moderate numbers of exposed and infected individuals in the provinces of Bergamo, Brescia, and Cremona (see Fig. 2). We used homogeneous Neumann boundary conditions for simplicity, mimicking a complete isolation of the region. We acknowledge that this is likely unrealistic (e.g. if lockdown orders are relaxed) and will be considered in future efforts. 4 A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 Fig. 1. Flow chart describing the dynamics of contagion between the population subgroups considered in our model. We assume σ = 1/7 day−1, φr = 1/24 day−1, φd = 1/160 day−1, and φe = 1/6 day−1. These values were based on available data from the literature regarding the mortality, incubation period, and recovery time for infected and asymptomatic patients [2,4,12,13,16]. Additionally, we do not consider births or non-COVID-19 mortality (i.e., we set α = 0 and μ = 0, respectively), given the time scale of months in our simulations. The remainder of the model parameters is estimated in a two-step approach. First, we fit a 0D SEIRD version of our model (i.e., consisting of a system of exclusively time-dependent ordinary differential equations and no diffusion terms) to match the temporal dynamics of the outbreak. Then, we iteratively refined these values by means of recursive simulations using Eqs. (1)–(5) to match the spatio-temporal epidemiological data. We use the R2 coefficient and the root mean squared error (RMSE) to assess the goodness-of-fit. Given the uncertainty in the currently available COVID-19 data, we think that parameter estimation aiming at matching the dynamics of all model compartments is not viable. As not every member of the population is tested for infection and asymptomatic cases are known to exist in possibly large numbers, we think that the available data of infected cases might lead to unrealistic parameter fitting. Conversely, the data reported for COVID-19 deaths offer more reliability to calibrate the model parameters. Therefore, we pursue quantitative agreement in the deceased compartment (i.e., d), and qualitative agreement for the rest of the model subgroups (i.e., s, e, i, r). In Section 4, our results will focus on the model forecasts of exposed and infected cases because these are key data for public health officials, e.g., in deciding resource allocation and measures to prevent contagion. We assume βi and βe to be equal, as precise estimates on the relative infectivity levels between the symptomatic and asymptomatic pools are unclear [12]. We define βi,e with decreasing piecewise constant values in time to model the escalation of the lockdown restrictions. Following the results of parameter calibration, we initially set1 βi,e = 3.3 · 10−4 persons−1·day−1 on 27 February 2020, reducing this to βi,e = 8.5 · 10−5 persons−1·day−1 after the first lockdown measures on 9 March 2020, to βi,e = 6.275 · 10−5 persons−1·day−1 after the additional restrictions on 22 March 2020, and to βi,e = 4.125 · 10−5 persons−1·day−1 following the final restrictions on 28 March 2020. Similarly, we assume νs,e,r = 0.0435, 1 We call “persons” the unit measure for the populations. In other papers, this is called “contacts”. A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 5 Fig. 2. Model forecast of COVID-19 spread in Lombardy. (A) Main areas affected by the pandemic in Lombardy. (B) Initially, the main affected areas are Lodi and Cremona and, to lesser extent, Bergamo and Brescia. (C–E) Our model predicts increasing exposures in Bergamo and Brescia. The outbreak in Lodi soon moves north into the Milan metro area, where it further spreads despite the lockdown restrictions. (F) The model also predicts that governmental restrictions eventually succeed in reducing the exposure to the disease, which is faster in Brescia and Bergamo than in Milan. (G) Cumulative curves of infections according to reported data (dots) and simulations (dashed lines) for the three main areas of contagion: Bergamo, Brescia, and Milan. The model has been calibrated to match the data reported for the deceased subgroup, resulting in a forecast of a larger number of infections. To highlight the qualitative agreement of our simulations, we also show the numerical results scaled to match the order of magnitude of the reported infectious data (solid lines). 0.0198, 0.0090, and 0.0075 km2· persons−1· day−1 over the respective phases. The Allee term A is set to 1000 persons, and we fix νi = 1.0 · 10−4 km2· persons−1· day−1 throughout, assuming that symptomatic individuals are largely immobile. 4. Results 4.1. Forecasting the spatio-temporal dynamics of the COVID-19 pandemic in Lombardy Fig. 2 shows the evolving spatial pattern of the COVID-19 outbreak in Lombardy, beginning with exposure in Bergamo, Brescia, Cremona and Lodi. The contagion moves north from Lodi into Milan via the southern suburbs and eventually reaches the city center. We note that although Lodi and Cremona are the most affected areas at the onset of the outbreak, they quickly improve and avoid the explosive growth found in Milan, Bergamo and Brescia. This is consistent with the reported data [15]. However, we found Cremona to be somewhat underpredicted by our model. This might be attributable to the presence of the neighboring city of Piacenza, which shares its metropolitan area with Cremona but is not included in our simulations because it belongs to the adjacent region of Emilia-Romagna. The inclusion of a possible inflow of infected 6 A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 Fig. 3. Effect of the boundary conditions: a spike of the outbreak in Cremona (black arrow) induced by an inflow of individuals from Piacenza, a city in Emilia Romagna southbound of Lombardy. through the southern border with Emilia-Romagna does show a small spike in Cremona, as reported in Fig. 3. In this case, in the absence of available data, we postulated with an educated guess an abnormal influx from the South (100 susceptible persons Km−1 day−1), 1 exposed person Km−1 day−1, 1 infected person Km−1 day−1). This result pinpoints the importance of an accurate calibration of the boundary conditions. In Fig. 2, we demonstrate the remarkable qualitative agreement in the outbreak dynamics between our model forecasts and data in the three main affected areas: Milan, Bergamo, and Brescia. The R2 between the model forecast and the data of infected cases are 0.997, 0.977, 0.976, and 0.998 for all Lombardy, Bergamo, Brescia, and Milan, respectively. We observe that the outbreak emerges in Milan later, where it grows more steadily, eventually becoming the most affected area in Lombardy. We also note that the lockdowns appear to have effectively halted the spread in Bergamo and Brescia. These restrictions notably reduce the spread in Milan, limiting the virus to a linear growth pattern, but fail to stop it. We observe that our simulations predict a larger number of infections than the reported data. This results from using the data for deceased cases for calibration, which is comparatively more accurate than infections (see Section 3). We obtained R2 = 0.972 and range-normalized RMSE = 7.6% for this subgroup. Thus, the difference in predicted and measured infections suggests a lack in the reporting of real cases, probably due to the deficiencies and difficulties of testing a significant sample of the whole living population. However, we also remark that COVID-19 mortality data depend on the currently unknown transmission rates, which emphasizes the importance of qualitative agreement to test novel modeling approaches. To this end, we show that it is possible to rescale our simulation results to accurately match the order of magnitude of the reported infected case data Fig. 2, though we emphasize that this is purely for the purposes of visualization. 4.2. Exploring alternative reopening scenarios We further use our model to assess four illustrative reopening scenarios over four months following 27 February 2020: maintenance of restrictions, relaxation of the lockdown everywhere on 3 May 2020 under two different sets of assumptions, and a combination of maintenance of restrictions in Milan and relaxation elsewhere in Lombardy. We still consider the changes in parameter values induced by the sequential restrictions (see Section 3) and the lockdown relaxation is modeled by setting νs,e,r = 2.175 · 10−2 km2· persons−1· day−1 and βi,e = 9.0 · 10−5 persons−1·day−1 (scenario A), and βi,e = 6.6 · 10−5 persons−1·day−1 (scenario B). Scenario A is a pessimistic scenario, which assumes that the population contact rate is similar to A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 7 Fig. 4. Results of simulations over alternative reopening cases. (A) Cumulative infected cases in Lombardy over the different reopening scenarios. Our simulations suggest that maintaining a strict lockdown outside of Milan offers little benefit; however, keeping lockdown restriction in Milan may prevent explosive growth. (B) Comparison of the cumulative infections in the three largest metropolitan areas (i.e., Bergamo, Brescia and Milan) for the global reopening A (GR A, dashed lines), global reopening B (GR B, dotted lines), and maintenance of lockdown (L, solid lines). early-outbreak levels. Scenario B is more optimistic and assumes that the generally greater public awareness of preventative measures (such as mask-wearing and social distancing) translates to greater success in limiting contact, despite increased mobility. Fig. 4 shows the resulting outbreak dynamics for these four reopening scenarios. Our simulations suggest that relaxing the lockdown restrictions in the entire region may cause severe and rapid growth in the Milan area. However, major urban zones far from Milan (e.g., Brescia and Bergamo) just experience a marginal increase in growth and still show a favorable trend in time. Conversely, if we maintain the lockdown restrictions in Milan and relax them elsewhere, the outbreak shows more favorable dynamics, similar to those obtained for Brescia and Bergamo. Thus, our results suggest that maintenance of lockdown measures in high-population, high-density areas like Milan may be necessary for longer times to effectively arrest the spread of contagious diseases like COVID-19. 5. Discussion We have introduced a compartmental PDE model describing spatio-temporal propagation of disease contagion and applied it to the 2020 outbreak of COVID-19 in Lombardy specifically. Our simulations are intended to provide a proof-of-concept of the potential of PDEs for regional modeling of the outbreak. Nonetheless, they show good qualitative agreement with reality, accurately predicting the outbreak dynamics in different areas and recreating the transmission path in time and space. We then used the model to examine some possible reopening scenarios, which suggested that reopening may be best determined based on local population and contagion dynamics, and not by a one-size-fits-all approach. Our model is in a very early stage, with ample room for improvement. We plan to consider non- constant model parameters and adaptively update them according to measured data usingyes’;font-family:’Arial Unicode MS’;mso-fareast-font-family:MicrosoftYaHei; font-size:11.0000pt;” > 8 A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 modeled by fractional operators can also be included [11]. These considerations may be crucial whenever using the model in larger geographical domains. We would also like to extend our framework into more sophisticated compartmental models including, e.g., hospitalizations, patients in intensive care units, or age and biological sex structures [3,4]. This would further increase the utility of the model, potentially helping decision-makers to determine the allocation of resources among different areas. The present results clearly pinpoint the current standpoints of virologists, emphasizing the need of restrictions. Finally, the socio-economical costs of lockdowns are not included here, but could ultimately be incorporated in future quantitative analyses, e.g., aiming at the comprehensive optimization of pandemic-arresting measures. CRediT authorship contribution statement Alex Viguerie: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Resources, Data curation, Writing – original draft, Visualization, Project administration. Guillermo Lorenzo: Conceptualization, Methodology, Investigation, Writing – original draft, Writing – review & editing, Visualization. Ferdinando Auricchio: Conceptualization, Methodology, Writing – original draft, Funding acquisition. Davide Baroli: Conceptualization, Methodology. Thomas J.R. Hughes: Writing – review & editing, Supervision. Alessia Patton: Writing – review & editing. Alessandro Reali: Writing – review & editing, Supervision. Thomas E. Yankeelov: Conceptualization, Writing – review & editing, Visualization. Alessandro Veneziani: Conceptualization, Methodology, Software, Investigation, Writing – original draft, Writing – review & editing, Supervision, Project administration. Acknowledgments The authors would like to acknowledge the work of Marco Demarziani, Luigi Greco, Isabella Atcha, Kasey Cervantes, Sanne Glastra, Shreya Rana, Stefano Minelli, Chiara Macchello, Simona Petralia, Anita de Franco, and Martina Moschella for their crucial help with data acquisition. References [1] A. Remuzzi, G. Remuzzi, COVID-19 and Italy: What next? Lancet (2020). [2] N. Ferguson, D. Laydon, G.G. Nedjati, N. Imai, K. Ainslie, M. Baguelin, S. Bhatia, A. Boonyasiri, P.Z. Cucunuba, G. Cuomo-Dannenburg, et al., Impact of non-pharmaceutical interventions (NPIs) to reduce COVID19 mortality and healthcare demand, Imperial College London, 2020. [3] M. Gatto, E. Bertuzzo, L. Mari, S. Miccoli, L. Carraro, R. Casagrandi, A. Rinaldo, Spread and dynamics of the COVID-19 epidemic in Italy: Effects of emergency containment measures, Proc. Natl. Acad. Sci. USA (2020) http: //dx.doi.org/10.1073/pnas.2004978117. [4] G. Giordano, F. Blanchini, R. Bruno, P. Colaneri, A. Di Filippo, A. Di Matteo, M. Colaneri, Modelling the COVID-19 epidemic and implementation of population-wide interventions in Italy, Nat. Med. (2020) 1–6. [5] X. Zhang, H. Rao, Y. Wu, Y. Huang, H. Dai, Comparison of the spatiotemporal characteristics of the COVID-19 and SARS outbreaks in mainland China, medRxiv (2020). [6] E.E. Holmes, M.A. Lewis, J.E. Banks, R.R. Veit, Partial differential equations in ecology: Spatial interactions and population dynamics, Ecology 75 (1) (1994) 17–29. [7] J. Keller, L. Gerardo-Giorda, A. Veneziani, Numerical simulation of a suceptible-exposed-infectious space-continuous model for the spread of rabies in raccoons across a realistic landscape , J. Biol. Dyn. 7 (1) (2013) 31–46. [8] M.Y. Kim, Galerkin Methods for a model of population dynamics with nonlinear diffusion, Numer. Methods Partial Differ. Equ. 12 (1996) 59–73. [9] M.Y. Kim, A numerical method for spatial diffusion in age-structured populations, Numer. Methods Partial Differ. Equ. 12 (2) (2010) 253–273. [10] S. Salsa, Partial Differential Equations in Action, from Modeling to Theory, Springer, 2009. [11] M.A. Khan, A. Atangana, Modeling the dynamics of novel coronavirus (2019-nCov) with fractional derivative, Alex. Eng. J. (2020). [12] Z. Du, X. Xu, Y. Wu, L. Wang, B.J. Cowling, Lauren A. Meyers, Serial interval of COVID-19 among publicly reported confirmed cases, Emerg. Infect. Diseases 26 (6) (2020). A. Viguerie, G. Lorenzo, F. Auricchio et al. / Applied Mathematics Letters 111 (2021) 106617 9 [13] H. Nishiura, N.M. Linton, A.R. Akhmetzhanov, Serial interval of novel coronavirus (COVID-19) infections, Int. J. Infect. Dis. 93 (2020) 284–286. [14] L. Lan, D. Xu, G. Ye, Postitive RT-PCR test results in patients recovered from COVID-19, JAMA 323 (15) (2020) 1502–1503. [15] Coronavirus in Italy: Updated map and case count, 2020, https://lab24.ilsole24ore.com/coronavirus/en/. (Accessed 26 April 2020). [16] J. Backer, D. Klinkeberg, J. Wallinga, Incubation period of 2019 novel coronavirus (2019-CoV) infections among travelers from wuhan, China, 20-28 January 2020, Euro. Surveill. 25 (5) (2020) 2000062. [17] G. Evensen, The Ensemble Kalman Filter, Springer, 2007. [18] Luca Formaggia, Fabio Nobile, Alfio Quarteroni, Alessandro Veneziani, Multiscale modelling of the circulatory system: A preliminary analysis, Comput. Vis. Sci. 2 (2–3) (1999) 75–83. [19] Irene E. Vignon-Clementel, C. Alberto Figueroa, Kenneth E. J