Significance A successful mitigation strategy for climate warming agents such as black carbon (BC) requires reliable source information from bottom-up emission inventory data, which can only be verified by observation. We measured BC in one of the fastest-warming and, at the same time, substantially understudied regions on our planet, the northeastern Siberian Arctic. Our observations, compared with an atmospheric transport model, imply that quantification and spatial allocation of emissions at high latitudes, specifically in the Russian Arctic, need improvement by reallocating emissions and significantly shifting source contributions for the transport, domestic, power plant, and gas flaring sectors. This strong shift in reported emissions has potentially considerable implications for climate modeling and BC mitigation efforts.

Abstract Black carbon (BC) in haze and deposited on snow and ice can have strong effects on the radiative balance of the Arctic. There is a geographic bias in Arctic BC studies toward the Atlantic sector, with lack of observational constraints for the extensive Russian Siberian Arctic, spanning nearly half of the circum-Arctic. Here, 2 y of observations at Tiksi (East Siberian Arctic) establish a strong seasonality in both BC concentrations (8 ng⋅m−3 to 302 ng⋅m−3) and dual-isotope–constrained sources (19 to 73% contribution from biomass burning). Comparisons between observations and a dispersion model, coupled to an anthropogenic emissions inventory and a fire emissions inventory, give mixed results. In the European Arctic, this model has proven to simulate BC concentrations and source contributions well. However, the model is less successful in reproducing BC concentrations and sources for the Russian Arctic. Using a Bayesian approach, we show that, in contrast to earlier studies, contributions from gas flaring (6%), power plants (9%), and open fires (12%) are relatively small, with the major sources instead being domestic (35%) and transport (38%). The observation-based evaluation of reported emissions identifies errors in spatial allocation of BC sources in the inventory and highlights the importance of improving emission distribution and source attribution, to develop reliable mitigation strategies for efficient reduction of BC impact on the Russian Arctic, one of the fastest-warming regions on Earth.

Black carbon (BC) is a short-lived climate pollutant, formed during incomplete combustion of biomass and fossil fuels and contributes to the amplified warming in the Arctic (1⇓⇓–4). However, estimates of the magnitude of added radiative forcing to the global atmosphere by BC span a large range (0.2 W⋅m−2 to 1 W⋅m−2) (1, 5). Due to its short atmospheric lifetime, BC is a potential target for climate change mitigation. Historically, BC concentrations have been decreasing in the Arctic air (6), but their future fate is unclear. Projections range from increasing concentrations due to a decrease in rainfall (wet scavenging) (7), changes in wind patterns (8), an increase in emissions from wildfires (9), and increased shipping and extraction of natural resources (10) to decreasing concentrations due to more efficient wet scavenging (8). Chemical transport and climate model predictions of BC in the Arctic were, until recently, unsatisfactory and failed to reproduce the observed magnitude and amplitude of BC concentrations (11, 12). However, developments in atmospheric transport models show increasing model skills (12), especially for the European Arctic (13). A key component to the model−observation offset is the large uncertainty connected to emission inventories (EIs) (14⇓–16) used by the models. Implications of these model uncertainties include challenges of accurately assessing the radiative forcing of BC (1, 17). Improvements in simulating BC concentrations have been shown to also improve simulation of aerosol optical depth (18). Several recent assessments urge for observation-based source apportionments to improve BC EIs (1, 19⇓⇓⇓–23).

EIs are based on a bottom-up approach and heavily influenced by assumptions on emission factors (i.e., amount of BC released per amount of burned fuel in a given technology/source sector) and activity (i.e., amount of burnt fuel). Further, the spatially disaggregated national or regional inventories often use (coarse) large-scale proxies, such as total population, for aggregated sector categories (e.g., land transport, residential combustion), leading to spatial misallocation of sources. Even the most recent EI for Russia has an estimated range of total BC emissions that spans one order of magnitude (95% confidence) for anthropogenic emissions alone (24). This emphasizes the need for further source apportionment studies in this key source region for the Arctic, to eventually better understand climate effects of BC in the Arctic.

Carbon isotope characterization of elemental carbon (EC) aerosols [a BC analog (25)] has proven to be an important tool to constrain different BC source contributions (biomass burning vs. fossil fuel) and diagnose causes for differences between 14C-based diagnostic source apportionment and EI models in different parts of the world (16, 26⇓⇓⇓⇓–31). Here, we present a 2-y study combining radiocarbon and stable carbon isotope analysis of EC in the vastly understudied Siberian Arctic (Tiksi). The total aerosol loading was collected continuously (16 April 2012 to 07 March 2014), with a high-volume sampler for total suspended particles (TSP), to meet the sample size requirements for natural abundance 14C of the EC fraction. Many BC sources for Russia are reported as TSP by Russia's Federal State Statistics Service (24). These observations were directly compared with corresponding simulations with the atmospheric transport model FLEXPART (flexible particle dispersion model) for the same site, coupled to the EI ECLIPSE (Evaluating the Climate and Air Quality Impacts of Short-Lived Pollutants) (23, 32), which was recently proven to emulate observations in the European Arctic well (13). The ECLIPSE EI can probably be considered the most suitable global EI for the Russian Arctic, because it also includes BC emissions from gas flaring, thought to be responsible for a significant portion of the Arctic BC burden (11, 24). Because ECLIPSE only covers anthropogenic emissions, open fire emissions (including agricultural waste burning and wildfires) were included by using the Global Fire Emissions Database (GFED), based on satellite data (33). The model results were then subject to forward isotope modeling to create a best-fit scenario for all modeled sources. This combination allowed for a direct comparison between top-down measurements and bottom-up simulations of BC concentrations and sources of biomass burning and fossil fuels, including biofuels, open fires, coal, liquid fossil fuels, and gas flaring.

SI Methods Carbon Isotope End-Member Determination. The end-member values (source-specific signatures) of the carbon isotopes were compiled from the literature (Table S2). The fossil Δ14C sources (liquid fossil, coal, and gas flaring) are well constrained at −1,000 ± 0‰. The main source of biomass in Russia is wood burning, estimated at +225 ± 60‰ (13, 31). The δ13C end-member values for coal (−23.4 ± 1.3‰) and biomass (C3 plants, −26.7 ± 1.8‰) were collected from the literature review conducted by Andersson et al. (28). The liquid fossil end-member was taken from typical liquid fossil sources used in Russia (−31.4 ± 1‰) (40). The perhaps most uncertain end-member is the δ13C value for gas flaring. Here a wide variability of 3‰ is used, but the average position may be dependent on the relative contribution from the different gas components, e.g., methane, ethane, and propane, which all have quite varying δ13C-value signatures, although methane is expected to dominate the mixtures (δ13C for methane ≈ –60‰). Notable is that BC formation, in general, is associated with an enrichment in the δ13C isotope (higher δ13C), because the lighter 12C is more prone to CO 2 formation. However, an enrichment of more than a few per mill is not expected. In general, δ13C lower than –38‰ is expected for certain situations (corresponding to an even lower posterior flaring contribution). Here, –38‰ was used in agreement with the (to our knowledge) only published δ13C characterization of BC from gas flaring (44). Domestic is a mixed source (60% biomass, 39% coal, and 1% liquid fossil) (24). The end-member values (for both δ13C and Δ14C) were estimated assuming normal distribution mixing, where the mean (μ) and SD (σ) were estimated as: μ D = 0.6 × μ B + 0.39 × μ C + 0.01 × μ L [S1a] σ D = ( 0.6 × σ B ) 2 + ( 0.39 × σ C ) 2 + ( 0.01 × σ L ) 2 , [S1b]where D = domestic, B = biomass, C = coal, and L = liquid fossil. Estimating the SD of the Priors. The fractional source contribution priors for the Bayesian analysis are assumed to be normally distributed with a mean (μ) equal to the FEG model results. The second parameter in the normal distribution is the SD (σ), which also needs to be estimated. BC bottom-up EIs typically report large uncertainties (e.g., 125 to 500%) in the flux (14, 15). However, the uncertainties for the fractional contributions from different sources need to be lower (<100%), as constrained by the mass balance criterion. To fully estimate the uncertainties of the relative source contributions to Tiksi BC levels from the FEG model would mean a vast modeling effort, i.e., coupled uncertainty propagation from FLEXPART, ECLIPSE, and GFED. Instead, given the comparably large expected uncertainties of the FEG modeling results, we estimated the uncertainties (σ) using the weakest assumptions of prior knowledge: Each possible fractional source combination is equally probable. The distribution that describes this scenario for n sources (here n = 5) is the n-dimensional Dirichlet distribution where the shape factors (α) are all equal to 1 (a.k.a., the n-dimensional standard simplex). The SD of the marginal distribution (the distribution of one of the fractional source contributions) of this symmetric Dirichlet distribution is given by: σ = n − 1 n 2 × ( n + 1 ) . [S2]Five sources thus give an uncertainty of roughly 16.3%. For comparison, the two-source case of this Dirichlet distribution is the standard uniform distribution, with σ ≈ 28.9%. MCMC Analysis. The MCMC analysis was conducted using in-house scripts written in Matlab (ver. 2014b). A Metropolis−Hastings algorithm (63, 64) was used for sampling the parameter space, using 1,000,000 iterations, a burn-in of 10,000, and a data thinning of 10. The jump size of the stochastic perturbation was adjusted to obtain an acceptance rate of ∼0.23. The end-member distributions were described with normal distributions (28). The different samples were given different weight in the fitting, depending on observational sampling duration. The prior relative source distributions before the Bayesian modeling were assigned with normal distributions with (mean) values from the FEG model, and with an SD of 16.3%. For the sources where a general shift was introduced, the prior was exchanged from one per sample to one parameter affecting all samples, where the shift (here exemplified by gas flaring, x F ) was allowed to take any positive number, exemplified as: f F ( j ) p o s t e r i o r = x F × f F ( j ) p r i o r x F × f F ( j ) p r i o r + f O F ( j ) p r i o r + f C ( j ) p r i o r + f L ( j ) p r i o r + f D ( j ) p r i o r [S3a] f O F ( j ) p o s t e r i o r = f O F ( j ) p r i o r x F × f F ( j ) p r i o r + f O F ( j ) p r i o r + f C ( j ) p r i o r + f L ( j ) p r i o r + f D ( j ) p r i o r [S3b] f C ( j ) p o s t e r i o r = f C ( j ) p r i o r x F × f F ( j ) p r i o r + f O F ( j ) p r i o r + f C ( j ) p r i o r + f L ( j ) p r i o r + f D ( j ) p r i o r [S3c] f L ( j ) p o s t e r i o r = f L ( j ) p r i o r x F × f F ( j ) p r i o r + f O F ( j ) p r i o r + f C ( j ) p r i o r + f L ( j ) p r i o r + f D ( j ) p r i o r [S3d] f D ( j ) p o s t e r i o r = f D ( j ) p r i o r x F × f F ( j ) p r i o r + f O F ( j ) p r i o r + f C ( j ) p r i o r + f L ( j ) p r i o r + f D ( j ) p r i o r , [S3e]where F = gas flaring, OF = open fires, C = coal, L = liquid fossil, D = domestic, and f denotes the corresponding fractional contributions. The index j denotes sample number. Thus, if x takes on a value larger than 1, the overall fractional contribution of source gas flaring increases (and all others decrease), and, if it is smaller than 1, it decreases (and all others increase). Bayesian Model Comparison. For comparing the significance of the data fit for different models, evaluation of the Bayes factor is a powerful approach. In the present case, the different models that were compared shifted the importance of one or more sources to better fit with the observational data, i.e., the Δ14C and δ13C data for EC. The Bayes factor (K) is defined as the conditional probability (P) of data (D) given model 1 (the evidence for model 1, M 1 ) divided by the conditional probability for data given model 2 (the evidence for model 2, M 2 ): K = P ( D | M 1 ) P ( D | M 2 ) . [S4]The evidence for a model was computed by integrating (marginalizing) over the parameters assigned to that model. As an example, the Bayes factor comparing shifting both flaring (F) and open fires (OF) vs. shifting only flaring is given by the ratio between the evidence for the first model divided by the evidence for the second: K = ∬ P ( D | x F , x O F , M F , O F ) × ​ P ( x F , x O F | M F , O F ) d x F d x O F ∫ ​ P ( D | x F , M F ) × P ( x F | M F ) d x F , [S5]where x denotes the shifting factor. For the present example, two shifting factors are associated with model 1 and one shifting factor for model 2. The evidence for each model was computed separately (within model approach) using MCMC sampling over the parameter space, and the Bayes factor was computed through combinatorial comparison of the different models. For this particular case (comparing the two models with the best fits), the Bayes factor is >100, which is a “decisive” or “very strong” in favor of the combination of combining gas flaring and open fires. Adding a third parameter, e.g., coal, does not significantly improve the fit, as the Bayes factor is ∼1. An advantage with using Bayes factors for model comparisons is that the degrees of freedom, e.g., the number of sources that are shifted, are naturally incorporated in the integral in Eq. S5.

Discussion Observations and model predictions of BC are disparate in the Russian Arctic and other high latitudes (12, 24, 47). At the same time, the Russian Arctic is one of the fastest-warming regions on our planet (48), especially in the spring season (high-temperature anomaly for March−April−May). The reasons for this anomaly are unclear. However, the BC burden for Tiksi is higher than for other remote sites, and it is especially high during the Arctic haze season in spring. Our 2-y observations showed clear and strong seasonality with respect to sources and concentrations of BC. The isotope-constrained source apportionment pinpointed that biomass burning was dominant during low BC concentration seasons in summer whereas fossil sources dominated the Arctic haze seasons. Although the FEG model overpredicted contributions of gas flaring and biomass burning to BC, the observed isotopic signature pointed to liquid fossil fuels and biomass burning as the major sources, or, in policy terms, the transport and domestic emission sectors. In this, and a previous independent isotope-based study of BC in the European Arctic (13), we found that gas flaring did not play a dominant role, contrary to expectations from current EIs. Both receptor sites (Abisko in Sweden and Tiksi in Russia) were at a similar, rather large, distance from the major gas flaring sources in the Northwestern and Ural Federals Districts, which could be one explanation of this finding. The applied model did not agree as well with observations in the northern Siberian Arctic as it did for the European Arctic. This was, in part, due to overestimation of biomass burning BC (fire emissions inventory) and underestimation of fossil fuel sources (anthropogenic EI). The latter is most likely caused by an inchoate EI, where entire regions and major sources, most likely of the liquid fossil fuel type, appear to be completely missing or significantly underestimated, e.g., road transport or shipping. The appearance of white spots in the EI and underestimation of observed BC in Tiksi could also be due to misallocation of sources, such as mining or nonroad transport, i.e., shipping on rivers, as well as transportation on frozen rivers. Smaller settlements could likely be misrepresented in the coarse spatial pattern of the inventory, and low local emissions could have a significant impact on observations. However, the overestimation of fossil (gas flaring) sources with occasional overestimation of biomass burning sources are the reason why predicted BC sources based on ∆14C appear to be quite accurate, looking at the whole campaign. Application of a Bayesian statistics-based model allowed creating a best-fit scenario (posterior) to the originally predicted values by the model (prior). This posterior result underlined our previous observational findings. The contributions of gas flaring were overestimated, whereas contributions of liquid fossil (transport) and (domestic and power plant) coal were underestimated. The few exceptions of biomass burning overestimation can potentially be explained by positive fire count anomalies of local fire patches (33). Although the discrepancy of the gas flaring estimates and the observations could be due to issues within FLEXPART (e.g., too little wet scavenging causing a too-long BC lifetime and thus a too-large contribution of the remote flaring emissions relative to closer other sources), a too-long BC lifetime in the model would also lead to a general overestimation of observed BC, which is not the case. In addition, the gas flaring end-member could deviate from the herein applied value. This end-member depends on the hydrocarbon composition of the flared gas and, to a certain degree, on the biogeochemical origin of the hydrocarbons (49, 50), i.e., the isotope end-member is most likely dependent on the local origin of the gas flaring BC. Currently, only one study on the end-members of gas flaring exists (44), illustrating the need for more measurements, especially close to sources. Lastly, the model−observation dichotomy could be due to difference in actual 2012–2014 BC emissions versus the emission estimates used in this work, i.e., from the year 2010.

Conclusions This continuous 2-y study provided an opportunity to compare observations and model predictions of BC in the vastly understudied northeastern Siberian Arctic. Our isotope-based observations showed that sources of pollution were, to a significant extent, due to liquid fossil fuel of Russian origin, identified by its distinct isotopic signature. A high seasonality was observed with regard to both BC concentrations and sources. Regardless of local pollution, the MCMC calculations showed that a best-fit model could be generated by reducing gas flaring (–84%) and open fire emissions (–53%) and increasing transport (+139%), domestic (+113%), and power plant (+109%) emissions. The emissions accounted for in the EI (which is coupled to the model) appear to be uncertain, and sources related to BC emissions received in northeastern Siberia seem to be poorly allocated in the inventory. Apparent misestimation of gas flaring and other sources could be due to the misallocation of distant or close sources, relative to the remote receptor site, or simply missing local sources. This uncertainty may also be higher due to, e.g., inaccurate scavenging coefficient and aerosol lifetime in the model, as well as unexpected end-member variations in the isotopic sources of gas flaring in the observation. BC concentrations at Tiksi are higher than at other remote sites in the Arctic, which suggests that local sources enhance the otherwise uniform Arctic background at Tiksi. However, there is a near-absence of BC emissions in northeastern Siberia in the EI inventory. This 2-y continuous record provides a robust baseline for assessing the current loadings, and for diagnosing areas for further improvement in modeling and observations—a necessity also for more reliable estimates of the climate effect of BC in the rapidly warming Russian Arctic.

Methods The Far-East Siberian Arctic Receptor Site. The high-volume filtration systems were continuously operated for 24 mo in a specially constructed new sampling cabin at the observatory of the Russian Academy of Sciences (RAS), situated ∼10 km southwest of the Tiksi settlement to be away from immediate influence (Fig. S4). Tiksi is an urban settlement (population ∼2,000) situated on the Lena River Delta (34). The RAS Polar Geocosmophysical Observatory site (71.4°N, 128.5°E) was founded in 1958. Permanent technical staff attend the station on a daily basis. This study also benefitted from measurements of absorption-based BC measurements taken at the Tiksi Hydrometeorological Research Observatory, which opened in 2010, based on the 1932-founded Polyarka Station, 7 km south of the Tiksi settlement (34, 51). Fig. S4. The study site. Elevation map of the RAS Polar Geocosmophysical Observatory (star) and Hydromet observatory sites (black circle) in relation to Tiksi village (red circle); a.s.l., above sea level. EC and OC Analysis. Samples were collected from April 2012 to March 2014 with continuous sample intervals of 15 d to 25 d, depending on the weather conditions. The coldest winter temperatures turned out to be too demanding for the tubing used for the TSP inlet, which is why there is a 1-mo gap in observations (14 January 2013 to 6 February 2013), during which the tubing had to be replaced with a (low) temperature-resistant type. Aerosols were collected on precombusted quartz fiber filters (Millipore) using high-volume sampling with a TSP inlet (custom-built at Stockholm University). Carbonaceous aerosol concentrations (EC and OC) were measured with a standard thermal−optical transmission (TOT) analyzer (Sunset Laboratory Inc.) using the National Institute for Occupational Safety and Health 5040 method (52). Parts of the OC could potentially char during the application of this method. The OC would end up as pyrogenic carbon in the EC fraction, leading to an overestimation of the fraction of biomass burning. This potential effect has been evaluated in earlier work by sensitivity analysis, where it was found that the fraction of biomass burning could, in extreme cases, be overestimated by up to 7% (31). A total of 10 field blanks has been analyzed, all having EC concentrations below detection limit. Carbon Isotope Analysis. The 33 samples were pooled into 17 composites, with emphasis on higher temporal resolution during the Arctic haze period. The isotopic analysis of EC was performed as described in previous work (13, 26, 31, 53). Briefly, the EC fraction was cryogenically trapped for further off-line isotopic analysis after regular Sunset TOT conversion to CO 2 . Total sample size was at least 40 μg C. Both carbon isotopes were analyzed using accelerator mass spectrometry (AMS) at the United States National Science Foundation National Ocean Sciences Accelerator Mass Spectrometry Facility (53⇓–55). Fractions of biomass burning (which includes anthropogenic biofuels or wood combustion as well as natural wild fires) of EC were calculated based on the radiocarbon result (13), with an isotopic mass balance equation (26). Contemporary atmospheric CO 2 and freshly produced biomass have a ∆14C end-member of ∼+25‰ (56). However, the most common source of biomass burning in the Arctic is wood, which has reported contemporary end-members between +90 and +282‰, depending on age and species (region) (56). Here, we apply a mean isotopic end-member for biomass BC of +225 ± 60‰, assuming a similar wood fuel turnover time and biota for the Eurasian boreal and northern temperate regions (31). The conservative variability of ±60‰ was restricted by MCMC simulation, leading to a mean variability (plus SD) in the fraction biomass burning of 3 ± 2%. Open Fire Estimate by Satellite-Based Fire Emissions Inventory. The Global Fire Emissions Database, version GFED4.1s, was applied to get an estimate of the biomass burning contribution to BC (13). This fire emissions inventory is based on satellite data and quantifies open fires, as well as fires from agricultural waste burning (33, 57), as burned area product (Collection 5.1 Moderate Resolution Imaging Spectroradiometer) (33, 57, 58). The original resolution of 0.25 was changed to 0.5°. The monthly dataset includes small fires (59). Transport Modeling with Bottom-Up Emission Inventory. To predict the BC concentrations at Tiksi, the particle dispersion model FLEXPART (60, 61), version 9.2, was applied in backward mode for the exact same time periods as the measurements (13). The simulations used meteorological operational analysis data from the European Centre for Medium-Range Weather Forecasts (ECMWF) at a resolution of 1°. One data point (in summer 2013) is missing due to ECMWF’s increase of vertical model resolution on 25 June 2013. Simulations extended over 20 d back in time, sufficient to include most aerosol emissions arriving at the station, given a typical BC lifetime (∼1 wk). A mean particle diameter of 250 nm was used, with a logarithmic size distribution and a logarithmic SD (sigma) of 1.25. To estimate the anthropogenic BC emissions, FLEXPART was coupled to the ECLIPSE, version 5, EI based on the Greenhouse gas–Air pollution Interactions and Synergies model (62), using the year 2010 baseline scenario (23, 32). All emissions were available at yearly resolution for various source types, which were split into monthly resolution using monthly disaggregation factors from the ECLIPSE data set. Emissions from agricultural waste burning were excluded, because those were included in the GFED. Additionally, all emissions were explicitly split between biofuels (modern; e.g., wood burning) and fossil fuel emissions (Table S4). Table S4. Emission partitioning of the ECLIPSE EI data Data Availability. The observational data that support the findings of this study are available on request from the corresponding author (Ö.G.) and will be available in the Bolin Centre Database (bolin.su.se/data/). EI data for GFED are freely available and can be found on the website www.globalfiredata.org/data.html. The FLEXPART model is freely available to the scientific community. It can be accessed under https://www.flexpart.eu/. For an ECLIPSE version with emissions split into fossil and biofuel, please contact Z.K. or C.H. directly. The data for total emissions of BC for different emission scenarios of ECLIPSE are freely available from IIASA: www.iiasa.ac.at/web/home/research/researchPrograms/air/Global_emissions.html.

Acknowledgments Eija Asmi and John Backman (Finnish Meteorological Institute) are acknowledged for numerous helpful discussions. We are also grateful to two anonymous reviewers for helpful comments on an earlier version of the manuscript. A.A., Ö.G., and P.W. acknowledge financial support from the Swedish Energy Agency (Contract 35450-2), the European Union under the 7th Framework Programme (FP7) [Compound Specific Isotopes (CSI):Environment, Contract PITN-GA-2010-264329], the Nordic Council of Ministries Defrost project as part of the Nordic Centre of Excellence, the Swedish Research Council Formas (Contract 942-2015-1070), and the European Research Council [ERC-Advanced Grant (AdG) Cryosphere-Carbon on Top of the Earth (CC-Top) Project 695331]. I.P.S. and A.C. acknowledge support from the Russian Government (Contract 14.Z50.31.0012/03.19.2014). N.S. and O.D. acknowledge the Russian Science Foundation (Contract 15-17-20032). ECMWF is acknowledged for meteorological data, and part of this study [Norwegian Institute for Air Research (NILU) and International Institute for Applied Systems Analysis (IIASA)] was funded under the FP7 ECLIPSE (Project 282688) and the Norwegian Research Council (NFR) project Emissions of Short-Lived Climate Forcers near and in the Arctic (SLICFONIA).