Start Submission Become a Reviewer

Reading: Analysis and optimal design of air quality monitoring networks using a variational approach

Download

A- A+
Alt. Display

Original Research Papers

Analysis and optimal design of air quality monitoring networks using a variational approach

Authors:

Adolfo Henriquez ,

Center for Climate and Resilience Research (CR2), Blanco Encalada 2002; Centro de Modelamiento Matemático, Universidad de Chile-CNRS, Blanco Encalada 2120, piso 7, CL
X close

Axel Osses,

Center for Climate and Resilience Research (CR2), Blanco Encalada 2002; Centro de Modelamiento Matemático, Universidad de Chile-CNRS, Blanco Encalada 2120, piso 7; Departamento de Ingeniería Matemática, Universidad de Chile, Blanco Encalada 2120, piso 5, CL
X close

Laura Gallardo,

Center for Climate and Resilience Research (CR2), Blanco Encalada 2002; Departamento de Geofísica, Universidad de Chile, Blanco Encalada 2002, CL
X close

Melisa Diaz Resquin

Center for Climate and Resilience Research (CR2), Blanco Encalada 2002, CL; Comisión Nacional de Energía Atómica, Gerencia Química, Av. Gral. Paz 1499 (B1650KNA), San Martín, Pcia; Departamento de Ingeniería Química, Facultad de Ingeniería, Universidad de Buenos Aires, AR
X close

Abstract

Air quality networks need revision and optimisation as instruments and network requirements, both scientific and societal, evolve over time. Assessing and optimising the information content of a monitoring network is a non-trivial problem. Here, we introduce a methodology formulated in a variational framework using an air quality model to simulate the dispersion of carbon monoxide (CO) as a passive tracer at the city scale. We address the specific case of adding or removing stations, and the more general situation of optimally distributing a given number of stations in a domain taking into account transport patterns and spatial factors such as population density and emission patterns. We consider three quality indicators: precision gain, information gain and degrees of freedom for a signal. These metrics are all functions of the singular values of the sensitivity matrix that links emissions and observations in the variational framework. We illustrate the application of the methodology in the case of Santiago (33.5°S, 70.5°W, 500 m a.s.l.), a city of ca. 7 million inhabitants with significant pollution levels. We deem information gain as the best of the above indicators for this case. We then quantify the actual evolution of Santiago's network and compare it with the optimal configuration suggested by our methodology and with results previously obtained using a statistical approach. The application is restricted to diurnal and summer conditions, for which the dispersion model shows a good agreement with observations. The current method offers advantages in that it allows extending a network to include new sites, and it explicitly considers the effects of dispersion patterns, and desired weighting functions such as emission fluxes and population density. We find that Santiago's air quality has improved two-fold since 1988, regarding CO under diurnal summer conditions. Still, according to our results, the current configuration could be improved by integrating more suburban stations in the southwest of the basin.

How to Cite: Henriquez, A., Osses, A., Gallardo, L. and Resquin, M.D., 2015. Analysis and optimal design of air quality monitoring networks using a variational approach. Tellus B: Chemical and Physical Meteorology, 67(1), p.25385. DOI: http://doi.org/10.3402/tellusb.v67.25385
  Published on 01 Jan 2015
 Accepted on 11 Aug 2015            Submitted on 4 Jul 2015

1. Introduction

Observing systems, in general, and air quality networks, in particular, need revision and optimisation as monitoring instruments, scientific and societal requirements change over time. The objectives of an air quality monitoring network are multiple, including compliance of air quality standards, evolution of air quality and efficiency of curbing measures, as well as quantifying impacts on human health, ecosystems and climate. The question of how to best sample airborne pollutants in a monitoring network is non-trivial. The difficulties arise from the questions of defining location, sampling frequency and choice of variable to measure, and so on, in an optimal sense according to different criteria, and the large number of potential sites and times to collect data in a given system. In general, optimal network design is framed as a complex optimisation or decision problem (Caselton and Zidek, 1984; Zidek and Zimmerman, 2010). Currently, urban sustainability in the framework of a changing climate, and the need for combining multiple observing platforms, poses new challenges to the optimal design of networks (Lahoz and Schneider, 2014). There is an increasing body of literature addressing optimal network design in connection with air quality problems. Several approaches have been suggested over the last decades. There are those using pure measurements and statistics, and those using pure numerical simulations, as well as combinations of both. Of course, many of these methodologies allow network reduction or extension and weighted criteria under suitable social, economic or technical restrictions (e.g. Langstaff et al., 1967; Husain and Khan, 1983; Shindo et al., 1990; Trujillo–Ventura and Ellis, 1991; Nychka and Saltzman, 1998; Tseng and Chang, 2001; Bocquet, 2005a, 2005b; Kanaroglou et al., 2005; Kao and Hsieh, 2006; Abida et al., 2008; Elkamel et al., 2008; Abida and Bocquet, 2009; Ainslie et al., 2009; Mofarrah and Husain, 2010; Wu et al., 2010; Romary et al., 2011; Saunier et al., 2011; Sturman et al., 2011; Wu and Bocquet, 2011; Ruiz-Cárdenas et al., 2012; Osses et al., 2013; Singh et al., 2013; Wang et al., 2013; Cioaca and Sandu, 2014; Gómez-Losada et al., 2014; Romary et al., 2014; Zidek et al., 2000, 2014).

In this work, we introduce a methodology formulated in a variational framework, that is, based on the optimisation of quality indicators as functions of the analysis covariance matrix in the context of data assimilation. Therefore, the methodology also incorporates transport and mixing in a city by means of the sensitivity matrix calculated using a chemistry transport model (CTM). In this way, we quantify the changes in quality of a monitoring network when adding or removing stations, and when searching for the optimal spatial distribution of a given number of stations in a domain. We use three quality indicators: precision gain, information gain and degrees of freedom for the signal (DFS) (Shannon, 1948; Rodgers, 2000; Fisher, 2003). These metrics are used in optimal network design, but they are seldom put in a common framework. In fact, we show these metrics in terms of singular values of the sensitivity matrix that links emissions and observations in the variational or data assimilation framework. This methodology seems to be new in the context of network design for atmospheric dispersion. Indeed, even if these indicators (degrees of freedom, information content, total precision) are widely used in meteorological data assimilation, they are not extensively used in atmospheric chemistry data assimilation or atmospheric network design. The methodology also allows to weight in emission fields, population density or other spatially varying characteristics in a straightforward manner.

In sum, the present study includes the intercomparison of several network quality indicators, the use of geographical weights to obtain realistic and feasible networks, and a clear methodology to test the robustness analysis of the proposed optimisation procedure with respect to the initial network.

Furthermore, we apply the proposed methodology at a city scale to the case of Santiago in Chile (33.5°S, 70.5°W, 500 m a.s.l.). The air quality monitoring network of Santiago de Chile, a city of ca. 6 million inhabitants with significant levels of pollution, has evolved as stations have been removed and relocated, and the spatial coverage has grown to include present-day suburban stations, as shown in Fig. 1. These changes have been motivated by the need to monitor the attainment of air quality standards in a city that has rapidly expanded (Romero and Ordenes, 2004), and changed its emission rates and patterns (Gallardo et al., 2012), similar to many other cities in the world (Zhu et al., 2012). For this application, we use a modelling realisation described by Saide et al. (2009) and Saide et al. (2011a). These simulations describe the dispersion of carbon monoxide (CO) over Santiago de Chile under diurnal conditions in a summer month. We further address the evolution of the air quality network in terms of its information gain by taking into account both the number and the location of CO stations. We discuss the pros and cons of the variational approach introduced here and compare it with the results of a statistical approach previously presented by Osses et al. (2013).

Fig. 1  

Santiago's air monitoring sites (black dots) and stations symbols (white letters) for the years 1988, 1997 and 2009. The main topographic features of the Santiago basin are also shown (shading). The approximate urban–rural limit for the years 1989, 1996 and 2010 is shown in a black line (www.ocuc.cl). The area in which stations can be placed in this study is indicated by a dashed line.

In the next section, we introduce various metrics to quantify changes in the information content of a network including precision gain (pg), information gain (ig) and DFS (ds). Section 3 covers how these metrics are applied to address the addition and removal of stations, and to the question of optimal design of a monitoring network, including the use of geographical weighting functions such as population density and emission patterns. In Section 4, we apply this methodology to the case of Santiago's air quality network and its evolution. Summary and conclusions are presented in Section 5.

2. Quantifying the quality of a network

In the context of optimal network design using a variational approach, one encounters the use of various metrics to assess the quality of an observing system. Here, we consider three indicators of relative changes in quality: precision gain, information gain and DFS (e.g. Rodgers, 2000). Along these lines, one can evaluate a set of stations in an air quality network to improve the precision of the estimate of an unknown parameter, for example, emissions of a given tracer, here CO. In this particular case, by assuming a linear relationship between emissions and observations of CO, one links emissions x with measured observations y at the sites of a given network, using the sensitivity matrix H (with rank r>0 and an associate error ɛ

(1 )
y=Hx+ε.

This is a reasonable assumption for CO at the city scale since it behaves as a transport tracer (e.g. Saide et al., 2011b). Note that the sensitivity matrix H is calculated for unitary emissions and evaluated at each observation site. Thus, H depends on mixing and transport, and the location and timing of the observation network.

Before specifying the selected quality metrics, let us introduce the basic definitions and notations of the data assimilation framework (e.g. Daley, 1993; Bouttier and Courtier, 1997; Kalnay, 2003). Let xt be the true n-dimensional state and xb its background estimate with error covariance B, that is, xbN(xt,B). Let yo be the m-dimensional measurement observation vector with error covariance R, that is, yoN(Hxt,R). Let us assume that B and R are positive definite matrices. The analysis (estimate) of the state xa is the unique solution of the optimisation problem

(2 )
f¯(Λ2)1

where xM=xTMx. The analysis is given by

(3 )
xa=xb+ΣaHTR-1d
(4 )
Σa=(B-1+HTR-1H)-1,

where Σa is the error covariance matrix of the analysis [i.e. xaN(xt,Σa)] and d=y0-Hxb is the so-called innovation vector. Note that Σa-B-1=HTR-1H0, indicating that the assimilation of observations always improves the precision of analysis Σa with respect to the precision of the background B−1.

Let us consider the transformation H~=R-12HB-12, where B12 and R12 stand for the square roots from the Cholesky decomposition (e.g. Horn and Johnson, 2012) and the singular value decomposition H~=U~Λ~V~T, with Λ~ the diagonal matrix of non-zero singular values λ~1,,λ~r, and U~ and V~ unitary matrices of singular vectors u~i and v~i, respectively. Also, transforming d~=R-12d and V¯=B12V~, eqs. (3) and (4) become

(5 )
xa=xb+V¯(In+Λ~2)-1Λ~U~Td~=xb+i=1r(λ~i21+λ~i2)u~iTd~λ~iv¯i
(6 )
Σa=B12(In+H~TH~)-1BT2.

The term in parenthesis is known as the Tychonoff filter factor (Johnson, 2003). This factor ponders the contribution of the ith component of the singular value decomposition in the total analysis. When the corresponding λ~i1, this contribution is relevant, whereas when λ~i1 the opposite occurs.

In the particular case when the covariances are scalar matrices B=σb2In and R=σo2Im, with ratio μ2=σo2/σb2, this expression simplifies to

(7 )
xa=xb+i=1r(λi2μ21+λi2μ2)uiTdλivi
(8 )
Σa=σb2V(I+Λ2μ2)-1VT,

using in this case the singular value decomposition H=UΛVT, with Λ a diagonal matrix of non-zero singular values λ1,,λr, and U and V unitary matrices of singular vectors ui and vi, respectively.

Following from this, let Z be the set of locations of the stations of a city's air quality network

(9 )
Z={z1,z2,,zN}

with zi=(xi,yi) and cardinality Z=N. To each set Z, we associate a network of observation sites, and a corresponding sensitivity matrix H and covariance matrix Σa with quality indicator given by Qa)

ZHΣaQ(Σa).

In order to simplify the notation, we will associate to each network Z the corresponding quality indicator Q(Z).

ZQ(Z).

2.1. Precision gain

We introduce total precision of a given network Z characterised by sensitivity matrix H, as the trace of the precision matrix of the corresponding analysis, that is

(10 )
p(Z)=Tr(Σa)=Tr(B-1+HTR-1H),

which can be simplified in the case of diagonal error covariance matrices B and R as follows:

(11 )
p(Z)=nσb2+1σb2i=1rλi2μ2.

Then, precision gain is obtained by subtracting the total precision after assimilation of the network observations from the total precision before assimilation:

(12 )
pg(Z)=Tr(Σa)-Tr(B-1)=1σb2i=1rλi2μ2.

In simple words, this indicator follows the notion that precision increases as errors decrease in a quadratic manner. It is worth noting that this metric is not sensitive to B.

2.2. Information gain

The total information of a random variable x with values in Rn normally distributed with covariance Σx is given by

(13 )
i(x)=12lnΣx-n2(1+ln(2π)).

Therefore, information gain for a network Z characterised by sensitivity H is obtained by subtracting the total information after and before the observations of the network are assimilated:

(14 )
ig(Z)=i(xa)-i(xb)=12lnΣa-12lnB-1=12i=1rln(1+λi2μ2).

An alternative but equivalent definition uses the Kullback-Liebler distance (Kullback, 1959) as described in Osses et al. (2013). Since µ ponders the relative importance of the errors in observations and emissions, the extreme cases correspond to no knowledge of observations relative to emissions, that is, µ→+∞ so ig=0, and perfect knowledge of observations relative to emissions, that is, µ→0 so ig=+∞.

2.3. Degrees of freedom for the signal

The DFS associated with a network Z characterised by the sensibility matrix H are given by

(15 )
ds(Z)=r-Tr(B-1Σa)
(15 )
=r-i=1r(1+λi2μ2)-1,

where r is the rank of H (note that rnm in the application discussed later). DFS represent the number of useful pieces of information present in the noisy measurements of a given network in order to effectively retrieve emissions. It can also be interpreted as the effective rank (but not integer in general) of the sensitivity matrix (see Rodgers, 2000). Extreme cases correspond to no knowledge of observation relative to emissions, that is, µ→+∞ so ds=0, and perfect knowledge of observations relative to emissions, that is, µ→0 so ds=r.

2.4. Comparison among indicators

All previous quality indicators are functions of the analysis covariance matrix Σa associated with network Z, thus non-linear functions of the associated sensitivity matrix H. In the particular case of scalar covariance matrices for observations (R) and emissions (B), these terms can be expressed as the ratio λ i/µ, where λi are the singular values of the sensitivity matrix, and µ is the ratio between observational and background deviations, in the general form if(λi/μ), where the corresponding values for f are indicated in Table 1. A similar analysis can be done in the non-diagonal case by using expressions (5) and (6). Figure 2 compares the behaviour of f for the quality indicators: pg, ig and ds.

Fig. 2  

Comparison of quality indicators as a function of singular values λ. The vertical continuous black line indicates the value of the parameter μ10=0.0075, and the vertical dashed grey line indicates the value of the parameter µ=0.075. The bottom of the figure shows the distribution of the 644 singular values of the sensitivity matrix for Santiago's network in 2009.

Calculating the singular values for the full sensitivity matrix can be cumbersome. For example, in the case of Santiago shown later, the observation operator for the whole domain has 154 560 rows and 9660 columns. Thus, to assert the dependence on singular values, we reviewed a random sample of 250 networks of 11 stations, and for which singular values were calculated. These calculations indicate that, on average, λ¯min10-5 and λ¯max10. Let us now consider the singular values in the regions Λ1=[λ¯min,μ10] and Λ2=[μ10,λ¯max]. In region Λ1, all quality metrics expressed as functions of singular values behave similarly, that is, λ2μ2, and show negligible differences if compared with respect to their values in region Λ2. In region Λ2, the functions behave very differently. In this region, and for the random set of 250 networks, precision gain has the steepest slope. Furthermore, an average of 70% of the spectrum is in region Λ1. The remaining 30% of the spectrum is highly determined by singular values in region Λ2, since f¯(Λ1)10-3 and f¯(Λ2)1. In order to compare the different quality indicators in the range of all scales of singular value, a logarithmic scale in singular values is needed. Precision gain varies quadratically with singular values in region Λ2, whereas DFS show an asymptotic growth, covering the full range of singular values. Information gain shows a linear growth, in between precision gain and DFS. For these reasons, we adopt information gain and DFS as quality indicators to be applied hereafter in the case of Santiago. Thus, for the simulated meteorological and dispersion conditions used for Santiago, the most robust indicators with respect to the leading singular values correspond to information gain and DFS, even if DFS is less sensitive to changes in singular values than information gain.

3. Methodology

In the previous section, we defined a number of metrics to address information content of a monitoring network. Each quality indicator Q(Z) is associated, via the covariance matrix of the analysis, with some particular network Z={z1,z2,,zN} characterised by the geographical location of each station zi=(xi,yi).

We first present the methodology for removing or adding one station, taking into account the influence of transport fields as well as the original geographical distribution of the network. Thereafter, we consider the problem of optimally distributing in space a given number of stations, taking into account circulation patterns, and weighting functions such as population density and emission fluxes. The inclusion of weighting functions allows a straightforward manner to ponder the importance of different factors, which in turn results in targeted optimal networks.

3.1. Reduction, extension and optimal distribution

Reduction. Consider a quality indicator Q. The reduction of a network Z by subtracting a station located at Zi={(xi,yi)} in the discrete spatial domain will correspond to a new network ZZi. For each station that can be extracted from the network, we calculate the decrease in quality by

(16 )
Q(Z)-Q(ZZi).

Thus, the greater the decrease in quality for a particular station, the lower its priority to be removed.

• Extension. Consider a quality indicator Q. The extension of a network Z by adding a new station located at Z0={(x0,y0)} in the discrete spatial domain will correspond to the new network ZZ0. For each possible new station, we calculate the increase in quality by

(17 )
Q(ZZ0)-Q(Z).

Thus, the greater the increase in quality for a potential new station, the greater its priority to be added.

• Optimal distribution. Consider a quality indicator Q and a network size N. We find the optimal network by maximising the network's quality over all networks Z of size N in the discrete spatial domain D, according to

(18 )
maxZDZ=NQ(Z)

This optimisation problem has two major difficulties. First, the feasible set is discrete, and second the direct evaluation of all possible values of the objective function is computationally cumbersome. Therefore, we adopt a heuristic optimisation algorithm known as simulated annealing (SA) (e.g. Bertsimas and Tsitsiklis, 1993). This algorithm gives an estimate of the optimal configuration in a manageable computation time even with a personal computer. In the Santiago case, for a network of 11 stations, and a sensitivity matrix of size 2640×644, corresponding to 11 observation sites over 10 d and 24 hours, and 28×23 emission grids. We simplified the emissions by averaging over time or by choosing a specific hour. We solved this problem in 10 hours using a system of two cores of 2.3 GHz.

3.2. Simulated annealing

The procedure of SA starts from an initial feasible network and, at each iteration, it chooses, according to a given probability, to move to another feasible network in the vicinity of the current network or stay in the current network. The mobility of the search algorithm depends on a global parameter T and it is called temperature. Initially SA uses large values of T, allowing a random search. Then, as the values of the objective function increase, the value of T decreases, restricting the search and ensuring convergence to a local optimum.

For a network quality indicator Q, a feasible network Z1, a test feasible network Z2, and δ=Q(Z2)-Q(Z1), we define the probabilities

(19 )
ρ={1ifδ0eδTifδ<0
(20 )
ρaU([0,1]).

If δ≥0, SA moves from Z1 to Z2, increasing the value of Q. If δ<0 and ρ>ρa, SA moves from Z1 to Z2, without increasing the value of Q, but allowing the function to increase its value in future iterations. If δ<0 and ρ<ρa, SA stays in Z1 and chooses another test network. As T decreases, the latter ensures the convergence of SA near the local optimum.

The max–min network, which is a solution of the problem

(21 )
maxZDZ=Nminz1,z2Zz1z2z1-z2,

is a network equally spaced along the observation domain. This is an adequate initial condition for SA.

3.3. Geographical weighting

When designing monitoring networks, in practice one has to consider different criteria and objectives, as well as the geographical distribution of information, such as emissions, population density, socioeconomic indexes, among others. To do so, we apply suitable geographical weights to the sensitivity matrix. These weights π=(πi)i are modulated by a penalisation parameter β≥0, according to

(22 )
Hβ=HΠβ

where Πβ=diag(π1β2,,πnβ2), so we apply a weight πjβ2 to each grid point j. Hence, the previous definitions are slightly modified by changing H by Hβ. In particular, this implies that the singular values of Hβ depend on the value of β. Recall that for a linear tracer such as CO at the city scale, the sensitivity matrix depends on meteorological parameters. Thus, the inclusion of parameter β alters the sensitivity matrix making it sensitive to geographical weights. This is a simple way to include geographical weights and it is equivalent to a change of variables in the parameter space with the same observation equation: y=Hβx+ε, where xβ=Πβ-1x.

Figure 3 shows two geographical weights used in this work, namely population density as in the census for 2002 (www.ine.cl) and CO emission fluxes by Saide et al. (2011b). The results are strongly modulated by the selection of the penalisation parameter β. This parameter controls the relative importance used to correct the sensitivity matrix given the minima and maxima of a given geographical field. Note that it is not possible to adjust β in order to balance the different quality indicators uniformly with respect to the network since the dependence with respect to β is logarithmic for information quantity and quadratic for total precision. The weight and the corresponding chosen value of β should be defined by the user and it will depend on the specific objective of the network. In our case, we will take a value of β that produces an even quality gain for urban and rural areas when adding a new station (see Section 4.3 and Fig. 5). In any case for β=0, the weighting field is not taken into account (uniform weight). In the case of population and emission distributions, we chose β as a trade-off between urban and suburban areas; this corresponds to β=0.7 for information gain and β=1.2 for DFS.

Fig. 3  

Weighting functions applied in this work. The upper panel shows logarithm of population density (hab/km2). The lower panel shows logarithm of normalised CO summer emission fluxes (mol km−2 h−1). The white contour indicates the approximate urban–rural limit of Santiago in 2010. White circles show the location of monitoring stations in 2009.

4. Results

As described earlier, we consider information gain to be the best metric for assessing the quality of a network for the case of Santiago de Chile under summer conditions. Therefore, hereafter we show results using this indicator. We performed several tests and found similar results using DFS (not shown).

4.1. Experimental setting

For simulating CO dispersion over Santiago, we adopt the same modelling tool used by Saide et al. (2009) and Saide et al. (2011a), currently used by the Chilean Weather Office for forecasting ozone in summer. This system uses an air quality platform consisting of the CTM Polyphemus (Mallet et al., 2007) fed with meteorological fields calculated with the MM5 model (Grell et al., 1995). The simulations were performed for the 15–25 January 2002 period using four nested grids with horizontal resolution of 54, 18, 6 and 2 km, respectively. The 2 km resolution domain, centred in the city of Santiago, was fed into Polyphemus. The Polyphemus domain is also centred at Santiago, with 70×63 grid cells of 2 km resolution, and 12 vertical levels: from 13 to 6000 m above the surface, with finer resolution in the boundary layer. These simulations used zero initial and boundary conditions. Given the relatively clean Pacific air that arrives in Santiago, these assumptions appear as reasonable (e.g. Allen et al., 2011).

The sensitivity matrix (H) was calculated using unitary emissions in a domain of 28×23 cells depicted in Fig. 1. For building (H), we followed the same procedure as described by Saide et al. (2011a) for diurnal conditions between 6:00 and 20:00 hours, local time. Following these authors, for simplicity and lack of better knowledge, B and R are treated as scalar matrices, being the error in emissions 20%, and the error in observations 5%. Note that the null boundary conditions for CO are far away from the emission grid, since the CTM simulations were performed in a bigger and nested mesh. We point out that for the geographical weights considered in this study, and for the different initial networks used, the final optimal networks sometimes presented few stations located near the boundary of the emission grid (see Figs. 6 and 9) but not affecting the main orientation or distribution of the networks identified in Fig. 7. Therefore, for the purposes of this study, that is, influence of meteorology in the orientation of the network, influence of the geographical weight, robustness of the optimisation procedure, and so on, the chosen emission grid was adequate.

4.2. Reduction of the network

According to our calculations, if one were to remove a monitoring station from Santiago's network as in 2009, using population density as the relevant weighting function, one would conclude that the removal of station O would result in the least decrease in quality of the network (Fig. 4, left panel). In contrast, one would not suggest removing stations Q, N, M, L or F, removal of which would result in a significant loss in network quality. If instead, the factor under consideration was CO emissions, one would still remove O and suggest keeping S, Q, P, V, M and T. Thus, depending on the network's objectives, one would prioritise to keep different sets of stations.

Fig. 4  

Decreases in network's quality after removing stations. The centre of the circles indicates the position of the stations and each radius is a linear factor of the decrease in quality. The left (right) panel corresponds to the use of population density (emission fluxes) as weighting factors.

It is interesting to note that the potential removal of stations O and R has the least impact on the quality of the network for assessing CO, either considering population density or emission fluxes. This is independent of their close geographical location. In fact, if we remove either station or if we merge the two stations by averaging their sensitivity values, the final quality decrease map is similar to the one shown in Fig. 4. This means that O and R share information and that the information provided by these stations is less relevant to the network's quality. This is firstly because the summer circulation in Santiago shows a clear southwest to northeast axis (following the predominant wind direction) and secondly because neither population nor emission density is distinct over this area. The former results in a relatively larger sensitivity of stations along the circulation axis and he latter results in too small a weighting factor to change the sensitivity for O and R. In the case of station S in the southeast of Santiago, the weighting factor is significant enough to counteract the fact that this station is off the predominant wind direction. Station T, which is located along the air circulation (SW to NE), increases its importance when emission criteria are prioritized because it is located close to a major highway with relatively large emissions. To verify these results, we calculated our metrics for the case β=0, finding that that stations located in areas with the lowest density of stations (e.g. T) are those that have the largest impact on the quality of the network. It must be pointed out that these results use a sensitivity matrix restricted to diurnal conditions, when circulation and stability patterns are better captured by the dispersion model (e.g. Saide et al., 2011a). Thus, by no means do we question the location of stations O and R in the west of the city, which in fact capture the nocturnal maxima in pollution conditions (e.g. Saide et al., 2011b).

All in all, our results for summer conditions suggest that downtown stations are important from the point of view of total information content regarding CO. A previous study that considered multiple species using statistical information indexes by Osses et al. (2013) also suggested retention of downtown stations. These stations were in fact deemed to be the most representative of dispersion patterns for the Santiago basin, considering all species and all hours of the day, and all seasons of the year.

4.3. Extension of the network

Here, we address the question of where to add stations so that the quality of the network increases by growing information content. In other words, we calculate the information gain provided by adding a monitoring station in every grid box lacking monitoring at present. According to the values of the sensitivity matrix and the weighting function used, one may determine the change in information gain. In this way, we can build an information gain map as shown in Fig. 5. Both weighting functions suggest the addition of stations in the urban north and northeast areas of the basin, and in the suburban areas to the west and south of the city. This follows from the coincident geographical patterns of population and emission fields. In fact, the main difference between the suggested new sites is in a location southwest of the city, which appears when weighting the sensitivity matrix by emissions along the highway, as was mentioned earlier. Interestingly, the extension of the network that took place in 2009 did add station T in the southwest of the basin, a conclusion that is supported by evaluating the information gain for the 2005 network. This is discussed in the following section.

Fig. 5  

Spatial distribution of increases in quality (shaded contours) for the 2009 network (white circles), after adding stations. White squares suggest the potential location of new stations coinciding with maxima in information gain. The left (right) panel shows the corresponding results using population density (emission fluxes) as weighting functions.

4.4. Optimal network

Next, we address the question of designing the optimal network for a given sensitivity matrix, that is, we search for the network that maximises its quality indicator, in this case information gain. To do so, we start with a network of 11 stations evenly distributed in the domain using a max–min approach. We then optimise the network using simulated annealing. As in the previous section, the resulting optimal network follows from both the dispersion patterns embedded in the sensitivity matrix and the weighting function used. If we use H with no weighting or a radial weighting centred in the middle of the domain, the optimal network shows an axis along the prevailing wind direction, with a more even distribution in the morning hours than in the afternoon hours (not shown). In the morning, the low wind speeds (≤3 m/s) and inhibited vertical mixing result in a larger sensitivity to the impact of neighbouring grid boxes compared to the impact of grid boxes farther away. Under these conditions, the resulting optimal distribution shows a larger spatial spread, that is, it privileges the coverage of the domain. In the afternoon, high wind speeds and vigorous vertical mixing result in a less geographically spread distribution. But since no weight is used, the stations tend to be located outside the city. Considering the modified sensitivity matrix by weighting functions, the corresponding optimal networks have similar characteristics, but the coverage of the domain is now complemented with stations in the city according to the weighting functions. These results are illustrated in Fig. 6.

Fig. 6  

Optimal networks according to wind pattern. Upper row: morning sensitivity weighted by population density (left), by CO emission fluxes (centre), and morning winds (right). Lower row: the same but for the afternoon sensitivity and afternoon winds.

4.5. Intercomparison of optimal networks

A network Z={z1,z2,,zN} can be characterised by two independent geometrical parameters. One corresponds to the Hausdorff distance with respect to a reference network Z′ (Huttenlocher et al., 1993), defined as:

(23 )
dH(Z,Z)=max{supzZinfzZz-z,supzZinfzZz-z},

which can be loosely defined as the maximum distance between two separated points but still close enough to be members of the same set of points (network). A second parameter is the spread or spatial dispersion of the network:

(24 )
1NzZz-z¯2,

where z¯ corresponds to the centroid of the network. Other shape indicators such as centroid, skewness, main direction and eccentricity may also be used. In Fig. 7, we show the corresponding Hausdorff distance and geometrical dispersion for a set of 18 optimal networks according to different hours of the day, initial networks for SA and quality indicators. These optimal networks were generated starting from three very different configurations, namely, the actual Santiago network in 2009, the maximum coverage configuration (max–min approach) and a synthetic network. The synthetic configuration has an orientation roughly perpendicular to Santiago's network. We determine optimal networks using both information gain and DFS as quality indicators. Furthermore, to characterise the influence of wind fields, we present the results for morning and afternoon hours, and also for the average diurnal conditions. The Hausdorff distance is calculated with respect to Santiago's network in 2009. In all cases, we use population density as weighting function. Independent of the starting network, hour of the day and quality indicator, all optimal networks have more spatial coverage (dispersion) and are distant from the current configuration of Santiago's network. In other words, according to our calculations, the current configuration is suboptimal. Depending on the wind and dispersion conditions, the optimal networks may require more (morning) or less (afternoon) spatial coverage as discussed in Section 4.4. Note that both information gain and DFS are robust quality criteria with respect to the choice of the initial network and converge to optimal networks with similar characteristics as shown in Fig. 7.

Fig. 7  

Comparison among optimal networks using Hausdorff distance (in °), and geometrical dispersion (in °) for a set of 18 optimal networks for different hours of the day, initial conditions for SA and quality indicators. Black diamonds indicate the initial configurations used in the optimisation procedure: 1) Santiago 2009, 2) max–min and 3) synthetic. A sketch of each network is also shown. Stars indicate the optimal networks for morning hours, circles correspond to afternoon hours and triangles to average diurnal conditions. Networks obtained using information gain (degrees of freedom for signal) are shown by grey shaded (white) symbols.

4.6. Evolution of a network

Finally, we compare the actual evolution of Santiago's air monitoring network with that of the optimal network according to our calculations for summer conditions. We use information gain weighted by emissions, starting each time from the actual distribution and adding/removing the number of stations actually added/removed over time in Santiago. Except in 2002, when stations were removed, the quality of the network has improved over time, with today's network being roughly twice as informative as the network of 1988 (Fig. 8).

Fig. 8  

Historical (dashed line) and optimal (continuous line) evolution of Santiago's air quality network using information gain and CO emission fluxes as weighting function. All values are normalised with respect to the optimal future scenario.

Note that the methodology to obtain an optimal network of 11 stations is different from the methodology to evolve a network by adding stations ‘one to one’ in an optimal way. In principle, one could obtain different results. Here we noted that the resulting networks are similar in shape and quality.

Figure 9 illustrates the spatial distribution of stations according to the actual evolution and the corresponding optimal evolution, considering information gain and a sensitivity matrix pondered by emissions. Each optimisation contemplates the same number of added stations as in the actual evolution for the transitions between 1988 and 1997, and between 1997 and 2009. We include a future scenario considering the addition of three stations, both starting with the 2009 configuration. Our results suggest adding stations in the southwest and northeast of the city. In fact, this example shows how to construct and justify a feasible change in the actual real network by mixing the different methodologies of the previous sections (adding, optimising, evolving).

Fig. 9  

Evolution of Santiago's air quality network. The upper (lower) panel shows the actual or historical (optimal) evolution. New stations are highlighted by a black square.

5. Summary and conclusions

We have presented a methodology to evaluate the quality of an air pollutant monitoring network. The methodology uses quality metrics calculated on the basis of a sensitivity matrix estimated by means of dispersion modelling, that is, it is formulated in a variational framework as emissions and observations are connected by means of a sensitivity matrix calculated with a CTM. The metrics evaluated are precision gain, information gain and DFS. These metrics are expressed in a common framework defined by the singular values of the sensitivity matrix and the covariance matrices that characterise the errors in emissions and observations. The methodology allows the assessment of adding or removing stations in a network and the more general situation of optimally distributing a given number of stations in a domain taking into account transport patterns in a city. Compared with statistical approaches addressing the optimal design of monitoring networks, the use of a variational framework is advantageous in that it explicitly factors in air circulation patterns and it allows the assessment of future scenarios. However, it requires the use of CTMs to calculate sensitivity matrices, albeit for linear tracers such as CO at the city scale.

We applied the methodology to evaluate the evolution and optimal design of Santiago de Chile's air quality network. We show evidence of the importance of considering multiple factors when designing and modifying a network. In particular, we address the impacts of transport patterns and of geographically varying characteristics such as population density and emissions. For instance, downtown stations in Santiago are key to the quality of the network, a finding that is consistent with previous studies (Osses et al., 2013). This follows from their location along the main axis of air circulation in a relatively flat area of the basin, and the high emission and population densities. Also, as the spatial coverage of the network has increased by adding peripheral stations (e.g. T, V, S), so has its quality increased as determined by information gain and other metrics. In fact, according to our results, the quality of the network has increased over time, with the present configuration being twice as informative for CO as it was in 1998. Still, the evolution has been suboptimal, and it can be significantly improved by integrating more suburban stations in the southwest of the basin. Again, these conclusions refer to diurnal conditions in summer. The purpose of this work was nevertheless to review the use of various quality metrics for characterising monitoring networks and to illustrate its use in a coherent framework. Future works should address the development of more accurate scenarios, consider nocturnal conditions and factor in the seasonal variability in circulation patterns, as well as other design objectives such as identifying maximum impact areas in terms of air pollution.

6. Acknowledgements

This work was initially carried out with the aid of a grant from the Inter-American Institute for Global Change Research (IAI) CRN II 2017 which is supported by the US National Science Foundation (Grant GEO-0452325), and it was completed with support by Conicyt/FONDAP/15110009. The authors greatly appreciate comments and suggestions from Dr. N. Huneeus, Dr. K. Pistone and Dr. P. Saide. The authors also appreciate the comments and suggestions by two anonymous referees. The second author acknowledges the partial support of Fondecyt 1151512 and ACT 1106 Conicyt grants.

References

  1. Abida R. , Bocquet M . Targeting of observations for accidental atmospheric release monitoring . Atmos. Environ . 2009 ; 43 : 6312 – 6327 .  

  2. Abida R. , Bocquet M. , Vercauteren N. , Isnard O . Design of a monitoring network over France in case of a radiological accidental release . Atmos. Environ . 2008 ; 42 ( 21 ): 5205 – 5219 .  

  3. Ainslie B. , Reuten C. , Steyn D. G. , Le N. D. , Zidek J. V . Application of an entropy-based Bayesian optimization technique to the redesign of an existing monitoring network for single air pollutants . J. Environ. Manage . 2009 ; 90 ( 8 ): 2715 – 2729 .  

  4. Allen G. , Coe H. , Clarke A. , Bretherton C. , Wood R . South East Pacific atmospheric composition and variability sampled along 20 S during VOCALS-Rex . Stat. Sci . 2011 ; 8 ( 1 ): 10 – 15 .  

  5. Bertsimas D. , Tsitsiklis J . Simulated annealing . Atmos. Chem. Phys . 1993 ; 11 : 5237 – 5262 .  

  6. Bocquet M . Reconstruction of an atmospheric tracer source using the principle of maximum entropy. I: theory . Q. J. Roy. Meteorol. Soc . 2005a ; 131 : 2191 – 2208 .  

  7. Bocquet M . Reconstruction of an atmospheric tracer source using the principle of maximum entropy. II: applications . Q. J. Roy. Meteorol. Soc . 2005b ; 131 : 2209 – 2223 .  

  8. Bouttier P. , Courtier F . Data assimilation: concepts and methods . Training Course Notes of ECMWF . 1997 ; UK : Reading . Vol. 53 .  

  9. Caselton W. , Zidek J . Optimal monitoring network designs . Stat. Probab. Lett . 1984 ; 2 : 223 – 227 .  

  10. Cioaca A. , Sandu A . Low-rank approximations for computing observation impact in 4D-Var data assimilation . Comp. Math. Appl . 2014 ; 67 : 2112 – 2126 .  

  11. Daley R . Atmospheric Data Analysis . 1993 ; Cambridge University Press, Cambridge .  

  12. Elkamel A. , Fatehifar E. , Taheri M. , Al-Rashidi M. S. , Lohi A . A heuristic optimization approach for Air Quality Monitoring Network design with the simultaneous consideration of multiple pollutants . J. Environ. Manage . 2008 ; 88 : 507 – 516 .  

  13. Fisher M . Estimation of entropy reduction and degrees of freedom for signal for large variational analysis systems . ECMWF Tech. Memo. No. 397 . 2003 ; 18 . Online at: http://old.ecmwf.int/publications/library/ecpublications/_pdf/tm/301-400/tm397.pdf .  

  14. Gallardo L. , Dawidowski L. , Rojas N. J. , Andrade M. F. , Osses M . Evaluation of vehicle emission inventories for carbon monoxide and nitrogen oxides for Bogo[tacute]a, Buenos Aires, Santiago, and Sao Paulo . Atmos. Environ . 2012 ; 47 : 12 – 19 .  

  15. Gómez-Losada A. , Lozano-García A. , Pino-Mejías R. , Contreras-González J . Finite mixture models to characterize and refine air quality monitoring networks . Sci. Total Environ . 2014 ; 485–486 : 292 – 299 .  

  16. Grell G. A. , Dudhiao J. , Stauffe D. R . A Description of the Fifth Generation Penn State/NCAR Mesoscale Model (MM5) . 1995 . NCAR Tech. Note TN-398+ STR, 1995 .  

  17. Horn R. , Johnson C . Matrix Analysis . 2012 ; Cambridge University Press, Cambridge .  

  18. Husain T. , Khan H . Shannon's entropy concept in optimum air monitoring network design . Sci. Total Environ . 1983 ; 30 : 181 – 190 .  

  19. Huttenlocher D. , Klanderman G. , Rucklidge W . Comparing images using the Hausdorff distance . IEEE Trans. Pattern Anal. Mach. Intell . 1993 ; 15 ( 9 ): 850 – 863 .  

  20. Johnson C . Information Content of Observations in Variational Data Assimilation . 2003 ; PhD Thesis. Department of Meteorology, The University of Reading, UK .  

  21. Kalnay E . Atmospheric Modeling, Data Assimilation and Predictability . 2003 ; Cambridge University Press, Cambridge .  

  22. Kanaroglou P. S. , Jerret M. , Morrison J. , Beckerman B. , Arain M. A. , co-authors . Establishing an air pollution monitoring network for intra-urban population exposure assessment: a location–allocation approach . Atmos. Environ . 2005 ; 39 : 2399 – 2409 .  

  23. Kao J. J. , Hsieh M. R . Utilizing multiobjective analysis to determine an air quality monitoring network in an industrial district . Atmos. Environ . 2006 ; 40 : 1092 – 1103 .  

  24. Kullback S . Information Theory and Statistics . 1959 ; Wiley, New York .  

  25. Lahoz W. A. , Schneider P . Data assimilation: making sense of Earth Observation . Front. Environ. Sci . 2014 ; 2 : 16 .  

  26. Langstaff J. , Seigneur C. , Mei-Kao L. , Behar J. , McElroy J. L . Design of an optimum air monitoring network for exposure assessments . Atmos. Environ . 1967 ; 21 : 1393 – 1410 .  

  27. Mallet V. , Quélo D. , Sportisse B. , Ahmed de Biasi M. , Debry E. , co-authors . Technical note: the air quality modeling system Polyphemus . Atmos. Chem. Phys. Discuss . 2007 ; 7 : 6459 – 6486 .  

  28. Mofarrah A. , Husain T . A holistic approach for optimal design of air quality monitoring network expansion in an urban area . Atmos. Environ . 2010 ; 44 ( 3 ): 432 – 440 .  

  29. Nychka D. , Saltzman N . Design of air-quality monitoring networks . Case Stud. Environ. Stat . 1998 ; 132 : 51 – 76 .  

  30. Osses A. , Gallardo L. , Faundez T . Analysis and evolution of air quality monitoring networks using combined statistical information indexes . Tellus B . 2013 ; 65 : 19822 . DOI: http://dx.doi.org/10.3402/tellusb.v65i0.19822 .  

  31. Rodgers C . Inverse methods of atmospheric sounding, theory and practice . Atmos. Ocean. Planet. Phys . 2000 ; 2 : 12 – 41 .  

  32. Romary T. , De Fouquet C. , Malherbe L . Sampling design for air quality measurement surveys: an optimization approach . Atmos. Environ . 2011 ; 45 : 3613 – 3620 .  

  33. Romary T. , Malherbe L. , Fouquet C . Optimal spatial design for air quality measurement surveys . Environmetrics . 2014 ; 25 : 16 – 28 .  

  34. Romero H. , Ordenes F . Emerging urbanization in the Southern Andes: environmental impacts of urban sprawl in Santiago de Chile on the Andean piedmont . Mt. Res. Dev . 2004 ; 24 : 197 – 201 .  

  35. Ruiz-Cárdenas R. , Ferreira M. , Schmidt A . Evolutionary Markov chain Monte Carlo algorithms for optimal monitoring network designs . Stat. Methodol . 2012 ; 9 ( 1–2 ): 185 – 194 .  

  36. Saide P. , Bocquet M. , Osses A. , Gallardo L . Constraining surface emissions of air pollutants using inverse modeling: method inter-comparison and a new two-step multiscale approach . Tellus B . 2011a ; 63 : 360 – 370 .  

  37. Saide P. , Carmichael G. , Spak S. , Gallardo L. , Osses A. , co-authors . Forecasting urban PM10 and PM2.5 pollution episodes in very stable nocturnal conditions and complex terrain using WRF-Chem CO tracer model . Atmos. Environ . 2011b ; 45 : 2769 – 2780 .  

  38. Saide P. , Osses A. , Gallardo L. , Osses M . Adjoint inverse modeling of a CO emission inventory at the city scale: Santiago de Chiles case . Atmos. Chem. Phys. Discuss . 2009 ; 9 : 6325 – 6361 .  

  39. Saunier O. , Bocquet M. , Mathieu A. , Isnard O . Model reduction via principal component truncation for the optimal design of atmospheric monitoring networks . Atmos. Environ . 2011 ; 43 : 4940 – 4950 .  

  40. Shannon C . A mathematical theory of communication . Bell Syst. Tech. J . 1948 ; 27 ( 3 ): 379 – 423 .  

  41. Shindo J. , Oi K. , Matsumoto Y . Considerations on air pollution monitoring network design in the light of spatio-temporal variations of data . Atmos. Environ. B Urban Atmos . 1990 ; 24 : 335 – 342 .  

  42. Singh K. , Sandu A. , Jardak M. , Bowman K. , Lee M . A practical method to estimate information content in the context of 4D-var data assimilation . SIAM/ASA J. Uncertain. Quantif . 2013 ; 1 : 106 – 138 .  

  43. Sturman A. , Titov M. , Zawar-Reza P . Selecting optimal monitoring site locations for peak ambient particulate material concentrations using the MM5-CAMx4 numerical modelling system . Sci. Total Environ . 2011 ; 409 : 810 – 821 .  

  44. Trujillo-Ventura A. , Ellis J. H . Multiobjective air pollution monitoring network design . Atmos. Environ. A Gen Topics . 1991 ; 25 : 469 – 479 .  

  45. Tseng C. , Chang N. B . Assessing relocation strategies of urban air quality monitoring stations by GA-based compromise programming . Environ. Int . 2001 ; 26 : 523 – 541 .  

  46. Wang K. , Ding Y. , Zhao H. , Hou L. , Sun F . Optimization of air pollutant monitoring stations based on genetic algorithm . Proceedings of the Emerging Intelligent Data and Web Technologies (EIDWT) Fourth International Conference on 2013 . 2013 . IEEE, 680–684. .  

  47. Wu L. , Bocquet M . Optimal redistribution of the background ozone monitoring stations over France . Atmos. Environ . 2011 ; 45 ( 3 ): 772 – 783 .  

  48. Wu L. , Bocquet M. , Chevallier M . Optimal reduction of the ozone monitoring network over France . Atmos. Environ . 2010 ; 44 ( 25 ): 3071 – 3083 .  

  49. Zhu T. , Melamed M. L. , Parrish D. , Gauss M. , Gallardo L. , co-authors . WMO/IGAC impacts of megacities on air pollution and climate. (GAW Report No. 205, 205) . 2012 ; Geneva : World Meteorological Organization .  

  50. Zidek J. V. , Shaddick G. , Taylor C. G . Reducing estimation bias in adaptively changing monitoring networks with preferential site selection . Ann. Appl. Stat . 2014 ; 8 : 1640 – 1670 .  

  51. Zidek J. V. , Sun W. , Le N. D . Designing and integrating composite networks for monitoring multivariate Gaussian pollution . Appl. Stat . 2000 ; 49 ( 1 ): 63 – 79 .  

  52. Zidek J. V. , Zimmerman D. L . Gelfand A. E. , Diggle P. J. , Fuentes M. , Guttorp P . Chapter 10. Monitoring network design . Handbook of Spatial Statistics . 2010 ; CRC Press, Taylor and Francis Group, New York . 131 – 148 .  

comments powered by Disqus