Oxidative potential apportionment of atmospheric PM₁: a new approach combining high-sensitive online analysers for chemical composition and offline OP measurement technique

Source apportionment models were widely used to successfully assign highly time-resolved aerosol data to specific emissions and/or atmospheric chemical processes. These techniques are necessary for targeting the sources affecting air quality and for designing effective mitigation strategies. Moreover, evaluation of the toxicity of airborne particulate matter is important since the classically measured particulate matter (PM) concentrations appear insufficient for characterizing the impact on human health. Oxidative potential (OP) measurement has recently been developed to quantify the capability of PM to induce an oxidative imbalance in the lungs. As a result, this measurement unit could be a better proxy than PM mass concentration to represent PM toxicity. In the present study, two source apportionment analyses were performed using positive matrix factorization (PMF) from organic aerosol (OA) mass spectra measured at a 15 min time resolution using a time-of-flight aerosol chemical speciation monitor (ToF-ACSM) and from 19 trace elements measured on an hourly basis using an online metal analyser (Xact 625i). The field measurements were carried out in summer 2018. While it is common to perform PMF studies individually on ACSMs and more recently on Xact datasets, here we used a two-step methodology leading to a complete PM₁ source apportionment. The outputs from both OA PMF and Xact PMF, the inorganic species concentrations from the ACSM, and the black carbon (BC) fractions (fossil fuel and wood burning) measured using an Aethalometer (AE33) were gathered into a single dataset and subjected to a combined PMF analysis. Overall, eight factors were identified, each of them corresponding to a more precise source than performing single PMF analyses. The results show that besides the high contribution of secondary ammonium sulfate (28 %) and organic nitrate (19 %), about 50 % of PM₁ originated from distinct combustion sources, including emissions from traffic, shipping, industrial activities, cooking, and biomass burning. Simultaneously, PM₁ filters were collected during the experimental period on a 4 h sampling basis. On these filters, two acellular OP assays were performed (dithiothreitol; OP^DTT and ascorbic acid; OP^AA) and an inversion method was applied on factors issued from all PMFs to assess the contribution of the PM sources to the OP. This work highlights the sensitivity of OP^AA to industrial and dust resuspension sources and those of OP^DTT to secondary ammonium sulfate, shipping, and biomass burning.

1 Introduction

Airborne particulate matter (PM) is a significant contributor to air pollution, leading to adverse effects on ecosystems, climate stability, and the environment (Beelen et al., 2014; Cohen, 2017; Jacob, 1999). Air quality is currently a major public health outcome and is responsible for 4.2 million premature deaths worldwide each year (WHO, 2021). More worryingly, this figure is expected to double by 2050 (Lelieveld et al., 2015). However, the links between air pollution and human health effects are not fully understood yet, but PM appears to be a key pollutant in aerosol toxicity (Medina et al., 2016; Zhang et al., 2016). Depending on the size and chemical composition of the particles, PM may cause significant damage to pulmonary cells due to the ability of particles to penetrate deeply into the organism and induce inflammatory responses in lung cells (Strak et al., 2012). There has been growing interest in recent years in submicron aerosol PM₁, which are deposited deeply in the respiratory system, reaching the alveoli of the lungs (Sturm, 2020). PM₁ are associated with the physicochemical processing of compounds also resulting from anthropogenic sources, mainly combustion sources, and are known to contribute to the health impacts of PM (Grigas et al., 2017; Manigrasso et al., 2020).

Although an increasing number of studies investigated the potential effects of submicron particles on the risk of respiratory diseases, this subject still represents a research line that needs further toxicological and epidemiological studies (Hu et al., 2022). It was recently shown that these particles are strongly linked to the occurrence of cardiovascular disease, perhaps due to the higher surface-to-volume ratio (Münzel et al., 2022). While the air quality guidelines regarding PM levels were recently updated for PM₁₀ and PM_2.5 (WHO, 2021), no regulation for PM₁ has yet been established, and several studies suggest their monitoring should be considered in this respect (Kumar et al., 2010).

While all mechanisms leading to the toxicity of airborne particles are not fully established to date, it is generally believed that the activity of reactive oxygen species (ROS) could play an important role. These chemical species that are carried or induced in the lung are suspected of disrupting the natural redox balance, causing oxidative stress, a key factor in the inflammatory response in the organism, ultimately leading to diseases such as asthma or chronic bronchitis (Abrams et al., 2017; Dellinger et al., 2001; Møller, 2014; Pope, 2004). Thus, new approaches have been investigated over the past two decades to better quantify oxidative stress and, ultimately, the impact that exposure to airborne particulate pollution has on the population. Indeed, the measurement the of oxidative potential (OP) of PM is seen as a new promising metric and probably a better representative proxy of health impacts than the PM mass concentration, as it quantifies the ability to generate ROS in vivo. OP measurements integrate many properties of PM that are important for the interactions with lung fluid, such as size distribution, specific surface area, or PM chemical composition (Ayres et al., 2008). Indeed, changes in the chemical composition of aerosol and consequently the type of sources from which they are emitted may lead to a significant difference in airborne particle toxicity (Boogaard et al., 2012). Furthermore, previous studies showed the value of fine-mode particles in understanding the oxidative effects of PM (Chen et al., 2017).

An active field of research for air quality is understanding the emission sources of PM through source apportionment techniques, using approaches such as statistical receptor models (e.g. chemical mass balance (CMB) or positive matrix factorization (PMF)) (Paatero and Tapper, 1994). PMF, which is widely used in the community, can be based on data from online analysers such as aerosol mass spectrometers (AMSs) and aerosol chemical speciation monitors (ACSMs), which determine the resolution of organic mass spectra with high sensitivity (Bozzetti et al., 2017a; Chen et al., 2022; Crippa et al., 2013b) or from off-line filters analyses (Borlaza et al., 2021, 2022).

Several methods have been explored to combined datasets from various online analysers to perform PMF. Some studies proposed the combination of high-resolution AMS and proton transfer reaction–mass spectrometry (PTR-MS) measurements to refine the links between particle-phase organics and their precursors (Crippa et al., 2013a; Slowik et al., 2010). Some others performed source apportionment by including both organic and inorganic fractions from the AMS (Äijälä et al., 2019; McGuire et al., 2014; Sun et al., 2012), improving the resolution of the factors and their chemical nature. More recently, Zografou et al. (2022) performed PMF analysis on a combined organic and inorganic year-long dataset from a ToF-ACSM. Tong et al. (2022) combined into a single dataset AMS and extractive electrospray ionization time-of-flight mass spectrometry (EESI-ToF) measurements providing an optimized identification and quantification of the organic factors, more specifically the SOA fraction. Nursanto et al. (2023) also tried a hybrid approach by combining organic aerosol concentrations from a ToF-ACSM with the particle size distribution from a scanning mobility particle sizer (SMPS) and resolved organic factors related to new particle formation and growth.

However, more efforts are needed to combine instrument datasets and apportion the sources of the total PM₁ fraction. A multi-time resolution approach was suggested by Via et al. (2023) mixing PM₁ data from online analysers and offline filters.

Belis et al. (2019) conducted three separated PMF analyses (on offline filter data, online organic data, and online inorganic species) which produced reference profiles to constrain a fourth PMF with combined online data. Petit et al. (2014) also followed a multi-step methodology (PMF²) comprising the use of results from the first PMF runs as inputs for their combined PMF using ACSM and Aethalometer (AE33) data. To our knowledge, no other study has performed this type of PMF² approach, in particular for all PM₁ components. A known drawback of performing PMF on OA mass spectra from ACSM or AMS is the resolution of the secondary organic aerosol (SOA) origin. SOA factors are usually reported as either a single factor or two factors separated by their degree of oxygenation rather than in terms of sources. A PMF² approach using previous OA factors combined with other species and/or PMF factors may enable a more accurate identification and quantification of the SOA fraction in the PM sources. The current study addresses this challenge by using the PMF² method for the PM₁ fraction measured with online analysers (i.e. ToF-ACSM, Xact 625i, and AE33) at a high time resolution (< 1 h).

While the study of the relationship between the OP and single compounds is a step forward in understanding the chemical mechanisms involved, it seems unrealistic to characterize the OP of the aerosol as a whole in this way, given the myriad of chemical elements that make up the aerosol. A more integrative approach is to consider OP sources since aerosol is a complex mixture where synergistic and antagonistic effects can occur (Yu et al., 2018). This enables an estimation of intrinsic OP of sources but also of their contribution to a population's health-related exposure, which is more practical for adapting air quality management policies (Borlaza et al., 2018; Calas et al., 2019; Fang et al., 2016; Grange et al., 2022; Verma et al., 2014; Weber et al., 2018; Yu et al., 2019).

Major insights have been gained from OP measurements on filters. Many studies, based on filters collected at different sites, have established associations between OP and PM sources, thus improving our understanding of the health impact of specific types of airborne particles (Borlaza et al., 2022; Daellenbach et al., 2020; Fang et al., 2016; Weber et al., 2021). However, the next step is to develop online OP measurements to understand one-time extreme events and ultimately enable near-real-time (NRT) OP acquisition for better air quality management.

To make progress in this direction, here we propose to investigate the submicron aerosol sources – highly relevant in health impacts (Lin et al., 2020) – contributing most to two types of OP (ascorbic acid assay (OP^AA) and dithiothreitol assay (OP^DTT)) in a challenging environment, the Marseille area. Effectively, the activity of the industrial port complex and several industrial areas close to this highly urbanized city, combined with specific meteorological conditions (wind regimes and significant photochemistry in summer), involves frequent pollution episodes (Chazeau et al., 2022, 2021; Salameh et al., 2018). A challenging double-PMF method using chemical online analysers (ToF-ACSM, Xact, and the Aethalometer) is proposed as a first step. High-frequency acquisition of OP observations with filter sampling every 4 h for 7 weeks was then coupled to these data. On this basis, an OP source apportionment using a multilinear regression approach is provided. This method enables an estimation of the oxidizing capacity of each µg of PM₁ from the emission sources identified but also the relative contribution of each source to OP^AA and OP^DTT on a 4 h basis. The complementarity of these assays gives a broader and more representative view of the PM₁ health impact.

2 Material and methods

2.1 Site and sampling

2.1.1 Site

The Marseille-Longchamp supersite is an urban background site of Marseille – the second most populated city in France (about 870 000 inhabitants, with a density of 3600 inhabitants per km² in 2019). The site is located in the heart of the Parc Longchamp, in the fourth district of Marseille (43°${\mathrm{18}}^{\prime }{\mathrm{20}}^{\prime \prime }$ N; 5°23^′41; m a.s.l.). Figure 1 shows the site location in the city and the main areas in the surroundings with the maritime port of Marseille within 2 km and industrial areas within 30 km with petroleum refining, coke production, and steel facility activities (Salameh et al., 2018). In addition to shipping and local industrial sources, Marseille suffers from the second largest traffic congestions in France, generating a year-round source of background traffic (Chazeau et al., 2021). The location of the city also leads to influences of natural and biogenic emissions such as marine aerosol, terrestrial vegetation, Saharan dust, or crustal dust. Finally, secondary organic aerosol (SOA) formation events and high ozone concentrations formed by intense photochemistry are frequent during warm periods in the area (El Haddad et al., 2013; Flaounas et al., 2009).

2.1.2 Sampling campaign

The sampling campaign and site have been detailed elsewhere (Chazeau et al., 2022). Briefly, the field work took place during the summer of 2018 over a 7-week period (from 11 July to 1 September 2018). The sampling station was equipped with a range of analytical instruments for characterization of submicron aerosols: a time-of-flight aerosol chemical speciation monitor (ToF-ACSM; Fröhlich et al., 2013) to measure in near real time (10 min resolution) the chemical composition of non-refractory PM₁ (organic aerosol, NH${}_{\mathrm{4}}^{+}$, NO${}_{\mathrm{3}}^{-}$, SO${}_{\mathrm{4}}^{\mathrm{2}-}$, and Cl⁻); the dual-spot Magee Scientific AE33 Aethalometer (Drinovec et al., 2015) equipped with a PM_2.5 cut-off inlet to measure the equivalent black carbon concentrations (BC: with a distinction between BC_FF and BC_WB origins) at a 1 min resolution; and a Xact625i (Cooper Environmental) to measure a user-defined list of 25 PM₁ trace elements with a time resolution of 60 min. PM₁, PM_2.5, and PM₁₀ mass concentrations were determined with an optical particle counter (FIDAS 200; PALAS). A 3D sonic anemometer for temperature, wind direction, and velocity measurements and O₃, NO_x, and SO₂ analysers are also amongst the permanent instruments of the station. Finally, PM₁ collection for OP analysis was performed for 15 d (from 11 July and 25 July 2018) every 4 h on 150 mm diameter quartz fibre filters (Whatman Tissuquartz; pre-heated at 500 °C for 8 h), using a high-volume aerosol sampler (HiVol, Digitel DA80) at a flow rate of 30 m³ h⁻¹. A total of 90 samples and four blank filters were collected with a time resolution of 4 h. Procedural care for filter handling, sampling, and storage was taken to avoid contamination and modification of the filter deposit after sampling (Weber et al., 2018).

All the instruments ran acquisitions during the whole campaign period (i.e. 11 July to 1 September 2018). However, it should be noted that ToF-ACSM stopped between 12 July 2018 19:00 and 13 July 2018 03:00 (UTC) and between 14 July 2018 15:00 and 15 July 2018 03:00 (UTC), which means the number of filter samples decreased to 83 instead of 90 for the OP deconvolution model discussed below.

Figure 1 Location of the Marseille-Longchamp supersite and localization of the main industrial areas around Marseille, France (© PlanetObserver, Geoportail/© Google Maps).

2.2 OP analysis

OP was assessed by the depletion rate of anti-oxidant compounds, using two methods with dithiothreitol (DTT) and ascorbic acid (AA) (Calas et al., 2019, 2018, 2017). DTT depletion in contact with PM extracts was determined by dosing the remaining amount of DTT with dithionitrobenzoic acid (DTNB) at different reaction times (0, 15, and 30 min), and absorbency was measured at 412 nm using a plate spectrophotometer (Tecan, M200 Infinite). Briefly, 25 µM of DTT and phosphate buffer react for 30 min with PM suspensions at 25 µg mL⁻¹ in a simulated lining fluid. The AA assay is a simplified version of the synthetic respiratory tract lining fluid (RTLF) assay (Kelly, 2003), where only AA is used. A mix of 80 µL of PM suspension with 24 nmol of AA (100 µL of 0.24 mM AA solution in Milli-Q water) is used, and AA depletion is read continuously for 30 min by absorbency at 265 nm (TECAN, M1000 Infinite). The depletion rate of AA was determined by linear regression of the linear section data. For both assays, a 96-well plate was auto-shaken for 3 s before each measurement and kept at physiological conditions (37.4 °C).

Three filter blanks (laboratory blanks) and three positive controls (1,4-naphthoquinone, 24.7 µmol L⁻¹) were included in each plate (AA and DTT). The average values of these blanks were then subtracted from the sample measurements of the given plate. The detection limit value was defined as 3 times the standard deviation of laboratory blank measurements (laboratory blank filters in Gamble + DPPC solution). Three replicates are performed with two absorbance measurements at each time. The short duration of the sampling time (4 h), the type of PM fraction (PM₁), and consequently the low mass recovered on each filter resulted in some OP replicate measurements being below the detection limit, suggesting the highest uncertainties in these results as is usual for such measurements (Calas et al., 2018; Weber et al., 2021). A propagation of these uncertainties was carried out (on average 15 % of the OP measurement) to deal with thus issue. Hereafter, the OP normalized in volume relative to the AA assay and the DTT assay is denoted OP${}_{\mathrm{v}}^{\mathrm{AA}}$ and OP${}_{\mathrm{v}}^{\mathrm{DTT}}$, respectively.

2.3 Source apportionment using positive matrix factorization

Source apportionment was performed through the positive matrix factorization (PMF; Paatero and Tapper, 1994) method using the multi-linear engine (ME-2) solver (Paatero, 1999) and run within the Source Finder Professional (SoFi Pro) software (Datalystica Ltd., Villigen, Switzerland; Canonaco et al., 2021, 2013). PMF is a bilinear unmixing model widely used to determine the atmospheric aerosol sources based on online measurements (Canonaco et al., 2021; Chazeau et al., 2022; Chen et al., 2022). The equation is described as follows:

$$\begin{array}{}\text{(1)}& x_{i,j}={\sum }_{k=\mathrm{1}}^n{g}_{i,k}.f_{k,j}+e_{i,j},\end{array}$$

with x_i,j a non-negative matrix of measurements which is factorized into g_i,k the factor time series, f_k,j the factor profiles, and e_i,j the model residuals. The indexes i, j, k, and n are the time, variables, discrete factor numbers, and total number of factors in a solution, respectively.

ME-2 solves the model solution by using a least squares algorithm to iteratively minimize the following object function Q defined as the sum of the squared model residuals weighted by their respective uncertainties (σ_ij):

$$\begin{array}{}\text{(2)}& Q={\sum }_i{\sum }_j\left({\displaystyle \frac{e_{ij}}{{\mathit{\sigma }}_{ij}}}\right){.}^{\mathrm{2}}\end{array}$$

Generally, the PMF model does not result in a mathematically unique solution as a multiple combination of f_k,j and g_i,k may provide a similar value of Q. ME-2 allows us to introduce a priori information in the model using known source profiles or time series to orient solutions towards environmentally meaningful rotations (Paatero, 1999; Paatero and Hopke, 2009). An advantage of SoFi Pro is that it makes it possible to control the rotational ambiguity of the solution by applying constraints with the a-value approach (Canonaco et al., 2013):

$$\begin{array}{}\text{(3)}& f_{kj}=f{{}^{\prime }}_{kj}\pm a\cdot f{{}^{\prime }}_{kj}\end{array}$$

$$\begin{array}{}\text{(4)}& g_{ik}={g^{\prime }}_{ik}\pm a\cdot {g^{\prime }}_{ik},\end{array}$$

where the scalar a defines the range (between 0 and 1) to which $f{{}^{\prime }}_{kj}$ and ${g^{\prime }}_{ik}$ can vary from the known input profile (f_kj) or time series (g_ik).

In the current study, source apportionment was conducted using the PMF model on three distinct datasets. In the first two analyses, PMF was applied separately to the OA dataset (PMF_organics) and the Xact dataset (PMF_metals). As a second step, the outputs from the previous source apportionment analyses, the inorganic species concentrations from the ToF-ACSM measurements (NO${}_{\mathrm{3}}^{-}$, NH${}_{\mathrm{4}}^{+}$, SO${}_{\mathrm{4}}^{\mathrm{2}-}$, and Cl⁻), and the deconvolved BC source (BC_FF and BC_WB) concentrations from the AE33 were combined into a single dataset to perform a total PM₁ source apportionment (PMF${}_{{\mathrm{PM}}_{\mathrm{1}}})$.

2.4 PMF preparation and optimization

2.4.1 OA dataset

ACSM data were acquired with the Igor-based Acquility v2.1.4 software and analysed with Tofware v3.2 also developed in Igor Pro (Wave Metric inc., Lake Oswego, Oregon, USA). The data treatment including ionization efficiency (IE) and relative ionization efficiency (RIE) calibrations, collection efficiency (CE) correction, and detection limit determination is detailed in Chazeau et al. (2021). The complete PMF methodology and optimization applied to OA mass spectra at MRS-LCP are described in a dedicated paper (Chazeau et al., 2022) and we will provide only a brief summary here. The model is performed on data input including 185 variables from $m/z$ 12 to $m/z$ 214 and 4422 time steps (15 min intervals) from 1 July to 1 September 2018. The error matrix was exported from the Tofware software and the calculation included ion counting statistics, background errors, electronic noise, and a minimum error from the measurement of a single ion.

A five-factor solution was resolved with three constrained and two unconstrained factors. Both hydrocarbon-like organic aerosol (HOA) and cooking-like organic aerosol (COA) factor profiles were constrained using the reference profiles from Ng et al. (2011) and Crippa et al. (2013b), respectively. Shipping and industrial organic aerosol (Sh-IndOA) factor time series were constrained with the SO₂ concentrations as this is a specific proxy for these emissions in the Marseille area (El Haddad et al., 2013). The a-value ranges for these constraints were optimized based on previous sensitivity analyses (Chazeau et al., 2022) and random a values of 0–0.6, 0–0.2 and 0–0.2 for HOA, COA, and Sh-IndOA, respectively, were retained. The two remaining factors corresponded to the secondary oxygenated organic fraction separated into two components: a less oxidized organic aerosol (LOOA) factor and a more oxidized organic aerosol (MOOA) factor.

In order to explore the rotational ambiguity and statistical uncertainties of the PMF solution, a bootstrap resampling strategy is applied (Efron, 1979), where 100 repeated runs are performed to test the stability of the solution. The inspection of the 100 generated runs was achieved based on a selection of predefined criteria customized within SoFi Pro (Canonaco et al., 2021; Chazeau et al., 2022; Chen et al., 2022, 2021). First, three criteria are defined as acceptance thresholds to evaluate the quality of the PMF runs. The r Pearson correlation (denoted as r hereafter) with BC_FF for HOA, the ratio between lunch hours (11:00 and 12:00 UTC) and the average background hours in the morning (06:00–08:00 UTC) for COA, and the r Pearson correlation with SO₂ for Sh-IndOA were used. Then, the monitoring of f₄₃ intensity for LOOA and _f44 intensity for MOOA are used as repositioning criteria to avoid mixing of the unconstrained factors since they are not sorted automatically among the different PMF iterations. The bootstrapped runs fulfilling the criteria list were then averaged into a unique PMF solution.

2.4.2 Xact dataset

The PMF_metals was performed on the Xact data matrix of hourly element concentrations from 11 July to 1 September 2018. First, some individual species were excluded to improve the quality of the analysis. Elements whose concentrations were below their respective minimum detection limit (MDL) more than 90 % of the time were not included in the inputs (Fig. S1 in the Supplement). The MDLs were provided by the manufacturer and are given in Table S1 in the Supplement. From this approach, the following 17 elements remained: As, Br, Ca, Cd, Cr, Cu, Fe, K, Mn, Ni, Pb, Sb, Se, Sn, Ti, V, and Zn. Co was also included in the inputs as it showed a good correlation with the Ni element (R²=0.5).

A PMF error matrix (σ_i, j) was estimated using Eq. 5 for concentrations greater than the MDL (Reff et al., 2007; Ryder et al., 2020):

$$\begin{array}{}\text{(5)}& {\mathit{\sigma }}_{i,j}=\sqrt{{\mathrm{MDL}}_i^{\mathrm{2}}+u_{i,j}^{\mathrm{2}}},\end{array}$$

where u_i,j is a specific analytical uncertainty for each data point provided by the Xact software. It includes both the uncertainties of the sampling air volume and the uncertainties of the mass spectra deconvolution calculated by the Xact software. For the concentrations below the MDL, the values were replaced by the MDL of metal divided by 2. The corresponding error is often set to 5/6 × MDL (Polissar et al., 1998). Since the relative error for each data point was mostly less than 50 %, Polissar et al. (1998) recommended applying a relative error between 100 % and 250 % for values below the MDLs. In our dataset, some relative errors for data greater than the MDLs were much larger than 50 %. Following Polissar et al. (2001), the elements were down-weighted by using larger error estimates for values below the MDLs. The methodology and the PMF tests panel applied are described in the Supplement Sect. S1. All the data points with a signal-to-noise ratio (S2N) below 1 were down-weighted by adding a penalty function of 1/S2N to the error (Rai et al., 2020; Visser et al., 2015a). The weighting is performed cell-wise as some variables had low average S2N but some high S2N periods.

Intense firework episodes were recorded during the French National Day celebration on 14 July. While it is common to exclude such episodes from the PMF analyses to reduce modelling uncertainties linked to very high concentrations (Ducret-Stich et al., 2013), some studies succeeded in identifying a firework factor profile. Rai et al. (2020) performed a constrained PMF analysis only on the firework hours and identified a firework factor based on the $K/S$ elemental concentration ratio in black powder. The factor profile of fireworks was then constrained in the final complete dataset PMF analysis. Manousakas et al. (2022) ran a PMF on the entire dataset constraining all the source profiles except the firework factor and the time series of all sources. The firework time series were set to 0 except during the corresponding events. In this study, we followed the same logic and the firework hours (13 July 20:00 UTC to 14 July 05:00 UTC; 14 July 20:00 UTC to 15 July 14:00 UTC) were removed from the dataset to inspect the remaining sources. The PMF inputs without the firework points (WFP) consist of 1201 time points with 1 h step and 18 elements. As a second step, PMF analyses were conducted only on firework day (13 and 14 July) points (FDP) to determine a specific profile. The dataset represented 60 time points and 19 variables. In addition to the elements previously selected, Bi was included as it was exclusively associated with the firework events. Bi in the form of bismuth trioxide (Bi₂O₃) is commonly used instead of the toxic lead forms for crackling fireworks, the so-called dragon's eggs (Mohan, 2010; Perrino et al., 2011). Finally, PMF was performed on the total dataset (1230 time points and 19 variables) by adding a constrained firework factor.

One important step is selecting the number of factors based on both mathematical diagnostics and the environmental meaning of the factors. Solutions with a range of one to eight factors were examined for the WFP dataset. The selection is made based on the changes in $Q/Q_{\exp }$ ($\mathrm{\Delta }Q/Q_{\exp })$ and on real and noisy unexplained variation (ΔUEV_real, ΔUEV_noisy) when increasing the number of factors (Fig. S2). A large reduction in these values would indicate limited improvements of the model residuals and explained variability. There were no significant changes in ΔUEV_real and ΔUEV_noisy between two and eight factors. However, $\mathrm{\Delta }Q/Q_{\exp }$ showed a decrease of up to five factors, meaning the changes in $Q/Q_{\exp }$ were very low. To relate the factors from the PMF to specific sources, the diurnal trends, the time series of the elements, and a comparison with some external tracers were examined. We could clearly identify five environmentally relevant factors: dust resuspension, shipping, industrial, tire and/or brake wear, and regional background factors. Selecting the six-factor solution results in an unresolved Br-rich factor which cannot be attributed to a specific source or aerosol processes. Therefore, the five-factor solution was chosen as the best representation of the data.

For the FDP dataset, the PMF analysis resolved the same five factors in addition to a firework factor. However, the analysis showed some mixing between the regional background and the firework factors due to the large contribution of K in both factors. To avoid this mixing, the regional background profile was tightly constrained using the profile resolved with the WFP dataset and an a value of 0.1. The runs were repeated 50 times with a bootstrap resampling strategy to test the stability of the solution. This time a well-defined firework factor was resolved (Fig. S3) with an elemental composition in agreement with other studies (see Sect. 3.2.2).

The averaged firework factor profile retrieved from the FDP dataset was used as a constraint for the complete dataset with a values randomly initialized between 0 and 0.5 with an increment of 0.1. This initialization is used to evaluate whether larger deviations could improve the results (Canonaco et al., 2021). The factor time series were also constrained with similar a values and set to 0 except during the firework events. The remaining factors were left unconstrained in the complete dataset solution.

Similarly to the PMF_organics methodology, 100 bootstrapped runs are conducted and a criteria-based selection is used to assess the quality and position of the PMF runs (Fig. S4). This statistical selection is described in the Supplement Sect. S2.

2.4.3 Combined PM₁ dataset

Following the methodology described by Petit et al. (2014), we combined PMF outputs with the remaining chemical species of PM₁. Thus, the factors from both PMF_organics and PMF_metals were combined with BC_FF, BC_WB, NO${}_{\mathrm{3}}^{-}$, NH${}_{\mathrm{4}}^{+}$, SO${}_{\mathrm{4}}^{\mathrm{2}-}$, and Cl⁻ concentrations. BC_WB and BC_FF were deconvolved based on the model of Sandradewi et al. (2008). We used the 470 and 950 nm wavelengths with a constant absorption Angström exponent of 1.68 and 1.02 for pure wood burning and traffic, respectively, as recommended by Zotter et al. (2017) and Chazeau et al. (2021).

The uncertainty matrix was constructed as follows: uncertainties for the ToF-ACSM inorganic species were exported with Tofware similarly to those for the organics; BC_FF and BC_WB uncertainties were estimated based on Eq. 5 and on Eq. S1 in the Supplement, with MDLs set to 0.1 for both species and u_i,j, the relative uncertainties, set to 40 % (Petit et al. 2014) multiplied by the species concentrations; PMF_organics and PMF_metals outputs errors were taken from the standard deviations of the factor time series from the bootstrap analyses, which can be used as statistical uncertainties (Canonaco et al., 2021).

All the variables were synchronized to a 1 h time resolution corresponding to the 4 h time resolution of the filter sampling, and the missing measurement periods of each instrument were removed from the PMF analysis. The firework_metals factor was excluded here since the ToF-ACSM was not running during the main active period of this event (i.e. 14–15 July). Finally, this leads to PMF being performed over 16 variables and 849 time steps.

Since we combined data from three instruments presenting different measurement uncertainty calculations, S2N ratio, and relative numbers of variables, it is necessary to ensure good representativeness of each group in the PMF analysis (Tong et al., 2022). Some studies suggested applying a relative instrument weight to balance the scaled residuals of each subgroup of data (Crippa et al., 2013a; Slowik et al., 2010; Tong et al., 2022; Via et al., 2023). The details about this instrument weighting process are given in the Supplement Sect. S3. As described in Sect. 2.4.2, a cell-wise down-weighting was applied to data points with a weak S2N ratio. PMF runs were performed from 1 to 12 factors to inspect and identify the most physically meaningful factors. Eight factors were clearly determined: biomass burning, cooking, industrial, dust resuspension, traffic, organic nitrate-rich (ON-rich), shipping, and ammonium sulfate-rich (AS-rich) factors. However, conducting different seed runs showed a high degree of rotational ambiguity in the solution, with some unstable factors which cannot be resolved systematically (Table S4).

While it is common to set some variables to 0 in the factor profiles based on prior chemical knowledge of the sources (Bozzetti et al., 2017b; Weber et al., 2019), this method did not enable a clear separation of factors. Here, the entire profiles were constrained for the biomass burning, cooking, and industrial factors using their profiles from the most interpretable solutions as the base case. Profile constraints were applied with an a value of 0.4, 0.1, and 0.05 for biomass burning, cooking, and industrial factors, respectively, leading to a more stable solution.

Further discussions on the factor identification, the rotational ambiguity, the a-value selection for the constrained profiles, and the acceptance criteria are provided in the Supplement Sect. S4. Similarly to the two previous PMF analyses, a bootstrap analysis was conducted over 100 runs and all the accepted runs were averaged into the reported solution.

2.5 OP apportionment

An inversion method is applied on factors issued from all PMFs to assess contributions of the PM sources to the OP. The dependent variable OP expressed in nmol min⁻¹ m⁻³ is explained by a linear combination of the mass contribution of PM sources (µg m⁻³) taken as independent variables as follows:

$$\begin{array}{}\text{(6)}& \mathit{O}\mathit{P}=\mathbf{H}\times {\mathit{\beta }}_n+\mathit{\epsilon },\end{array}$$

where the OP vector (p×1) is the measured OP (p observations) with a 4 h time resolution, H is a matrix (nx(p+1)) of n sources resolved from the PMF analyses plus the intercept (data of the source PMF with a time step of 1 h have been averaged over the 4 h time step of the OP data), and the ε vector (p×1) accounts for the misfit between the observations and the model. The regression coefficient β provided by the model (Eq. 6) is interpreted as an intrinsic OP of the n sources (nmol min⁻¹ µg⁻¹). Basically, it expresses how much the OP would increase if we increase 1 µg m³ of the given source. The source-specific OP contribution is calculated by multiplying the regression coefficient β of each source by the respective mass contribution of the source to PM. This methodology is essentially based on previous works detailed in Borlaza et al. (2022) and Weber et al. (2021, 2018).

Three scenarios in the construction of the matrix H (Eq. 6), i.e. the contribution of various source factors of PM identified by each of the three PMFs, have been considered to make the best use of the results from the different PMFs. Three models were tested in each scenario (i.e. nine solutions): weighted least squares linear regression (WLS), weighted robust multiple linear regression with an iterative M-estimator, and partial least square regression (PLS). The description of the three scenarios (Eqs. S4, S5, and S6) and the three models (Eqs. S7, S8, and S9) are shown in the Supplement. The best model solutions with the lowest RSR (RMSE standard deviation ratio) are presented: M-estimator in scenarios 1 and 2 (see Table S3) and WLS in scenario 3. Source factor contributions exhibiting a Pearson's association with OP of less than 0.1 were discarded from the predictor variables. Finally, to provide robust estimates of the model output terms in scenario 3, the process was performed 500 times with bootstrapped inputs by following the method of Canonaco et al. (2021). Briefly, 15 % of the input samples were randomly removed or duplicated before each run. Output runs with $R_{\mathrm{adjusted}}^{\mathrm{2}}$ values below the value found in the first inversion model without bootstrapping (Eq. 6) (i.e. $R_{\mathrm{adjusted}}^{\mathrm{2}}<\mathrm{0.3}$ for OP^DTT and $R_{\mathrm{adjusted}}^{\mathrm{2}}<\mathrm{0.4}$ for OP^AA in scenario 3) were removed. Factors with a variance inflation factor (VIF) of > 10 were also removed, because this suggests strong collinearity between them (Calas et al., 2019). Robust linear regression was performed with the MASS package developed in R (Grange et al., 2022; Venables and Ripley, 1997) and weighted linear regression was performed using the stats package developed in R.

3 Results

3.1 OP results

Meteorological conditions typical of those prevailing during summertime in the region occurred during the period of the study, with land and sea breeze cycles (except between 20 and 24 July), associated with stable conditions characterized by ozone episodes (six regulatory exceedances, with a maximum at 166 µg m⁻³). NO_x, O₃, and SO₂ average concentrations during the OP apportionment period were, respectively, 20, 80, and 2.5 µg m⁻³ (Fig. S5). The overall period (n=83 samples) was characterized by an average PM₁ concentration of 13.2±3.4 µg m⁻³.

This study is the first to characterize OP in France with a 4 h time step, providing an overview of the daily OP_v variation. Figure 2a shows the difference between night and day and Fig. 2b presents the typical daily OP_v variation (without firework episodes), associated with organic fractions of aerosol quantified by ToF-ACSM, metallic fractions of aerosol quantified by Xact, and PM₁ variations. The averages were calculated using 15 d during the period. The mass of PM₁, the metallic elements of organic aerosol, and the OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ are quite higher during the day than at night, while there is no significant variation in OP${}_{\mathrm{v}}^{\mathrm{AA}}$ between night and day. A t test demonstrates no significant difference (p<0.05) between the OP_v measured during the day (07:00–23:00 UTC) and the OP_v measured at night (23:00–07:00 UTC), for both OP assays. Two ANOVA variance tests were separately performed on OP^AA and OP^DTT 4 h series, and the result showed no significant difference between the two assays.

Figure S6b presents the composition in major chemical components of PM₁ measured by ToF-ACSM, Xact, and AE33 (organic fraction, metallic fraction, NH${}_{\mathrm{4}}^{+}$, Cl⁻, NO${}_{\mathrm{3}}^{-}$, SO${}_{\mathrm{4}}^{\mathrm{2}-}$, BC_FF, and BC_WB), together with the comparison of the reconstructed mass with these chemical components and the PM₁ concentration measured with the FIDAS. The Spearman association between these two time series is r_s=0.47, p<0.001. Figure S6a shows periods when the reconstructed mass fits well with the mass provided by the FIDAS, and periods when the reconstructed mass is overestimated and/or underestimated.

OP^AA and OP^DTT median values are, respectively, 0.62 and 1.47 nmol min⁻¹ m⁻³. Figure 3 shows variations in both OP assays set against PM₁ mass. These OP values are characteristic of the coastal environment in warm periods (Calas et al., 2019), but we note that they are rather low compared to many other series (for PM₁₀ or PM_2.5) measured with the same methods in other environments (Weber et al., 2018).

Spearman coefficients (r_s) between PM₁ mass measured by FIDAS and OP display some differences (r_s PM₁ vs. OP${}_{\mathrm{v}}^{\mathrm{AA}}=\mathrm{0.23}$ (p<0.01) and r_s PM₁ vs. OP${}_{\mathrm{v}}^{\mathrm{DTT}}=\mathrm{0.63}$ (p<0.001)) where PM₁ has a greater association with OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ than with OP${}_{\mathrm{v}}^{\mathrm{AA}}$. These Spearman coefficients are close to those found by in't Veld et al. (2023) on PM₁ all year long in a similar urban coastal environment (Barcelona) (r_s PM₁ vs. OP${}_{\mathrm{v}}^{\mathrm{AA}}=\mathrm{0.29}$ (p<0.001) and r_s PM₁ vs. OP${}_{\mathrm{v}}^{\mathrm{DTT}}=\mathrm{0.73}$ (p<0.001)). The higher association between OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ and PM₁ compared to OP${}_{\mathrm{v}}^{\mathrm{AA}}$ and PM₁ has already been observed in other studies conducted on PM₁₀ (Calas et al., 2019; Weber et al., 2021; Janssen et al., 2014). This phenomenon is attributed to the heightened sensitivity of AA to chemical composition, exhibiting robust specificity. Moreover, DTT demonstrates superior sensitivity to aerosol concentration owing to its more balanced sensitivity to chemical constituents (Gao et al., 2020a).

Figure 2 Profiles of OP${}_{\mathrm{v}}^{\mathrm{AA}}$, OP${}_{\mathrm{v}}^{\mathrm{DTT}}$, and organic and metal fractions of submicron aerosol and PM₁ during (a) night (23:00–07:00 UTC) and day (07:00–23:00 UTC) and during (b) different times of the day following the time step of OP (4 h).

Figure 3 (a) OP${}_{\mathrm{v}}^{\mathrm{AA}}$, OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ and PM₁ mass time series over the campaign at the Marseille-Longchamp site. (b) Distribution of OP${}_{\mathrm{v}}^{\mathrm{AA}}$ and OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ in mean and median values.

3.2 Interpretation of PMF factors

3.2.1 PMF analysis of OA

The PMF_organics results are presented in the Supplement section with the profiles of the five factors (Fig. S7a), the time series (Fig. S7b), the relative fractions (Fig. S7c), and the diurnal cycles (Fig. S7d). As expected in this period, the secondary fraction including both LOOA and MOOA contributed to the highest part of OA with 33.4 % each. HOA represented 16.3 % of the total OA, followed by COA (13.7 %) and Sh-IndOA (3.2 %). HOA is assumed to be related to traffic exhaust emissions and displayed a distinct bimodal pattern with significant peaks during the morning and evening rush hours. COA also showed a bimodal pattern with increased concentrations during lunchtime and in the evening. Sh-IndOA concentrations accounted for combined plumes from the industrial area of Fos-sur-Mer and from the shipping activity of Marseille harbour. These emissions are advected on-site by sea breeze in the morning and slowly decrease throughout the day. LOOA and MOOA are distinguished based on their ${}_{f\mathrm{44}}/f_{\mathrm{43}}$ ratio, with the higher ratio for the more oxidized part. According to the diurnal cycles, MOOA showed a flat pattern suggesting a long-range transported origin and a slight increase at mid-day potentially attributed to more local photochemical activity. By contrast, LOOA concentrations were higher at night related to some night-time chemistry. Further descriptions of OA sources in MRS-LCP are provided in Chazeau et al. (2022).

3.2.2 PMF analysis of metals

The PMF_metals solution is investigated with the factor profiles and time series presented in Fig. 4, along with the factor-relative diurnal cycles and contributions shown in Fig. S8. Non-parametric wind regressions (NWR) were also performed to determine the source concentrations attributed to their geographical origins using the wind direction and velocity (Henry et al., 2009; Petit et al., 2017). The results are displayed in Fig. S9.

Fireworks: This factor was only resolved for a short time from 13 to 15 July (see Sect. 2.4.2). Over its activity period, the factor represented a major fraction of the total element mass (up to 80 %) and contributed to K (73 %), Bi (100 %), Ti (71 %), and Cu (68 %). These elements are usually found in the composition of fireworks (Manousakas et al., 2022; Perrino et al., 2011; Rai et al., 2020; Vecchi et al., 2008). K is both a component of gunpowder (Drewnick et al., 2006) and is used as an oxidizer for firework bangs, while Cu and Ti produce blue and white flames when ignited.

Figure 4 (a) Average factor profiles with the coloured sticks indicating the normalized contribution of the element to the factor (left axis) and markers showing the normalized factor contribution to each element (right axis). Error bars represent the standard deviation of each profile. (b) Time series for the six factors resolved by the PMF_metals analysis.

Dust resuspension: This factor accounted for the largest contribution to the total elemental composition (53.7 %). Its profile presented the main fraction to Ca (98 %) and significant contributions to Ti (76 %), Cu (51 %), and Fe (48 %). These elements are major constituents of crustal soils and can be considered here to be urban dust (Almeida et al., 2020; Rai et al., 2021). Ca is also a compound often used for construction materials (Manousakas et al., 2022). The influence of construction work is supported by the factor diurnal profile which displayed increasing concentrations at 07:00 UTC followed by stable levels during the day and very low concentrations at night. The NWR plot showed a clear geographical origin from the northwest associated with a strong wind. Dust resuspension might occur under the mistral conditions (a regional wind for the Rhone Valley) with enhanced concentrations under high-velocity wind and dry conditions (Fig. S11). The factor could also include resuspension from non-exhaust traffic emissions due to the significant contribution of Cu.

Tire and/or brake wear: This factor showed high contributions to Zn (76 %) and Sb (29 %). The factor diurnal profile did not exhibit a distinct pattern and the time series were often affected by some spikes (Fig. 4b), which can be attributed to local short events. The factor time series showed a moderate correlation with exhaust traffic proxies such as BC_FF (r=0.33). However, Zn and Sb elements are known tracers of the brake wear emissions (Grigoratos and Martini, 2015; and references therein). While Zn is one of the most abundant metals in the brake linings, Sb is contained in the form of stibnite (Sb2S3) and employed as a lubricant to reduce vibration (Roubicek et al., 2008). Moreover, Zn is considered as a marker of tire wear particles (Panko et al., 2018). It should be noted that the tire and/or brake wear factor contributes to a low extent to the total metal composition (5.3 %). In a previous study in Zurich, Bukowiecki et al. (2009) showed that the contribution to non-exhaust emissions from both light and heavy-duty vehicles was very low in the submicrometre mode. Visser et al. (2015b) demonstrated that elements usually assigned to brake lining and tire wear emissions (e.g. Cu, Sb, Fe or Sn) are mainly found in the coarse mode at the “Marylebone Road” kerbside site, and Hays et al. (2011) reported similar trends for a near-highway site in Raleigh, with Zn being the only element significantly present in the fine mode. Such results suggest the existence of significant alternative source for these elements, potentially mixed in the regional-scale background factor.

Regional background: This factor contributed to a large range of elements: K (71 %), Cd (78 %), Br (74 %), Sn (74 %) and to some extent to As (58 %), Sb (56 %), Se (54 %), and Pb (43 %). The factor diurnal profile was mostly flat suggesting long-range transport of aged background compounds. Since most trace elements in the fine mode are non-volatile, they can undergo long-range atmospheric transport (Morawska and Zhang, 2002). This is supported by a strong correlation with the MOOA factor (r=0.6) resolved during the PMF_organics analysis. Furthermore, the NWR analysis displayed a regional geographical origin, with enhanced concentrations from the southerly sector (the Mediterranean Sea) and from the northeast sector with land breeze advecting aged air masses back to the site.

Shipping and industrial: The shipping factor accounted for the main fractions of V (97 %), Ni (88 %), and Co (43 %). The $\mathrm{V}/\mathrm{Ni}$ ratio has often been suggested as a proxy of heavy fuel combustion (Pandolfi et al., 2011; Viana et al., 2014). Here, we found a ratio of ∼ 2, which is in agreement with the typical range for shipping emissions (between 2 and 4) and with a ratio found in a previous study in Marseille (2.35; Salameh et al., 2018). The industrial factor contributed to Fe (47 %), Mn (63 %), Pb (38 %), and Se (35 %). This factor profile showed similarities with profiles from several industrial areas (Fig. S15). The contributions of some major elements (i.e. Fe, Ca, Mn, As, Zn) to this factor were in the same range as those of two iron converter zones and two storage zones measured by ICP-MS for the PM_2.5 fraction (Sylvestre et al., 2017). Combining these two factors (shipping and industrial factors) showed a strong correlation (r=0.74) and a similar diurnal pattern with Sh-IndOA (Fig. S16), which accounted for both the industrial emissions from Fos-sur-Mer and the shipping activity from the harbour (Chazeau et al., 2022). The shipping diurnal profile in Fig. S8b displayed a bimodal pattern which is linked to the diurnal trend of the ship departures and arrivals at the harbour (Chazeau et al., 2021). Once the sea breeze sets in, the first peak related to the ship arrivals is observed, followed by a second peak at 17:00 UTC due to ship departures. The diurnal profile of the industrial factor also exhibits increasing concentrations once the sea breeze occurs, which then gradually decreases throughout the day. The NWR analyses support the notion that the industrial and shipping factors were advected on-site by breeze from the Mediterranean Sea as they revealed clear hotspots from the southwesterly sector. The shipping factor showed high concentrations at lower wind speed compared with the industrial factor, highlighting more local emissions that were expected due to the proximity of the harbour.

3.2.3 PMF analysis of PM₁

The full PM₁ source apportionment solution is explored in this section with the average factor profiles (Fig. 5a), the time series (Fig. 5b), the pie chart of mass contributions (Fig. 5c), the average diurnal profiles (Fig. 5d), and the NWR analyses (Fig. S17). The comparison of the time series over the OP sample period of these factor profiles with those of the two OP assays (OP^AA and OP^DTT) is shown in Fig. 6.

Biomass burning contributed greatly to BC_WB (81 %) and to a lower extent to regional background metals (32 %) and MOOA (18 %). This factor accounted for 5.1 % of the total PM₁ concentration. While no primary biomass burning organic aerosol (BBOA) factor was resolved with the PMF_organics analysis in summer, the presence of a significant MOOA contribution reflects the influence of a secondary process in this biomass burning factor. The low concentration of this factor is in agreement with minor regional emissions linked to agricultural activities, wildfires, and cooking practices such as BBQs, transformed through oxidation processes during regional transport and ageing (Chazeau et al., 2022; Cubison et al., 2011). The NWR analysis in Fig. S17 showed biomass burning concentrations associated with higher wind speed compared with sources with a local origin (traffic, shipping, cooking, and ON-rich), corresponding to southwesterly winds from the Mediterranean Sea. Additionally, the northeast land breeze advected these aged emissions back to the sampling site.

Cooking mainly contributes to COA (92 %) and MOOA (34 %) and represents a 14.2 % contribution to the total PM₁ mass. Even if its diurnal pattern is similar to that from primary COA, the factor is mixed with a secondary organic aerosol fraction. Moreover, the cooking source included an unexpected contribution to Cl⁻ (19 %), which was already observed in a rural environment in the Po Valley (Dall'Osto et al., 2015). While this study showed high contributions of Cl⁻ and oxygenated organic aerosol with COA in the cooking source that may be associated with some additional emissions from agricultural activities and waste disposal, these sources are not expected at our urban site.

As presumed, the industrial factor was characterized by high contributions to industrial metals (74 %) and Sh-IndOA (46 %). The factor contributes little to the PM₁ composition (3.2 %), which is expected as the size of the industrial particles generally belongs to the ultrafine mode (< 100 nm) (Riffault et al., 2015). Chazeau et al. (2021) and El Haddad et al. (2013) already described that plumes originating from the main industrial area of Fos-Berre are advected on-site by sea breeze conditions and are mainly attributed to ultrafine particles, influencing the mass concentrations only to a minor extent.

The dust resuspension factor was exclusively driven by metals and showed some identical features to the dust resuspension factor from the PMF_metals (100 % of the variable). It should be noted that the factor also included some weak contributions to the brake and/or tire wear and industrial metals (∼ 10 % each), which can be affected by the resuspension processes.

The traffic source was mostly composed of HOA (32 %), BC_FF (44 %), and to a lesser extent of LOOA (15 %) and also showed a large contribution of these variables (83 %, 79 %, and 18 %, respectively). Tire and/or brake wear metals were also important contributors to this factor (29 % of the total variable). However, this component presented the highest unexplained variation (> 30 %) from all of the dataset (Fig. S14) and revealed some over-splitting in several factors. It should be emphasized that 23 % of the traffic source was composed of SOA (LOOA and MOOA), meaning that the primary traffic contribution is mixed with secondary aerosol concentrations attributed to fast oxidation of freshly emitted particles (Chirico et al., 2011). The factor represented the second highest fraction of the PM₁ mass (21.4 %).

Figure 5 (a) Average factor profiles, with the coloured sticks representing the normalized contribution of the variable to the factor (left axis) and the markers showing the normalized factor contribution to each variable (right axis) for the eight factors from the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ solution. Error bars are the standard deviation of each profile. (b) Time series, (c) pie chart contributions, and (d) diurnal cycles (solid lines indicate the median, red circles the mean, and shaded areas the 25th–75th percentile range) for each factor of the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ solution.

The organic nitrate-rich (ON-rich) factor was resolved based on the high contribution to NO${}_{\mathrm{3}}^{-}$ (81 %) and LOOA (69 %). This factor accounted for 19.4 % of the PM₁. The diurnal trend of this factor suggested contributions from night-time chemistry. A significant pathway might be the oxidation of biogenic VOCs by NO₃ radicals to produce organic nitrate particles (Kiendler-Scharr et al., 2016; Xu et al., 2015). This factor displayed an origin from north to east within the land.

The shipping source showed an expected contribution from shipping metals (87 %) and Sh-IndOA (56 %) and accounted for 5.9 % of the PM₁. This factor further accounts for a noticeable variation in sulfate (11.6 % of the total sulfate concentration). This is in agreement with the results reported by Chazeau et al. (2021), indicating that during 25 % of the days in summer 2017, sulfate concentrations were prominently influenced by the nearby harbour. Shipping emissions from the Marseille harbour are further described in Sect. 3.3.2.

The AS-rich factor represented the largest fraction of the PM₁ (27.9 %) and was dominated by the ammonium sulfate contribution (80 % of SO${}_{\mathrm{4}}^{\mathrm{2}-}$ and 72 % of NH${}_{\mathrm{4}}^{+})$. The SO${}_{\mathrm{4}}^{\mathrm{2}-}$ $/$ NH${}_{\mathrm{4}}^{+}$ ratio in this profile was 2.57, which is consistent with the mass ratio of ammonium sulfate in the (NH₄)₂SO₄ form (2.66). The remaining fraction of SO${}_{\mathrm{4}}^{\mathrm{2}-}$ was attributed to the shipping (12 %) and industrial (4 %) factors. The diurnal profile of the AS-rich factor showed enhanced concentrations in the afternoon linked to the photochemical production of sulfate from its precursor SO₂ (Zhuang et al., 1999). The presence of the SOA contribution in this factor (16 %) suggests its formation by a similar process (i.e. photo-oxidation) (Bozzetti et al., 2017a; Salameh et al., 2018; Waked et al., 2014). This ammonium sulfate pattern was already pointed out in summer in Marseille and was attributed to mixed Mediterranean sulfate sources (including industrial and shipping emissions) from a regional origin transported by processed air masses (Chazeau et al., 2021). This interpretation is supported by the NWR analysis presented in Fig. S17. It should be noted that the AS-rich factor might also include some other anthropogenic influence due to its moderate composition of BC_FF (17 %).

To assess the robustness of the PMF² solution, the results were compared with a PMF solution utilizing the OA factors from PMF_organics, ACSM inorganic species (SO${}_{\mathrm{4}}^{\mathrm{2}-}$, NO${}_{\mathrm{3}}^{-}$, NH${}_{\mathrm{4}}^{+}$, Cl⁻), BC sources, and metal concentrations as the input dataset. Consistent with the PMF² method, constraints, instrument weighting, criteria selection, and bootstrap analysis were applied and are reported in the Supplement Sect. S5. This alternative approach successfully identified the same eight factors (Fig. S18), exhibiting comparable mass contributions and very high correlations with the PMF² factor time series (Table S5), all exceeding an R² of 0.9, except for shipping (R²=0.81).

The biomass burning and shipping factors accounted for slightly higher concentrations in the PMF² solution, due to the slightly elevated contribution of SO${}_{\mathrm{4}}^{\mathrm{2}-}$, NH${}_{\mathrm{4}}^{+}$, and MOOA concentrations which dominate the PM₁ mass. The metal composition found in the factors from this alternative PMF approach is in agreement with the metal profiles from the PMF_metals solution. Note that Zn and Sb, the most prominent elements in the tire and/or brake metal factor, were mainly present in the traffic source. However, they again displayed some mixing with other factors (dust resuspension, AS-rich, and cooking), suggesting additional sources unresolved by the current PMF solutions. Previous studies suggested that Zn may originate from waste incineration or other industrial processes (Belis et al., 2019; Manousakas et al., 2022; Visser et al., 2015a). Comparable results in terms of explained variability were observed, emphasizing the suitability of both methods for such study.

A PMF analysis was also conducted on all instrument datasets (i.e. organic fragments from $m/z\mathrm{12}$ to 100, ACSM inorganic species, BC fractions, and metals) merged into a unique input matrix and did not result in a satisfactory solution (see the Supplement Sect. S6 and Fig. S19). Overall, the present PMF approach successfully identified various sources of PM₁ during the summer season, consistent with previous studies in Marseille. These sources include traffic (El Haddad et al., 2013; Bozzetti et al., 2017a; Salameh et al., 2018), cooking (Bozzetti et al., 2017a), and a minor contribution from biomass burning (Bozzetti et al., 2017a; Salameh et al., 2018). However, this study marks the first identification of an ON-rich factor. A previous source apportionment of PM_2.5 markers by Salameh et al. (2018) highlighted the dominant contribution of ammonium sulfate in summer (35 %) and identified a dust factor with a metal composition similar to the current study (Cu, Fe, Ca). While they identified a fossil fuel factor attributed to mixed harbour and industrial emissions, our results provide new insights by distinctly separating industrial and shipping emissions simultaneously advected on-site by sea breeze.

Figure 6 (a) Time series of both OP assays during the OP sampling campaign. (b) Contribution of source factors provided by the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ to PM₁ over time.

3.3 Results of OP inversion for the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ sources

Associations between each of the sources provided by the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ and the OP measurements are shown in Table S6. The observation of these direct correlations shows that none of the sources identified is dominant on its own in explaining the changes in both OP assays. It is the combination of sources that ultimately leads to the observed OP^AA and OP^DTT.

3.3.1 Model accuracy

The M-estimator inversion model results issued from PMF_organics (scenario 1) or PMF_metals (scenario 2) alone are, respectively, presented in Tables S3a and S3b and are discussed in the Supplement. The results obtained with the WLS inversion applied to the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ (scenario 3) are the most robust and founded on a geochemical base. The cooking source factor was not considered based on its anti-correlation with OP${}_{\mathrm{v}}^{\mathrm{AA}}$. All other source factors included in the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ were considered, since they did not show any multicollinearity (VIF < 5). The accuracy of this model was estimated by its robustness (validity of the error model generated by the bootstrap method) and by a residual analysis between the OP observed and the OP reconstructed by the model. The Breusch–Pagan test performed for each assay assessed the absence of heteroscedasticity in the model residuals ($p<{\mathrm{10}}^{-\mathrm{5}})$ (see Fig. S20). The observed OP and reconstructed OP showed a fair correlation for both assays (OP${}_{\mathrm{v}}^{\mathrm{AA}}$: r=0.44, $R_{\mathrm{adjusted}}^{\mathrm{2}}=\mathrm{0.4}$ – OP${}_{\mathrm{v}}^{\mathrm{DTT}}$: r=0.54, $R_{\mathrm{adjusted}}^{\mathrm{2}}=\mathrm{0.3}$, p<0.001 in both cases).

3.3.2 Intrinsic OP^DTT and OP^AA

Intrinsic OP (i.e β coefficients provided by WLS regression models in scenario 3 (see Sect. 2.5), thereafter denoted OP_m) of the source factor contributions identified by PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ are shown in Table 1 and are discussed here.

Dust resuspension and industrial factors are the main reactive sources to the AA assay, with OP${}_{\mathrm{m}}^{\mathrm{AA}}$ mean values of 0.26±0.03 and 0.22±0.05 nmol min⁻¹ µg⁻¹, respectively. Since 51 % of copper is found in the dust resuspension factor and as various metals (Fe, Cr, Ti, Mn, Pb, or Se) were found in these two factors, the AA assay confirms its metal-sensitivity especially to Cu (Daellenbach et al., 2020; Grange et al., 2022; Pant et al., 2015; Pant and Harrison, 2013). In parallel, Calas et al. (2019) and Gao et al. (2020b) pointed out the role of organic species in OP^AA response, which in this study is sensitive to both less and more oxidized organic aerosols (LOOA and MOOA factors from the PMF_organics) and organic compounds from harbour and plant activities (Sh-IndOA factor also from the PMF_organics). Thus, 27 % and 25 % of the industrial source factor is constituted by LOOA and Sh-IndOA factors, respectively, and 17 % of the dust resuspension source factor is constituted by MOOA factor. As already reported by Weber et al. (2019) in the OP apportionment study on PM₁₀ at the same site, DTT appears to be sensitive to a wide range of sources. In this study, AS-rich, shipping, and biomass burning sources are the main drivers of OP^DTT with, respectively, OP${}_{\mathrm{m}}^{\mathrm{DTT}}$ values of 0.18±0.02, 0.16±0.03, and 0.15±0.06 nmol min⁻¹ µg⁻¹.

In addition to the inherent reactivity of the chemical species, high levels of SO${}_{\mathrm{4}}^{\mathrm{2}-}$ in AS-rich and shipping sources (respectively, 54 % and 36 % of the source factor) may increase OP activity through the dissolution of some metallic elements under acidic conditions (Fang et al., 2017). In the same way, toxicological studies highlighted the role of BC_WB (81 % of the quantification of BC_WB is founded on the biomass burning source) as an indicator of co-transported high DTT-reactive species such as metals or quinones (Niranjan and Thakur, 2017; Shang et al., 2016).

Although several studies emphasized the role of road traffic in OP (Daellenbach et al., 2020; Fang et al., 2016; Saffari et al., 2015), OP_m of the traffic source is surprisingly very low for both OP assays in our case (OP${}_{\mathrm{m}}^{\mathrm{AA}}$: 0.01±0.02 nmol min⁻¹ µg⁻¹–OP${}_{\mathrm{m}}^{\mathrm{DTT}}=\mathrm{0.02}\pm \mathrm{0.02}$ nmol min⁻¹ µg⁻¹). This result might be explained by the fact that non-exhaust traffic emissions traditionally associated with OP are mainly found in a coarser mode than PM₁ (Piscitello et al., 2021).

Interestingly, the traffic source is well correlated with OP${}_{\mathrm{v}}^{\mathrm{AA}}$ and OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ (r=0.40–r=0.34, p<0.01) but associated with the lowest OP_m values. Although Pearson's correlation between OP values and source factor contributions is the first indication of the OP sensitivity towards certain sources, it is most likely preferable to be used in association with an MLR-like model, as already underlined by Weber et al. (2018) and in't Veld (2022).

Table 1 Intrinsic OP^AA and OP^DTT (OP_m) expressed in nmol min⁻¹ µg⁻¹ of sources provided by the PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ method over the OP sampling campaign without the firework episode (n=78 samples) for both OP assays. Values are the mean ± standard deviation from selected bootstraps of the optimal solution.

3.3.3 Population exposure: median contribution of OP

Source-specific contributions to OP^AA, OP^DTT, and PM₁ mass are presented here, ranked in decreasing order, and reported as median values in Fig. 7. The same results gathered in mean values are shown in Fig. S21. These two metrics do not address the same issue: mean value is generally used in the atmospheric community while epidemiological studies prefer to rely on median value, excluding outlier events which are not representative of a chronic exposure of the population. We observe little difference in the ranking of sources between mean and median values due to the overall low variability of the observed OP during the sampling campaign. However, the mean value and the median value of the contribution of the AS-rich source factor to OP^DTT are significantly different, with the mean value being 4 times higher than the median value. The median values of the contribution of the AS-rich source are close to those of biomass burning, shipping, and ON-rich source factor contributions (respectively, with OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ values of 0.1 ± 0.01, 0.09 ± 0.01, 0.08 ± 0.01, 0.08 ± 0.02 nmol min⁻¹ m⁻³).

As already observed in other studies (Borlaza et al., 2022; Weber et al., 2021, 2018), the main observation in Fig. 7 is the clear difference in the factor source contributions when considering the OP activity or PM₁ mass. This highlights the fact that the sources driving OP activity are not the same as those driving PM mass. While the industrial source contributes little to the PM₁ mass (0.26 µg m⁻³), it has the highest OP${}_{\mathrm{v}}^{\mathrm{AA}}$ value (0.06 ± 0.01 nmol min⁻¹ m⁻³) with ON-rich factor (0.05±0.02 nmol min⁻¹ m⁻³) and followed by the dust resuspension (0.04 ± 0.00 nmol min⁻¹ m⁻³) source. In the same way, shipping and biomass burning sources contribute greatly to OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ (0.08±0.03 and 0.09±0.02 nmol min⁻¹ m⁻³) while each source contributes less than 13 % of PM₁ total mass. On the contrary, AS-rich and cooking sources display negative contributions to OP${}_{\mathrm{v}}^{\mathrm{AA}}$ while they contribute significantly to the PM₁ mass. The traffic source does not appear to be a main driver in both types of OP despite its high PM₁ mass contribution (21.4 %). However, we note a very large standard deviation in the contribution of this factor to OP${}_{\mathrm{v}}^{\mathrm{DTT}}$.

A previous study was conducted over 1 year on PM₁₀ at the Marseille-Longchamp site, and already emphasizes the contribution of biomass burning and heavy fuel oil (HFO; related to shipping activity) to OP${}_{\mathrm{v}}^{\mathrm{AA}}$ and OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ (Weber et al., 2021). Although similar results can be found, they should be treated with caution as different fractions of PM are being studied.

Even though OA and metals are found in all sources, the results suggest that only a fraction of these compounds have a substantial impact on OP activity of PM₁. Overall, this draws our attention to the contribution of multiple sources (local and regional) with low PM₁ mass loading in the chronic exposure of PM pollutant.

Figure 7 Median contribution of the sources identified by PMF${}_{{\mathrm{PM}}_{\mathrm{1}}}$ over the OP sampling campaign without the firework episode (n=78 samples) to (a) OP${}_{\mathrm{v}}^{\mathrm{AA}}$, (b) OP${}_{\mathrm{v}}^{\mathrm{DTT}}$, and (c) PM₁. Error bars represent the standard deviation of the data distribution.

3.3.4 Discussion

As shown in Table 1 and Fig. 7, OP${}_{\mathrm{v}}^{\mathrm{AA}}$ and OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ display various sensitivities to sources, and considering both OP assays does not reveal the influence of any particular source. For example dust resuspension and industrial emissions display significant positive OP${}_{\mathrm{v}}^{\mathrm{AA}}$ values and negative or low OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ values. Conversely the shipping emission source and the AS-rich source show significant positive OP${}_{\mathrm{v}}^{\mathrm{DTT}}$ and negative or low OP${}_{\mathrm{v}}^{\mathrm{AA}}$ values. This disparity is most likely associated with different oxidative pathways of the two probes, which account for the diversity of defensive mechanisms operating in the pulmonary environment (Bates et al., 2019). AA is naturally present in the lungs, and its predominant anionic form in solution (HA⁻) is oxidized by various mechanisms facilitated by OH^•, O${}_{\mathrm{2}}^{•-}$, HO${}_{\mathrm{2}}^{•}$, and other radicals, and by transition metals such as Cu (II) or Fe (III) (Campbell et al., 2019). DTT has a disulfide bond and is considered a chemical substitute for cellular reducing agents such as nicotinamide adenine dinucleotide (NADH) or protein thiols (Verma et al., 2015; Borlaza et al., 2018). Protein thiols play an important role in major oxidative stress, restoring the redox balance by eliminating free radicals (Baba and Bhatnagar, 2018). Many studies have linked these two probes (AA and DTT) to transition metals (Cu, Fe, Mn, Zn), EC, and OC (Gao et al., 2020b). In addition, the different sensitivity of AA and DTT to both organic compounds and transition metals has been evidenced by Calas et al. (2018), Gao et al., 2020b and Pietrogrande et al. (2022).

Today, no consensus has yet been reached on which OP test is most representative of health impact, and the community still recommends the complementary use of OP tests, in particular the association of both AA and thiol-based (DTT or GSH) assays (Moufarrej et al., 2020). This association is currently the only way of assessing the full panel of the most oxidizing compounds of PM. However, recent studies have shown positive associations between OP^DTT and various acute cardiac (myocardial infarction) and respiratory endpoints, supporting the relevance of the OP^DTT assay for this purpose (Abrams et al., 2017; Weichenthal et al., 2016; He and Zhang, 2023). On the contrary, several studies did not find an association between OP^AA and health endpoints including early-life outcomes, respiratory and cardiovascular mortality, cardiorespiratory emergencies, and lung cancer mortality (Borlaza et al., 2023; Marsal et al., 2023b). Nonetheless, a recent study has associated OP^AA with oxidative damage to DNA (Marsal et al., 2023a). These results so far may suggest that OP^AA provides partial information on the link between OP and adverse health effects, and further epidemiological studies are needed to determine whether OP^AA should be considered as a proxy for health impact.

4 Limitations of the study

The robust method used in this study is promising for further practical applications. However, there are some limitations that can explain the moderate OP reconstruction. Mainly, little variability in OP measurements and in the chemical composition was observed at the Marseille-Longchamp site over a 15 d period, implying difficulties for models such as MLR to accurately reconstruct OP peaks and thus lowering the Pearson correlation coefficient r between observed OP and modelled OP. Then, averaging the online data which have a time step of 15 min and 1 h, respectively, for ToF-ACSM and Xact over the time step of the OP (4 h) led to a levelling of the pollution peaks visible on the online analysers. Also, the times series of OP and chemistry originates from offline and online methodologies that may introduce additional uncertainties to the results of the study, related to some extent to the sampling artefact. Furthermore, recent studies suspected non-linear relationships for the source apportionment of OP, especially when increasing PM mass (Calas et al., 2017; Charrier et al., 2016; Grange et al., 2022; Weber et al., 2021). Another main limitation is the lack of data between 12 July 2018 21:00–13 July 2018 05:00 UTC (3 points) and 14 July 2018 17:00–15 July 2018 05:00 UTC (4 points), resulting in difficulty for the model to reconstitute certain OP peaks. Lastly, the propagation of uncertainties can be a limitation as we increase the imprecision of the actual measured OP. The combination of all these limitations could explain the unrecovered percentage of the OP variance by the model (i.e. intercept).

5 Conclusions

To the best of our knowledge, this study is the first to apportion OP from sources provided by a two-step PMF approach using high time resolution online Xact and ToF-ACSM data. The PMF² approach successfully identified eight well-resolved sources (AS-rich, traffic, ON-rich, cooking, shipping, biomass burning, industrial, and dust resuspension), a solution not achievable through single PMFs conducted separately on OA and metal datasets. The method enabled the assignment of OA factors, which typically described components arising from a mixture of sources and chemical processes rather than a single emission source, to more specific PM₁ sources. Additionally, this approach allowed us to assess both the primary and secondary origin of anthropogenic sources, such as traffic and cooking. However, a limitation of this method is that non-explained variability and uncertainties of the factors from the first-step PMFs will propagate into the PMF² results and therefore need to be carefully assessed. The inclusion of additional element measurements, such as Ba, S, Cl, and Si, to the PMF_metals could be an interesting feature for refining some sources and addressing this limitation.

A redistribution of the sources between mass and OP contributions was observed in both OP^AA and OP^DTT assays, highlighting the limiting capacity of the mass concentration alone in understanding the redox activity of PM. While the PM₁ mass concentrations were dominated by AS-rich, traffic, and ON-rich sources, factors with an anthropogenic origin such as industry, shipping, biomass burning, and dust resuspension, represented only a few percent. However, these factors with a low mass concentration showed a high OP contribution. Therefore, making an effort to reduce these specific anthropogenic sources, regardless of the PM mass regulations, might result in a significant reduction of OP in the submicron mode.

We demonstrated that OP apportionment models performed on source factors deconvolved from partial aerosol composition (only metals or organics) are not robust. The two-step PMF² approach appears to be a substantial method for integrating a wider range of markers (organic and inorganic) for a better identification of PM sources. The OP apportionment via a WLS inversion model led to comprehensive and realistic OP sources.

Several improvements can be suggested to address the limitations mentioned in Sect. 4, such as the subsequently introduced source uncertainties from the PMF model into the OP deconvolution model. Furthermore, these findings pave the way for the coupling of recent prototypes allowing for in situ OP data with online chemical analysers aiming to understand short-lived processes.