**1. Introduction**

Air quality monitoring is one of the major challenges that European institutions jointly with national and local administrations are facing in terms of environmental protection. In particular, the 2008 European Air Quality Directive (AQD) 2008/50/EC [1] requires EU Member States to design appropriate air quality plans for zones where the air quality does not comply with the AQD limit values. In the last few decades, European countries implemented various modeling methods to assess the effects of local and regional emission abatement policy options on air quality and human health [2]. On the one side, they include scenario approaches, in which running a chemical-physical simulation model with and without a specific emission source allows for quantifying the impact on air quality levels [3,4]. On the other side, they also include more comprehensive and multidisciplinary approaches, such as Integrated Assessment Models (IAM), which combine simultaneously many features of the economy, society, and scientific findings. These models are based on the combination of multiple mathematical tools and allow for assessing the impact of environmental policies or to improve the air quality control system. Typical tools are the full cost-benefit analyses [5], in which abatement measures, costs, and benefits are expressed in monetary units, optimization, and spatial analysis [6,7].

In areas such as Northern Italy, where the industrial transition in the 1990s reduced coal burning and sulphur concentration, the large majority of environmental studies focus their attention on toxic pollutants that are produced by thermic vehicle engines and house heating plants. These are known to generate serious health effects [8]. Total nitrogen oxides (NOx), nitrogen dioxide (NO2), and particulates matters (PM10 and PM2.5) belong to this class.

According to the above EU rules, governments adopted standards and quantitative limits for pollutant emissions to make economic agents responsible and implement abatement policies. In particular, the maximum concentration for NOx and NO2 is set to 40 μg/m<sup>3</sup> annual average and 200 μg/m<sup>3</sup> hourly not to be exceeded more than 18 times in a single year. Figure 1 represents the average concentration levels of NO2 in Europe for the year 2018. The Po basin in Northern Italy stands out as a heavily polluted area with difficulties in pollution management. The negative impact on society is not limited to health only. There is increasing evidence showing that bad air quality in general, and high NO2 concentrations in particular, impact the economy, including finance [9] and tourism [10].

**Figure 1.** Annual average NO2 concentrations (μg/m3) in Europe during 2018. Levels are expressed in μg/m3. Source: European Environmental Agency

The present paper analyses the introduction of the first phase of an air quality control policy in the municipality of Milan, which started on 25 February 2019 and directly acts on traffic rules. The administration defined an extended limited traffic zone, named Area B (https://www.comune. milano.it/aree-tematiche/mobilita/area-b), where the access and circulation for the most polluting vehicles, as well as those longer than 12 meters, have some partial restrictions, enforced by a monitoring system of entrance gates controlling each license plate and imposing a fine on unauthorized vehicles. The access prohibition concerns Euro 0 petrol vehicles and a large part of Euro 0, 1, 2, and 3 diesel vehicles, with specific exemptions for public transport, itinerant traders, and residents, and it is active from Monday to Friday during business hours (from 7:30 a.m. to 7:30 p.m.), except holidays. According to the municipality of Milan, the share of cars registered in the Milan metropolitan area and involved in the restrictions in 2019 is close to 17%, while the share of freight transport vehicles is around 53% [11]. Area B is a progressive policy divided into various phases, which will concern an increasing number of vehicle classes. In terms of NOx emissions, the administration expects a reduction of 4–5% per year until 2022 and a reduction of 11% between 2023 and 2026. The policy will be fully operative within October 2030.

Area B extends the previously existing limited traffic zone, Area C, which covers just the historical city centre. The physical coverage of the two restriction zones is represented in Figure 2, which highlights the arrangemen<sup>t</sup> of both within the city borders. Area B covers almost the entire area of the city, excluding extreme peripheral districts.

Statistical literature on air quality grew up sharply in the last decades. Two main statistical modeling directions have been developed. One has a focus on pollutants concentration and the other on human exposure. Regarding the latter, recent advances are based on crowdsourced data, such as smartphone data modeling [12]. Regarding pollutants concentrations, increasing attention is being given to latent component models; see, as an example [13] and for the problem of misalignment. In particular, the use of the INLA-SPDE approach for misalignment between pollutant concentration and epidemiological data [14] and PCA based methods with missing data [15].

When the territory under study is large and spatial correlation is important, spatio-temporal models are appropriate. See, for example, the multivariate state space approach of Calculli et al. [16], which is capable of handling jointly PM10, NO2 and weather variables, the approach of Menezes et al. [17] for modeling daily NO2 trends in Portugal. Moreover, the land-use regression model (LUR) under a state space approach has been used for modeling air pollution in Tehran [18]. Despite this growing spatial literature, time series analysis methods have been recently developed to understand the effect of meteorology on pollutant concentration [19], which will be the main focus of this paper.

The previous Milan limited traffic zone, known as Area C, has already been treated in literature by Fassò [20], who analyzed its introduction through spatio-temporal models, by Invernizzi et al. [21], who considered its impact on black carbon, and by Percoco [22] who considered its effect on traffic. Moreover, similar problems have been studied for London "sulphur-free zone" [23] and the "low emission zone" in Munich [24]. In Fassò [20], the author considered both particulates and nitrogen oxides and observed the presence of a more pronounced permanent reduction of the latter within the restricted area, despite the data showing a strong spatial variability depending on the type of pollutant. This is consistent with the known emissions pattern of particulate matters and nitrogen oxides. The latter are mainly primary gaseous pollutants and can be directly attributed to anthropogenic sources, such as car traffic and house heating. Moreover, from the so-called INEMAR emission inventory [25], in Milan province, 68% of NOx and only 41% of PM10 are due to road traffic. Hence, in this first study of Area B, we will take into account NOx and NO2 and postpone the analysis of PM10 and PM2.5 to further research. To adjust for confounding factors, we will consider weather conditions in Milan, the main calendar events, and the concentration levels of oxides observed in neighbouring towns, as in a pseudo-treatment-control approach.

The study aims to identify and quantify variations in pollutant levels due to the above described Area B. Hence, the present paper will try to investigate and test the following two scientific hypotheses:

**Hypothesis 1.** *The introduction of Area B achieved significant changes in pollution concentration for the city of Milan;*

**Hypothesis 2.** *The variation occurred homogeneously on the territory and the stations do not show spatial variability of the effects.*

The first hypothesis aims to quantify the impact of the policy on pollution levels measured by several air quality stations scattered around the city and to assess whether this evidence is significantly supported by the data. The impact is evaluated both regarding the statistical significance of the estimates, the absolute magnitudes of the coefficients, and their signs. From the policy maker perspective, the expected coefficients should be negative, indicating a reduction effect on concentrations due to the car traffic restrictions. However, given the complexity of the phenomenon, a change of opposite sign cannot be ruled out either. The second research hypothesis is dedicated to the comparison of the estimates for the considered stations: the effect can be considered homogeneous when the sign and the magnitude of the coefficients for all the stations are similar.

The paper is structured as follows. Section 2 describes the dataset and the methodologies used for the analyses. In particular, we briefly explain the composition of both weather and air quality monitoring systems in Milan, available data sources, and metadata information. Then, we present the methodologies implemented for the preliminary analysis and the state space approach to time series analysis for air quality data. Section 3 reports and discusses the empirical results of the estimated models and their implications. Section 4 concludes the paper discussing the two research hypotheses in light of what emerged from the data analysis and gives some hints for future research developments.

**Figure 2.** Monitoring system in Milan. Air quality stations (blue points): Marche (501), Verziere (528), Senato (548), Liguria (539) and Città Studi (705). Weather stations (red points): Lambrate (100), Zavattari (503), Brera (620), Feltre (869), Rosellini (1327), Juvara (502), and Marche (501).

#### **2. Materials and Methods**

In this section, we present the structure of the ARPA dataset and briefly introduce the methodologies implemented for the analyses. Section 2.1 introduces the data source for the Milan case study and the spatio-temporal structure of the data and provides a brief description of the variables taken into consideration. Section 2.2 designs the preliminary analysis, which introduces a temporal treatment-control experiment to highlighting the differences in concentration levels before and after the policy intervention. Section 2.3 gives a detailed overview of the use of state space models in time series analysis for the study of air quality data, including also a specific subsection for model selection and policy intervention.
