NO2 pollution over selected cities in the Po valley in 2018-2021 and its possible effects on boosting COVID-19 deaths.

This work analyzes nitrogen dioxide (NO2) pollution over a set of cities in the Po Valley in northern Italy, using satellite and in situ observations. The cities include Milan, Bergamo, and Brescia, the first area of the COVID-19 outbreak and diffusion in Italy, with a higher mortality rate than in other parts of Italy and Europe. The analysis was performed for three years, from May 2018 to April 2021, including the period of first-wave diffusion of COVID-19 over the Po Valley, that is, January 2020-April 2020. The study aimed at giving a more general picture of the NO2 temporal and spatial variation, possibly due to the lockdown adopted for the pandemic crisis containment and other factors, such as the meteorological conditions and the seasonal cycle. We have mainly investigated two effects: first, the correlation of NO2 pollution with atmospheric parameters such as air and dew point temperature, and second the possible correlation between air quality and COVID-19 deaths, which could explain the high mortality rate. We have found a good relationship between air quality and temperature. In light of this relationship, we can conclude that the air quality improvement in March 2020 was primarily because of the lockdown adopted to prevent and limit virus diffusion. We also report a good correlation between NO2 pollution and COVID-19 deaths, which is not seen when considering a reference city in the South of Italy. The critical factor in explaining the difference is the persistence of air pollution in the Po Valley in wintertime. We found that NO2 pollution shows a seasonal cycle, yielding a non-causal correlation with the COVID-19 deaths. However, causality comes in once we read the correlation in the context of current and recent epidemiological evidence and leads us to conclude that air pollution may have acted as a significant risk factor in boosting COVID-19 fatalities.

Introduction

The Corona Virus Disease 2019 (COVID-19) is a severe acute respiratory syndrome coronavirus 2(SARS-CoV-2), which was first detected in Wuhan, China, in December 2019, and successively qualified as a global health emergency on Mar 11, 2020, by the World Health Organization (WHO), given the fatality rate and the rapid diffusion outside China (Ghebreyesus, 2020). With its first COVID-19 disease case on Feb 21, 2020, Italy was, after China, the country most affected by the pandemic, counting a large number of deaths (Migliaccio et al., 2021). Therefore, to contain the spread of the virus, the Italian Government imposed an extended period (Mar 9, 2020–May 3, 2020) of restriction on human activities, indicated with lockdown, characterized by curfews, quarantines, travel restrictions, smart working, closure of restaurants and non-essential shops and activities.

Following these measures, the relationship between air pollution mitigation and COVID-19 has attracted the attention of different researchers. In particular, the focus was on the concentration of nitrogen dioxide (NO2), which, together with nitrogen oxide (NO), is one of the major pollutants in Earth's atmosphere, present in both the stratosphere and troposphere, emitted because of anthropogenic processes and activities, such as vehicular traffic, biomass, and fuel combustion (van Geffen et al., 2019). Furthermore, the tropospheric NO2 can severely affect human health, causing, among others, severe respiratory diseases (World Health Organization, 2003; Thurston et al., 2017; Khomenko et al., 2021). Nitrogen dioxide is also a precursor of PM2.5, PM10 (e.g.see https://cfpub.epa.gov/roe/indicator.cfm?i=19). Nitrogen dioxide, among many other environmental gases, contributes to the secondary de-novo formation of PM2.5 and PM10 in the atmosphere. Some studies show how a high concentration of NO2 is associated with a relevant load of aerosol (e.g. (Veefkind et al., 2011),). In this respect, NO2 concentration can be considered a good air quality indicator.

The most recent literature has investigated different impacts of air pollution on the COVID-19 spread both at the global (e.g., Sannigrahi et al., 2021; Sarmadi et al., 2021; Cooper et al., 2022) and regional (e.g., Biswal et al., 2020; Dantas et al., 2020; Zangari et al., 2020; Lonati et al., 2021) scales. More specifically (Zhu et al., 2020), investigated the associations of six air pollutants (PM2.5, PM10, SO2, CO, NO2, and O3) with COVID-19 confirmed cases in China by applying a generalized additive model. They observed significantly positive associations, the same as found by (Fattorini and Regoli, 2020, Magazzino et al., 2020), which provided further evidence that chronic exposure to atmospheric contamination (NO2, O3, PM2.5, and PM10) represents a favorable context for the virus spread. Similar analyses have also been made in (Zoran et al., 2020) that assessed the relationship between levels of NO2 and O3, measured by the ground station, and COVID-19 infections in Milan. In particular, for the January–April 2020 period, time series of daily average inhalable gaseous pollutants, such as O3 and NO2, were analyzed together with climate variables. The results show a positive correlation of daily averaged O3 with air temperature and inversely correlations with relative humidity and precipitation rates.

Other studies have dealt with the link between air pollution or NO2 concentration and COVID-19-related deaths (e.g., Comunian et al., 2020; Filippini et al., 2021). For example, in (Ogen, 2020), the relationship between long-term exposure to NO2 and coronavirus fatality has been analyzed, considering the number of death cases taken from 66 administrative regions in Italy, Spain, France, and Germany. Results show that 78% of the fatality cases were in areas with the highest NO concentrations combined with downwards airflow, preventing efficient air pollution dispersion. In (Travaglio et al., 2021), it was demonstrated that a slight increase in air pollution leads to a significant increase in mortality from COVID-19 in Great Britain. In (Conticini et al., 2020) and (Sasidharan et al., 2020), the correlation between air pollution and virus lethality has been explored in Northern Italy and London, respectively. In (Pansini and Fornacca, 2020) also, the authors found a significant correlation between levels of air pollution measured by satellite and COVID-19 mortality and spread. In particular, they observed more viral infections in the areas affected by high NO2 and PM 2.5 values. Similar results were also found in (Coker et al., 2020), who analyzed the relationship between COVID-19 mortality and PM2.5 concentration.

On the other hand, different studies are based on the lockdown effects on pollution. For example (Donzelli et al., 2020), is a study based on the Italian spring lockdown impact on emissions in three cities in Tuscany that discovered a significant decrease in NO2 concentrations. In (Sasidharan et al., 2020), it was shown that the severe limitation of the lockdown period determined a substantial reduction in pollutants concentration mainly due to vehicular traffic (PM10, PM2.5, BC, benzene, CO, and NOx). In (Wang and Li, 2021) also, the authors investigate the nonlinear impact of COVID-19 lockdown on four pollutants, such as NO2, PM10, O3, and SO2, in eight selected cities using the Spearman correlation function model. They noticed that only NO2 and particles have decreased due to the lockdown, but not O3.

The paper wants to present more general aspects of the NO2 temporal variation and its eventual connection with the COVID-19 pandemic. For this reason, we will use TROPOMI (Copernicus Sentinel 5 Precursor Tropospheric Control Systems) retrievals of NO2 tropospheric column amounts.

Unlike previous studies, our analysis covers an extensive period of 36 months, ranging from May 2018 to April 2021, which also comprehends the period January 2020–April 2020 of COVID-19 diffusion over cities in the Po Valley hit by COVID-19, beginning March 2020. The cities include Milan, Bergamo, and Brescia, the first area of the COVID-19 outbreak and diffusion in Italy, with a higher mortality rate than in other parts of Italy and Europe.

A re-gridding technique extensively discussed and explained in (Cersosimo et al., 2020) was applied to the TROPOMI NO2 observations to have monthly NO2 data with a regular spatial resolution of 1 km. The higher spatial sampling NO2 dataset allowed us to compare the punctual surface observations better. In this respect, the re-sampling tool allowed us to perform an assessment of the correlation between Level 3 TROPOMI NO2 and NO2 concentrations observed measured at nine air quality monitoring stations.

Furthermore, the monthly time series of meteorological variables, such as temperature and relative humidity, always measured by the same ground stations, were investigated to study possible relationships with the NO2 time series and identify a seasonal cycle followed by the different variables. In this way, it was possible to consider a more complex scenario concerning the variation of pollution from NO2, not only concentrated in the lockdown period adopted to the pandemic crisis containment but also linked to other factors, such as the meteorological conditions and the seasonal cycles. In addition, a further aspect analyzed in this paper was the possible link between COVID-19-related deaths and the distribution of tropospheric NO2 pollution, which can severely affect the human respiratory system and thus make people more susceptible to developing severe symptoms if the virus is contracted.

The paper is organized as follows: Section 2 describes the data and the re-gridding methodology applied to TROPOMI NO2 observations. The re-sampling of Level 2 TROPOMI NO2 retrievals and their comparison with ground-based observations are shown in Section 3; the section is mainly focused on the effects of the epidemic crisis on air pollution and the link between NO2 pollution and COVID-19-related deaths. Conclusions are drawn in Section 4.

Data and methodology

Data collection

Different freely available data sources are used in this work, such as remote sensing observations, ground-station measurements, and Istituto Nazionale di Statistica (ISTAT) data.

The Copernicus Sentinel 5 Precursor Tropospheric Monitoring Instrument (S5P/TROPOMI) provided the satellite observations to get the NO2 tropospheric column Level 2 products. These data descriptions have already been detailed (Cersosimo et al., 2020). Therefore, here we limit ourselves to summarizing the essential aspects.

TROPOMI is a single payload onboard the Copernicus S5P, a near-polar sun-synchronous orbit satellite launched on Oct 13, 2017. The satellite flies at an altitude of 817 km, with an ascending node equatorial crossing at 13:30 h, Mean Local Solar time; the repeat cycle is 17 days. The instrument works in a push-broom configuration; the swath across-track is about 2600 km, and that along-track of 7 km improved to 5.6 km from Aug 6, 2019. The instrument achieves global coverage in a day. The instrument is a four-spectrometer system, each electronically split into two bands (2 in UV, 2 in VIS, 2 in NIR, and 2 in SWIR) (Babic et al., 2019; Griffin et al., 2019). The satellite footprint at the ground is 7 km × 3.5 km.

In this work, the Level 2 NO2 satellite data, with a qa_valuve, a quality assurance value greater than 0.75, have been used for May 2018–April 2021. According to (Griffin et al., 2019), qa_valuve have been used as a pixel filter to remove cloudy observations over regions covered by snow/ice, errors, and unreliable retrievals.

The Po Valley is the target area considered in this study. The site is located in the North of Italy; for the present analysis, it is the region extending in longitude from 7° E to 13° E and in latitude from 44° N to 46° N, as shown in Figure 1.

Figure 1

Map of the Po Valley (yellow box) analyzed in this study. Ground stations whose observations are used for validation are shown with blue symbols.

Level 2 NO2 TROPOMI data were re-gridded at a regular spatial grid mesh of step 1 km to obtain Level 3 NO2 total tropospheric column. Then, the data were accumulated and averaged monthly to filter out random fluctuations and, eventually, to fill gaps because of clouds. The methodology used to obtain the monthly maps of NO2 is extensively described in (Cersosimo et al., 2020) and is here briefly outlined for the reader's benefit.

(Cersosimo et al., 2020) has designed and implemented a re-sampling algorithm, which uses the ordinary kriging method (Cressie 2015) to go from Level 2 satellite data at a spatial resolution of 7 km × 3.5 km to Level 3 data on a regular grid with a mesh of 1 km. Let us assume Y(s) be a spatial field whose single values are given a set of spatial points, say, , where sj is the spatial coordinate, e.g., in the case of the 2-D field, the (longitude, latitude) of every single data point. Let us assume we need to estimate Y(s) at a location where there are no data points. According to (Cressie 2015), the estimated value of is given by

In Eq. (1), is a weight-vector whose elements are real numbers and is the field variable at spatial coordinates, s, where it is supposed to be known. We consider the additional condition of weights normalized to unity, i.e.,to get an unbiased estimate. The weights in Eq. (2) are computed, which minimizes the variance of the estimate using a procedure that exploits the correlation properties of the process . Again, according to (Cressie 2015), these properties are described through the so-called semivariogram,where stands for variance. For the analysis at hand, Eq. (3) is the variance of the variable at two different points of coordinates, respectively. Exploiting the given data, we can compute an empirical semivariogram, say Towards this objective, we consider a fitting model that extrapolates the observed points' spatial behavior to the area of interest. We could use several theoretical semivariogram models with known analytical properties and the physical meaning of parameters for our objective. Our methodology uses a set of models: linear, exponential, Gaussian, wave, circular and specular semivariograms (Cersosimo et al., 2020). We use a Least Squares criterion procedure to select which best fits the given data.

Figure 2 exemplifies the monthly spatial maps of Level 3 TROPOMI NO2 over the Po Valley area for a winter month (December 2020).

Figure 2

TROPOMI nitrogen dioxide tropospheric column (Level 3 data in units of Dobson) over the Po Valley for December 2020.

In situ measurements of nitrogen dioxide were obtained on a network of 9 monitoring stations distributed in the target region (blue symbols in Figure 1). The stations belong to Agenzia Regionale per la Protezione dell'Ambiente (ARPA) of Lombardia Region (ARPA Lombardia, 2022). Ground-based observations also included a 2-m air temperature and relative humidity, available at a time slot of 1 h. Therefore, hourly values have been extracted. The daily/monthly averages have been calculated, as well, for comparison with coherent satellite observations. In addition, the proper functioning of the air quality equipment is checked regularly and inspected over time to ensure that the data meet the quality standards.

For the Po Valley orography favors, in wintertime, air subsidence, with a poor atmospheric circulation, which prevents pollution from diffusion, hence dispersion, in the free troposphere (Cersosimo et al., 2020). To provide proper comparisons with an area that, although polluted, is not subject to air subsidence in winter, we have also used the data from the city of Naples in the South of Italy. In contrast, to the Po Valley, Naples is influenced by a marine environment, where the air is continuously moved because of atmospheric processes in the boundary layer.

Finally, to analyze the possible link between COVID-19 related deaths and the distribution of tropospheric NO2 polluting, it was necessary to procure the number of deaths in correspondence with the cities shown in Figure 1. The data consists of the monthly death number provided by ISTAT (ISTAT, 2022).

As shown in Figure 2, the cities considered in the present analysis are the bulk of the NO2 pollution. In particular, they include the cities of Bergamo and Brescia, where the highest rate of COVID-19 mortality was recorded during the first wave of virus diffusion in the winter of 2020 (Casti and Adobati, 2020a; 2020b).

Results

NO2 TROPOMI data vs. in situ measurements

To begin with, we will first compare ground-based measurements vs. Level 3 TROPOMI NO2 data to show the capability of TROPOMI to sense down to the surface and discriminate between low and high pollution loads. Next, the data were collocated adequately in space-time to compare satellite and in situ observations. In situ observations have been compared with the spatially closest level 3 satellite data. For temporal colocation, each TROPOMI observation has the satellite overpass time; therefore, we used the temporally closest in situ measurement.

The Pearson correlation coefficient, or simply correlation coefficient, , between Level 3 NO2 TROPOMI and NO2 concentration at every single station is summarized in Table 1. The Pearson correlation coefficient is defined according towhere cov and var stand for covariance and variance, respectively, and are two generic vector-valued parameters. In practice, and its 95% confidence interval is computed with the Matlab script corrcoef version 2020b.

Table 1

Correlation coefficient (R2) between Level 3 NO2 TROPOMI and NO2 concentration at the stations shown in Figure 1. For comparison, we also show the correlation between Naples (yellow), a marine city in South Italy. The last column is the Confidence Interval at the level of 95% of estimated.

Ground station	Longitude	Latitude	R2	95% CI
Bergamo	9.64	45.69	0.86	[0.74, 0.93]
Brescia	10.22	45.54	0.81	[0.66, 0.90]
Como	9.08	45.80	0.80	[0.64, 0.89]
Cremona	10.02	45.13	0.87	[0.76, 0.93]
Lecco	9.39	45.85	0.71	[0.50, 0.74]
Lodi	9.50	45.31	0.82	[0.67, 0.90]
Mantova	10.79	45.16	0.82	[0.67, 0.90]
Milan	9.20	45.46	0.82	[0.67, 0.90]
Pavia	9.15	45.18	0.86	[0.74, 0.93]
Naples	14.32	40.89	0.58	[0.31, 0.76]

The linear correlation coefficient, , computed according to Eq. (4), and its 95% confidence interval are both shown in Table 1, which covers the period May 2018 and April 2021 and considers monthly averages; therefore, the number of couples () is , where indicate the monthly average of TROPOMI level 3 data (in Dobson) and in situ measurements (in units of μg/m3), respectively. We also stress that the correlation refers to the single station whose spatial coordinates are shown in Table 1.

The analysis shown in Table 1 allows us to conclude that Level 3 NO2 TROPOMI data are representative of the air quality in the lower troposphere because they are well correlated with the surface, in situ, measurements of NO2 concentration. Therefore, satellite data can give better insights into understanding the evolution of spatial patterns. It is also interesting to note that the correlation falls to 0.58 when we consider Naples, a marine environment characterized by a crucial diurnal cycle in the boundary layer, which helps the diffusion of pollution in the free atmosphere. These processes introduce spatial and temporal variability at scales that are not entirely captured from the satellite (Cersosimo et al., 2020).

The Level 3 NO2 TROPOMI data at the Po Valley stations in Table 1 have been merged and correlated with the monthly air temperature, and dew point temperature, . The results are shown in Figure 3, and it is seen that the correlation with the dew point temperature is not better than that observed with air temperature. In other words, including information coming from humidity does not improve the correlation. For this reason, hereafter in, this paper will consider air temperature alone to represent the meteorological conditions.

Figure 3

Scatterplot of TROPOMI NO2 vs. ground-based measurements of temperature (a) and dew point temperature (b). The data correspond to May 2018–April 2021. All the data from the 9 Po Valley stations shown in Table 1 have been merged.

The good correlation with air temperature is further evidenced when considering the three-year time series of NO2 and temperature, as shown in Figure 4. Again, it is seen that the NO2 and temperature have an evident seasonal cycle, with peaks of NO2 pollution correlated with troughs of temperature.

Figure 4

Time series (top to bottom) of for the station of Bergamo showing the trend-seasonal model fitted to the data.

The example shown in Figure 4 corresponds to the Bergamo station; however, similar results are obtained at the other stations analyzed in this paper. The anti-correlation NO2 vs. temperature is mainly linked to colder temperatures, typical of the winter months, leading to more massive domestic and commercial heating and more extensive use of private cars because of the winter school season.

To have a quantitative assessment of the amplitude and phase of the seasonal cycle, we have fitted to the time series the mathematical model,which is made up of a linear trend and a seasonal component of 12 months. In Eq. (5) is the generic observation () at time (in units of months), is the basic period in the series ( =12 months), and the coefficients a, b, c, d are fit parameters computed, e.g., through the Least Squares criterion. Figure 4 also compares the model and data points. We see that most of the variability of the time series is expressed by the seasonal cycle, as it is summarized in Table 2 to Table 4, where we show the fit parameters a, b, c, d and the ratio of the variance of the cycle to the total variance of the series. It is also seen that the trend coefficient, is zero within the error bars for most of the cases. Notably, the NO2 concentration does not show an evident seasonal cycle for the Naples station and is characterized by a decreasing trend. However, focusing on the in situ NO2 concentration, the bulk pollution (the coefficient in Eq. (5)) for Naples is even more extensive than in Milan. Naples has levels of air pollution comparable with those of the Po Valley. The different climate and meteorological conditions make the difference. The Naples marine environment and its orography generate a turbulent boundary layer capable of mixing air and dispersing pollution in the free troposphere (e.g., (Cersosimo et al., 2020)). The good correlation we observe for the Po Valley among the various parameters, exemplify the presence of a poor dynamic in the boundary layer and the presence of characteristic time scales much larger than the diurnal cycle. TROPOMI is a polar satellite; hence, it cannot resolve the daily cycle and has 3–4 km spatial scales. Therefore, the fact that the parameters, are well correlated for the Po Valley confirms that air pollution has a high space-time persistence.

Table 2

Fit parameters a, b, c, d, correlation coefficient (R2), and variance ratio of the cycle to the total variance of the Level 3 TROPOMI NO2 series for all reference stations.

Level 3 TROPOMI NO2
Station	Fit coefficients				R^2	Variance Ratio
	a (Dobson) [95% CI]	b (Dobson/month) [95% CI]	c (Dobson)c [95% CI]	d (radians) [95% CI]
Bergamo	0.33 [0.26, 0.39]	-0.00055 [-0.00372, 0.00261]	-0.23 (-0.27, -0.18)	-1.28 [-1.48, -1.08]	0.76	0.75
Brescia	0.33 [0.28, 0.37]	-0.00169 [-0.00375, 0.00037]	0.21 [0.18, 0.24]	1.78 [1.64, 1.92]	0.87	0.88
Como	0.27 [0.21, 0.33]	-0.0012 [-0.0041, 0.0017]	0.17 [0.13, 0.21]	1.75 (1.50, 1.99]	0.68	0.69
Cremona	0.28 [0.24, 0.33]	-0.0021 [-0.004, 0.00004]	0.20 [0.17, 0.23]	-4.46 (-4.61, -4.3]	0.84	0.86
Lecco	0.21 [0.16, 0.26]	0.00018 [-0.00240, 0.00274]	-0.15 (-0.19, -0.11]	4.95 (4.71, 5.20]	0.69	0.67
Lodi	0.30 [0.24, 0.37]	-0.0013 [-0.0043, 0.0017]	0.23 [0.19, 0.28]	-10.72 (-10.90, -10.53]	0.78	0.79
Mantova	0.27 [0.23, 0.31]	-0.0023 [-0.0043, -0.0003]	0.18 [0.15, 0.21]	1.72 [1.56, 1.88]	0.83	0.85
Milano	0.41 [0.34, 0.48]	-0.0016 [-0.0050, 0.0016]	0.28 [0.23, 0.33]	1.84 [1.67, 2.01]	0.82	0.82
Pavia	0.29 [0.23, 0.35	-0.0017 [-0.0045, 0.0009]	0.21 [0.19, 0.25]	-4.47 [-4.67, -4.28]	0.78	0.79
Naples	0.17 [0.15, 0.19]	-0.0037 [-0.0048, -0.0027]	-0.031 [-0.046, -0.016]	-0.5388 [-1.05, -0.032]	0.68	0.17

Table 4

Fit parameters a, b, c, d, correlation coefficient (R2), and cycle variance-ratio to the Air Temperature series total variance for all reference stations.

Air Temperature
Station	Fit coefficients				R^2	Variance Ratio
	a (°C) [95% CI]	b (°C/month) [95% CI]	c (°C) [95% CI]	d (radians) [95% CI]
Bergamo	14.19 (12.01, 16.37)	-0.010 (-0.11, 0.09)	10.40 (8.89, 11.91)	-7.59 (-7.73, -7.44)	0.87	0.84
Brescia	16.08 (13.77, 18.39)	0.004 (-0.106, 0.114)	11.73 (10.13, 13.33)	-7.59 (-7.72, -7.45)	0.88	0.86
Como	14.97 (12.81, 17.13)	-0.003 (-0.105, 0.099)	10.37 (8.88, 11.87)	-13.87 (-14.01, -13.73)	0.87	0.84
Cremona	15.13 (13.03, 17.22)	-0.0018 (-0.1016, 0.0978)	-12.01 (-13.46, -10.56)	-4.44 (-4.56, -4.32)	0.90	0.88
Lecco	14.86 (12.79, 16.92)	0.003 (-0.095, 0.102)	11.57 (10.14, 13)	-13.86 (-13.98, -13.73)	0.90	0.88
Lodi	14.00 (12.01, 15.99)	0.0099 (-0.0846, 0.1046)	12.05 (10.68, 13.43)	-13.87 (-13.98, -13.75)	0.91	0.89
Mantova	14.86 (12.79, 16.92)	0.0033 (-0.0950, 0.1016)	-11.57 (-13, -10.14)	-4.43 (-4.55, -4.31)	0.90	0.88
Milano	14.6 (12.52, 16.69)	0.04 (-0.05, 0.14)	11.55 (10.11, 13)	-7.60 (-7.73, -7.48)	0.89	0.89
Pavia	14.86 (12.79, 16.92)	0.0032 (-0.0950, 0.1016)	11.57 (10.14, 13)	-7.57 (-7.69, -7.45)	0.90	0.88
Naples	17.77 [16.61, 18.93]	0.041 [-0.014, 0.096]	8.09 [7.28, 8.90]	-8.12 [-8.21, -8.02]	0.93	0.93

TROPOMI spatial patterns of NO2 pollution

Figure 5 shows the monthly maps of Level 3 TROPOMI NO2 over the Po Valley area for the winter months, December to March, and from 2018 to 2021.

Figure 5

NO2 TROPOMI monthly maps (in units of Dobson) over the Po Valley area between December 2018 and March 2021. Same months appear on the same line.

These maps help us better understand the evolution of nitrogen dioxide pollution over the Po Valley when COVID-19 spread, leading Italian authorities to declare a lockdown period. To allow for the comparison of no COVID-19 with COVID-19 time, we show the monthly maps of tropospheric nitrogen dioxide over the target area from December to March for 2018–2019 (no COVID-19), 2019–2020, and 2020–2021. The first feature that arises from comparing the corresponding maps for the 2018–19 (no COVID-19) and 2019–20 (the first wave COVID-19) winter period is the sharp decrease in air pollution in 2020 compared to 2019. However, successively, we can see a new increase if we compare the maps corresponding to the 2020–21 winter period (the second wave of COVID-19) to that of 2019–20.

The relatively lower pollution in winter 2019–2020 (December to March) was most likely due to the warmer 2019–20 winter, limiting heating power use for residential and commercial sectors (e.g., see https://climate.copernicus.eu/boreal-winter-season-1920-was-far-warmest-winter-season-ever-recorded-europe-0). From Figure 5, which refers to the different winters, it is seen that only January 2019 and 2020 exhibits comparable NO2 pollution. In contrast, February 2020 shows a dramatic decrease in NO2 pollution over Po Valley, and it has been the month of maximum virus spread in the Po valley. In March 2020, when the lockdown was active, we saw a 50% decrease in NO2 pollution compared to the same month of 2019.

To understand how atmospheric conditions drive the patterns shown in Figure 5, Figure 6 exemplifies the monthly mean temperature for the city of Milan in the three winters of interest. We see that for winter 2019–2020 (December to March), the air temperature has been higher or comparable than in 2018–2019 and 2020–2021, but for March, we see an average temperature lower than 5–6 °C in 2020. Because of the good correlation between NO2 pollution and air temperature, the comparison for March shows an anomaly because March 2020 has been harder than the two in 2019 and 2021. In effect, in March 2020, the Italian Government issued a rigid and restrictive lockdown, which stopped many industrial activities and limited the circulation of people.

Figure 6

Air temperature for the Milan station and the three winters 2018–2019, 2019–2020, 2020–2021.

Because of the good correlation NO2 vs. air temperature, we expected March 2020 to be more polluted than March 2019 and 2021 or at least to show comparable pollution. We see the reverse means that the lockdown effectively reduced pollution. Also, the decrease in pollution is likely because of the stop to vehicular traffic. Conversely, domestic heating increased because people were confined at home, and March 2020 has been colder. The lockdown experiment shows that vehicular traffic is the primary source of NO2 pollution in the PO Valley. Figure 6 refers to Milan; however, the same results as those shown in Figure 6 have been found for the other Po Valley stations in Table 1. They are not shown for the sake of brevity.

Effect of air pollution on deaths in the COVID-19 period

According to (Casti and Adobati, 2020a; 2020b), three factors have mainly explained the Po valley fragility and, particularly, the Lombard Region. These are.

Persistence of pollution

Health care system

Rhizomatic commuting

We see that one of the main factors is environmental pollution, or rather its persistence. As already mentioned, the particular local weather conditions determine the pollution persistence in the Po Valley. Wintertime weather conditions are primarily governed by air subsidence, an effect of the Po valley orography, which prevents dispersal. Rhizomatic commuting refers to the intricate net of the transportation system and the many people who use it to move around the area and get to work. The welfare system in Northern Italy is considered the best in Italy, although it has transformed from public to private in recent years. The risk factor is linked to the cultural-social environment that isolates the elderly in private nursing homes more than in the hospital system. In fact, with public hospitals closed due to Covid-19, such facilities have become hotbeds of infection during the pandemic. However, risk factors b) and c) while justifying a greater spread of disease do not explain the higher incidence of deaths from COVID-19 compared to other regions. For example, the mortality rate in Lombardy is 345 per 100,000 inhabitants, while, e.g., Campania, which has a less efficient welfare system (e.g., see https://careonline.it/2018/09/la-misura-della-performance-dei-servizi-sanitari-regionali/) than the Lombard one, it is only 145, or less than half (e.g., see https://lab24.ilsole24ore.com/coronavirus/). According to Osservatorio Nazionale sulla Salute nelle Regioni Italiane (osservatoriosullasalute.it) Covid-19 does not cause the same mortality everywhere but manifests itself with extreme variability in the Italian Regions, ranging from a maximum of 5.4% of positives in Lombardy to a minimum of 1.3% in Campania, with an average of 3.5% at the national level. Focusing on the period October–December (2020), in particular on the data from Oct 12 to Dec 6 (2020), it should be noted that the mortality levels for Covid-19 in the Italian regions vary significantly, with the same prevalence of new infections and regardless of the age structure of the resident population.

While risk factors s b) and c) can have increased the diffusion of the virus, air pollution is a direct cause of mortality (e.g. (Thurston et al., 2017),), and according to current statistics (Khomenko et al., 2021), the cities in Europe with the most significant risk of premature deaths from pollution are Bergamo and Brescia while Milan is still in eleventh place.

With this in mind, we analyzed for 2020 year the number of fatalities affecting the ten municipalities listed in Table 1, and we have correlated it with the corresponding in situ NO2 concentration. We obtained the monthly number of deaths per municipality thanks to the data made available by ISTAT. From the same source, we also have the monthly ends for cities averaged over the period 2015–2019, before the COVID-19 event. These background records can filter out the number of total deaths in 2020 to have the excess mortality because of COVID-19.

The background deaths say Dbg, are shown for Milan and Naples for comparison in Figure 7. They have a clear seasonal cycle, with more deaths occurring in the wintertime as expected. In effect, Dbg is well correlated with the climatology air temperature, which shows the same seasonality. For both stations was computed according to the model of Eq. (5) corresponding to Milan and Naples, respectively. In other words, for January to December was obtained by considering as defined in Eq. (5), corresponding to the air temperature variable and with ranging from January to December. To consider only the cycle, we set . For Milan city, the correlation is -0.67, whereas for Naples is -0.70. However, if we consider a delay between air temperature and deaths of one month, that is, we correlate (t) and , with month, the correlation rise to -0.79 and -0.82 for Milan and Naples, respectively. A delay between air temperature and death is expected because people first get hospitalized and eventually die because of the disease.

Figure 7

Background mortality (average 2015–2019) for the cities of Milan (Po Valley) and Naples (South of Italy)'.

The negative correlation says that the lower the temperature, the higher the ends, which is well understood because lower temperatures correspond to wintertime when, e.g., the elderly population is exposed to flu viruses and pathologies linked to respiratory diseases. Outside COVID-19 times, there is no important difference between Milan in the Po Valley and Naples in the South regarding the annual cycle of mortality. However, based on the population of Naples (∼1 M inhabitants) and that of Milan (∼1.35 inhabitants), the yearly mortality rate of the two cities is 10.52 ends per 1000 inhabitants for Naples and 12.09 for Milan. The larger rate for Milan is thought to be because of bad air quality and related mortality increases (Khomenko et al., 2021; Thurston et al., 2017).

The effect of air quality can be analyzed by considering the correlation of with the seasonal cycle, say , of NO2 pollution, estimated according to Eq. (1), the same we did for . For Milan, we have 0.74, whereas for Napoli 0.05. Once again, for Milan, if we consider a delay between air pollution and deaths of one month, that is, we correlate (t) and , with month, the correlation rises to 0.76. The analysis shows us to conclude that for Milan, the persistence of air pollution is a risk factor independent of the COVID-19 pandemic. In effect, Naples is polluted the same as Milan, if even worse, as evidenced in section 3.1 by comparing the background coefficient, a in Table 3 for Milan and Naples. Yet, for Naples, we have 0.05, as far as the correlation between air pollution and is concerned. The result for Naples shows that the correlation is not universal but depends on local effects, which for the case at hand, is likely driven by weather conditions, which favors the persistence of poor air quality.

Table 3

Fit parameters a, b, c, d, correlation coefficient (R2), and variance ratio of the cycle to the total variance of the in situ NO2 concentration series for all reference stations.

In situ NO2 concentration
Station	Fit coefficients				R^2	Variance Ratio
	a (μg/m3)[95% CI]	b (μg/m3-month) [95% CI]	c (μg/m3) [95% CI]	d (radians) [95% CI]
Bergamo	29.75 [25.52, 33.98]	-0.29 [-0.49, -0.09]	14.81 [11.87, 17.75]	7.99 [7.80, 8.18]	0.77	0.79
Brescia	34.39 [30.19, 38.59]	-0.12 [-0.32, 0.08]	-9.03 [-11.93, -6.12]	-1.24 [-1.56, -0.93]	0.56	0.57
Como	41.85 [38.65, 45.04]	-0.31 [-0.45, -0.15]	8.12 [5.90, 10.34]	1.67 [1.41, 1.94]	0.66	0.62
Cremona	25.58 [21.71, 29.46]	-0.13 [-0.31, 0.06]	14.37 [11.68, 17.07]	7.93 [7.75, 8.12]	0.79	0.81
Lecco	18.23 [14.26, 22.19]	0.23 [0.038, 0.41]	10.5 [7.75, 13.24]	8.10 [7.84, 8.34]	0.71	0.56
Lodi	31.46 [28.16, 34.76]	-0.21 [-0.37, -0.06]	11.72 [9.44, 14.01]	-35.92 [-36.11, -35.73]	0.77	0.79
Mantova	20.19 [16.58, 23.79]	0.09 [-0.08, 0.26]	10.71 [8.21, 13.20]	1.85 [1.62, 2.08]	0.73	0.66
Milano	50.08 [42.9, 57.26]	-0.40 [-0.74, -0.06]	19.55 [14.56, 24.54]	7.97 [7.72, 8.21]	0.67	0.69
Pavia	28.93 [25.92, 31.94]	-0.18 [-0.32, -0.03]	-13.90 [-15.99, -11.81]	11.13 [10.99, 11.28]	0.85	0.88
Naples	61.71 [53.21, 70.22]	-0.5784 [-0.98, -0.17]	2.27 [-3.59, 8.12]	-2.54 [-5.10, 0.03]	0.25	0.015

With this in mind, we try to establish if there is any correlation between air pollution and excess mortality, defined according to where represents the actual total deaths. The analysis is performed by using the in situ observations and because they are taken simultaneously and have the same time-spatial scales. As before, we consider monthly averages alone. The analysis will be performed by distinguishing noCOVID-19 from COVID-19 times.

To begin with, we analyze the correlation in 2019, that is, before the pandemic crisis. Because we expect a delay between mortality and the seasonal cycle of and we analyzed the correlation as a function of the delay, . The results are shown in Figure 8 for the two stations of Milan and Naples. We have explored the range, = [0,3] months, because larger delays are physically unrealistic and also because the degrees of freedom scale linearly according to , with months.

Figure 8

Milan and Naples station. The correlation coefficient between vs and vs as a function of the delay, . The two magenta lines give the 95% Confidence Interval for zero correlation; therefore, only outside the CI range are statistically and significantly different from zero.

Figure 8 shows that in noCOVID-19 time, there is no statistically significant correlation between vs and vs . For Naples, for the case vs , we observe 0.54 with = 0. This result represents a borderline case because it falls close to the 95% CI upper bound for zero correlation. In conclusion, the low values of observed for both Milan and Naples lead us to conclude that there is no evidence of statistically significant correlation once the deaths have been deseasonalized according to the background .

Conversely, when we consider the same correlation analysis for 2020 (see Figure 9), we have for Milan the values 0.85, if we consider a delay, , in the case of vs . Therefore, the excess of mortality , concerning the background, is correlated with air pollution. For temperature, we observe a peak 0.58, which is not statistically significant at a confidence level of 95%. For Naples, the correlation is not statistically significant for the case vs . However, although borderline, we observe a positive, counterintuitive correlation with air temperature. We have 0.70, with months. The second wave of COVID-19 occurred in the autumn of 2020, following the summer holidays when most of the lockdown measures were removed.

Figure 9

As Figure 8, but for 2020.

The result for Milan yields a significant correlation between vs . The finding is also substantiated by the analysis of the remaining eight stations in the Po Valley. The results are shown in Table 5 for the case months, and we see that the correlation is significantly diverse from zero for all Po Valley stations, but Lecco is borderline. To complete the analysis, Table 5 also shows the result for Naples, which offers a not statistically significant correlation (note that for Naples, the CI crosses the zero line, which is not the case for the Po Valley stations, but Lecco as already said).

Table 5

Correlation index between NO2 pollution and as a station function and delay months. The CI is related directly to the observed value of shown in the second column. If CI includes a zero value, the correlation is not statistically significant.

Station	R^2	95% Confidence Interval for R^2
Bergamo	0.92	[0.69, 0.98]
Brescia	0.86	[0.50, 0.97]
Como	0.89	[0.59, 097]
Cremona	0.89	[0.59, 0.97]
Lecco	0.62	[-0.02, 0.90]
Lodi	0.92	[0.69, 0.98]
Mantova	0.84	[0.45, 0.96]
Milan	0.85	[0.47, 0.96]
Pavia	0.83	[0.42, 0.96]
Naples	0.46	[-0.24, 0.85]

It is noteworthy that reaches its maximum for the Bergamo station. In effect, the Bergamo municipality has been that with the highest mortality rate in the Po Valley.

Finally, we stress that a delay between severe air pollution conditions and deaths is expected. COVID-19 is not killing people instantaneously. Instead, people get hospitalized, eventually succumbing to the disease later. Figure 10 shows the normalized series of infection, hospitalization, intensive care, and deaths in Italy from Jul 1, 2020, to Jun 30, 2021. The time lag between detected illness and hospitalization and death is estimated at 17 and 27 days. These have to be understood as the lower bound of the delay because it does not include the delay between the time infection developed and the time it was detected. Also, the data shown in Figure 10 refers to the period between summer 2020 and summer 2021, when the COVID-19 monitoring system was active and working, which was not the case in winter 2020. Therefore, Figure 10 boosts the presence of a delay months we have revealed with our analysis because it is consistent with the present information about the spread, infections, and deaths of COVID-19.

Figure 10

Normalized time series for detected infection (blue), deaths (black), hospitalization (yellow), and intensive care (green) in Italy in the period between Jul 1, 2020, and Jun 30, 2021.

Discussion and conclusions

We have analyzed nitrogen dioxide (NO2) pollution over the Po Valley in Northern Italy. The analysis was performed for three years, from May 2018 to April 2021, covering January 2020–April 2020 of the first COVID-19 diffusion over the Po Valley. The study was twofold: first, we wanted to analyze the effect of the lockdown on air pollution, and second, the possible correlation between pollution and the excess of COVID-19 deaths observed on the first wave of the pandemic spread in the Po Valley compared to the other densely populated part of Italy but with climate and orography which favor a dynamic, turbulent, boundary layer. Therefore, we chose the city of Naples.

The decrease of air pollution during COVID-19 has been variously claimed. However, the concomitance of weather effects has made it very difficult to assess the lockdown impact. Therefore, in this study, we have performed an in-depth analysis using satellite and in situ observations of NO2 spanning from noCOVID-19 to COVID-19 times and air temperature records.

For the Po Valley, we have shown a good correlation between satellite and in situ observations of NO2. The good correlation gives us confidence in using satellite data to analyze spatial patterns. The analysis of these spatial patterns has revealed that the bulk of NO2 pollution extends to Milan, Bergamo, and Brescia, which COVID-19 severely hit in terms of mortality.

Another significant result is the good, negative correlation of NO2 concentration with air temperature. According to ECMWF, winter 2020 has been the warmest on record. However, March 2020 has been colder over the Po Valley than, e.g., those in 2019 and 2021. Yet, the air pollution in March 2020 has been some 50% lower than that observed in March 2019. The only possible explanation is the lockdown issued in March 2020. The lockdown was an experiment that stopped traffic and many industrial activities.

On the other hand, domestic heating was seemingly at its peak because of the colder temperature in March 2020. As a result, the lockdown experiment has shown that NO2 pollution's primary source is traffic and non-electric transportation. This fact supports the current EU policies for the green deal (https://ec.europa.eu/info/strategy/priorities-2019-2024/european-green-deal_it) is a helpful tool to improve air quality and reduce CO2 emissions simultaneously.

Table 6 summarizes the different cases we have examined to establish a correlation between pollution and deaths. Regarding Table 6, it is seen that irrespectively to COVID-19, there is a good background correlation between mortality and air temperature, which share the same seasonal cycle. Furthermore, heavily polluted cities such as Milan correlate mortality and air pollution. However, this correlation seems peculiar to the Po Valley and primarily due to air subsidence during winter. In effect, heavily polluted cities, such as Naples, in a more favorable environment, such as a marine setting, do not show any critical correlation between air quality and mortality.

Table 6

Summary of the higher correlation index as a function of the station and delay for the cases and couples examined in this study, the delay is units of months. Green boxes show the matter for which the correlation is statistically significant; yellow boxes refer to borderline cases (the 95% confidence interval). White boxes refer to values that are not statistically significant.

Analysis of the monthly series of background values
D_bg vs Y_T_a		D_bg vs Y_NO_2
Milan	Naples	Milan	Naples
R^2=-0.79, τ=1	R^2=-0.82, τ=1	R^2=-0.74, τ=1	R^2= 0.05, τ=0
Analysis for the 2019 monthly series (noCOVID)
ΔD vs T_a		ΔD vs q_G
Milan	Naples	Milan	Naples
R^2=-0.40, τ=0	R^2=-0.55, τ=3	R^2= 0.01, τ=2	R^2= 0.54, τ=0
Analysis for the 2020 monthly series (COVID)
ΔD vs T_a		ΔD vs q_G
Milan	Naples	Milan	Naples
R^2=-0.58, τ=1	R^2= 0.70, τ=3	R^2= 0.85, τ=2	R^2= 0.47, τ=2

When we consider noCOVID-19 times, the deseasonalized deaths do not correlate with temperature and NO2 pollution. In other words, the correlation is an effect of the seasonal cycle alone.

COVID-19 has imposed an additional cycle (e.g., see Figure 10), which boosts the correlation, especially with air pollution, which is a new fact. For the Po Valley, the COVID-19 deaths correlate better with air quality than the air temperature. However, we know that correlation does not (necessarily) imply causation. To check if the correlation underlies a mere statistical effect alone, we have also used a reference city outside the Po Valley, Naples, polluted than Milan. The comparison helped us establish that the critical element that introduces correlation is persistence, peculiar to the Po Valley. The persistence introduces a seasonal component to the NO2 pollution, which is not seen in Naples.

COVID-19 diffusion and NO2 pollution are primarily driven by weather, although the sustaining mechanisms are different. NO2 load is sustained by a temperature inversion in the boundary layer, likely to occur in wintertime. The mechanism traps polluting agents close to the surface. On the other hand, COVID-19 diffusion is sustained by indoor living because of adverse weather conditions. The expected weather forcing makes the two events develop a similar seasonal cycle: the larger the COVID-19 diffusion, the larger the NO2 pollution. The correlation we see is just an effect of the seasonal process, and as such, it could be a simple mathematical consequence with no physical causation. However, although independent, the two cycles can interact through an epidemiological mode, which boosts the COVID-19 fatalities. In fact, according to (Khomenko et al., 2021), the relationship between respiratory diseases and air pollution is no more a statistical correlation but a piece of epidemiological evidence. People have been exposed to the most prominent air pollution in the Po Valley during the most extensive diffusion of COVID-19, which can explain why the disease has been more lethal in Milan than in Naples. A mere non-causal effect would not justify the difference in mortality rate because of the pandemic.

Epidemiology could reasonably explain a causal effect between COVID-19 mortality excess and NO2 pollution. Air pollution leaves people more exposed to a respiratory virus, such as COVID-19, because poor air quality is a concomitant agent, increasing the risk of premature deaths (Khomenko et al., 2021; Thurston et al., 2017). It is also noteworthy that the same conclusion does not hold for Naples. Again, a key risk factor seems to be the persistent exposure to a high concentration of pollutants during the most intense virus spreading. This condition is missing for Naples, which shows a poor correlation between NO2 pollution and COVID-19 deaths. It could be argued that air quality is assessed in terms of PM2.5 and PM10, rather than NO2. However, we stress that NO2 is a good indicator of air quality because it is a precursor of PM2.5 and PM10 (Veefkind et al., 2011).

The causal relationship between COVID-19 mortality and air pollution has also been addressed by (Magazzino et al., 2021), who demonstrated, in line with our causal findings, that there exists a unidirectional causal effect from PM2.5 to Deaths, NO2 to Deaths, and economic growth to both PM2.5 and NO2. The analysis has been performed for the New York state and is limited to the first COVID-19 wave. In contrast, our study includes the second wave and elucidates the interrelationship with meteorological conditions.

Our analysis has been shown using the excess of mortality ΔD, which, based on its definition, is the difference between the total deaths and the background. The excess of mortality has also been used by (Coker et al., 2020), who analyzed the correlation between COVID-19 mortality and PM2.5. An analysis for Italy based directly on the COVID-19 deaths or case fatality rate (CFR) has been performed by (Timelli et al., 2021). The excess of mortality as defined in this study has been proposed and used by ISTAT and ISS (the National Department of Health) to address the impact of the COVID-19 pandemic on the total fatality rate in Italy (e.g., https://www.istat.it/it/files//2021/06/Report_ISS_Istat_2021_10_giugno.pdf). In the same report, they show that and CFR exhibit the same time behaviour, which demonstrates that the excess of mortality is dominated by COVID-19 deaths. This is also exemplified with the help of Figure 11, which compares for the city of Milan, the background, and its 95% confidence interval to the number of deaths in 2020. We remember that the background was the average number of total ends in 2015–2019. It is seen that outside the COVID-19 outbreaks, the background and 2020 deaths agree within the error bars. Conversely, during COVID-19 diffusion peaks, the number of deaths tends to double concerning the background and is outside the natural statistical variability of the background itself.

Figure 11

Milan station. The figure compares the background and its 95% confidence interval to the number of deaths in 2020. The 95% Confidence Interval is shown with error bars whose magnitude is two times the monthly sample standard deviation.

In conclusion, we have demonstrated that in Po Valley, the weather conditions develop a seasonal cycle of NO2 load, which correlates with the COVID-19 diffusion. Combining this result with the evidence that poor air quality increases the risk of premature deaths (Khomenko et al., 2021; Thurston et al., 2017), we conclude that NO2 pollution may be considered a significant risk factor. Therefore, it demands actions to mitigate its effects and new approaches to sustainable mobility policies, as suggested, e.g., in (Migliaccio et al., 2021).

Declarations

Author contribution statement

Carmine Serio: Conceived and designed the experiments; Analyzed and interpreted the data; Wrote the paper.

Angela Cersosimo: Contributed reagents, materials, analysis tools or data; Wrote the paper.

Guido Masiello: Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data.

Funding statement

The authors were supported by MIUR, Italy [D.D. 2261 del 6.9.2018, PON R&I 2014–2020 and FSC), project ARS01_00405].

Data availability statement

Data will be made available on request.

Declaration of interests statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.