Responses of wildfire-induced global black carbon pollution and radiative forcing to climate change

The impacts of climate change on wildfires have been studied extensively. Along with declining emissions from fossil fuel combustion due to anthropogenic emission control, black carbon (BC) released from wildfires is expected to contribute a more significant portion to its atmospheric burden. However, from a global perspective, little is known about the BC burden and radiative forcing caused by wildfires. Here, we report the results from the long-term wildfire-induced BC concentration and direct radiative forcing (DRF) from 1981 to 2010 globally simulated by an Earth System Model using an updated wildfire BC emission inventory. We show that wildfire-induced BC concentration and DRF varied significantly spatially and temporarily, with the highest in sub-Saharan Africa, attributable to its highest level of wildfire BC emission worldwide. The temporal trends of near-surface air temperature, precipitation, and evapotranspiration and their association with wildfire-induced BC concentration are explored using the multidimensional ensemble empirical mode decomposition. A statistically significant relation between changes in climate parameters and wildfire-induced BC concentration was found for 53% of the land grid cells, explaining on average 33% of the concentration variations. The result suggests that the wildfire-induced BC DRF, with an increasing trend, would be an important contributor to climate change, especially in sub-Saharan Africa. Occurrences of wildfires in the Amazon Basin respond most strongly to climate change and play an increasingly important role in changing the global climate.


Introduction
Wildfire exerts a substantial impact on ecosystems and human societies.While anthropogenic burning activities of forests and other vegetation are conducted for human benefits (for example, the preparation for cropland), wildfire is considered a natural hazard that can cause adverse outcomes, including the loss of lives and damage to the public infrastructures (Moritz et al 2014, Bowman et al 2017).Wildfire is prominent in boreal forest and rainforest zones, such as in sub-Saharan Africa, South America, and Australia (Fried et al 2004, Chapin et al 2008, Macias Fauria and Johnson 2008, Wulder et al 2009, Wotton et al 2017, Kraaij et al 2018), and plays a vital role in land cover changes (Ehrlich et al 1997, Eva and Lambin 2000, Hoscilo et al 2011).Due to the nonuniformity of fuel types, weather, and terrain, wildfires can be spatially and temporally heterogeneous (Cochrane 2009).The flammability of various vegetation on different land cover varies.For example, deciduous broadleaved forests are less fire-prone than shrublands (Moreira et al 2012).In forest regions where the amount of fuel is unlimited, the aridity of fuel becomes the major constraint of wildfire activities (Abatzoglou and Williams 2016).Meteorological conditions like temperature, humidity, and wind speed can influence the flammability of biomass fuel, thus modulating the occurrence of wildfires (Pinol et al 1998).Under global climate change, the area burned and the duration of wildfire activities have been altered significantly (Marion et al 2009, Liu et al 2013, Abatzoglou and Williams 2016, Bowman et al 2017).Wildfires can influence extreme rainfall postwildfire (Touma et al 2022) and drought (Whitman et al 2019).The linkages between wildfire and anthropogenic climate change have been widely studied using well-established weather and fire danger indices like vapor pressure deficit, Palmer drought severity index, and McArthur forest fire danger index (Kraaij et al 2018, Williams et al 2019, Canadell et al 2021).It has been reported that anthropogenic climate change is a significant contributor to growing wildfire activities (Westerling and Bryant 2007, Marion et al 2009, Liu et al 2010, 2013, Westerling et al 2011, Abatzoglou and Williams 2016, Williams et al 2019).
As a significant component of biomass burning, wildfire is a natural emission source of greenhouse gases (GHGs) and aerosols, such as black carbon (BC).Efforts have been made to develop BC wildfire emission inventories, such as the fourth generation of the Global Fire Emissions Database (GFED4) (Giglio et al 2013) and Peking University (PKU) inventory (Wang et al 2014).However, due to the risks and unpredictability, it is still difficult to conduct observational field campaigns to quantify the emissions of wildfire (Alves et al 2011).Attempts have been made to measure the concentration of many air pollutants in wildfire plumes through aircraft and satellite remote sensing and to estimate the wildfire emissions (Cook et al 2007, Alvarado et al 2010, Akagi et al 2011, Urbanski 2013, Liu et al 2017).The evaluation of wildfire emissions often solely relied on fuel consumption data because of the lack of measurements of wildfire emission factors (van der Werf et al 2017) or adopted emission factors estimated from temperate prescribed burning beyond the wildfire season or laboratory combustion (Yokelson et al 2007, Akagi et al 2011, Alves et al 2011, Mebust et al 2011, Urbanski et al 2011, Urbanski 2013, Liu et al 2017).The injection height of wildfire emissions is also crucial to air quality modeling, which has been thoroughly studied using semi-empirical formulas and satellite observations to constrain the vertical profile of wildfire emissions (Rio et al 2010, Sofiev et al 2013, Paugam et al 2016).And it has been revealed that wildfire plume top height showed enhanced trend throughout the Western US (Wilmot et al 2022).
GHGs and some particulate climate forcers like BC released by wildfire can alter regional and global climate significantly (Phillips et al 2022).The significance of wildfire-induced BC emissions has received increasing attention in recent years, along with declining BC emissions from anthropogenic sources due to worldwide BC emission control.Many studies have revealed that aerosols emitted by wildfires tend to yield a negative radiative forcing at the surface and lower the surface temperature (Ward et al 2012, Athanasopoulou et al 2014, Chang et al 2021).However, during an intense wildfire event in Europe, smoke aerosols caused an increase of 10-35 W m −2 atmospheric radiative forcing and leaded to decreasing photolysis rates (Hodzic et al 2007).Escalating carbon emissions from North American boreal forest wildfires were considered a threat to the goal of limiting temperature rise within 1.5 • C (Phillips et al 2022).
Originating from incomplete combustion, BC can absorb solar radiation, thus warming the atmosphere.Compared with BC from fossil fuel combustion, BC from biomass burning has larger sizes, more coated particles, thicker coatings, and more absorption per unit mass (Schwarz et al 2008).The rapid increase in BC levels can be observed during wildfire events compared to the pre-and post-wildfire periods (Athanasopoulou et al 2014).It is estimated that BC accounts for about 10% of wildfire plume mass and is a most critical climate agent subject to aerosolradiative forcing introduced by wildfires (Veira et al 2016).Yu et al reported that BC heated the atmosphere and facilitated the rise and spread of the plume during a persistent wildfire event in 2017 in western Canada (2019).
Wildfire emission varies significantly among different places and years (Mieville et al 2010).To assess the air pollution and climate change induced by wildfire BC emission, we simulated the longterm temporal and global distribution of BC concentration (ng kg −1 ) and its direct radiative forcing (DRF, mW m −2 ) from 1980 to 2010 using the Community Earth System Model version 2 (CESM2).Peking University PKU-BC-v2 BC emission inventory (Wang et al 2014) was employed in the present modeling investigation.We explore the relationship between BC concentration and two meteorological indices to paint a new picture of the nexus between wildfire, air quality, and climate change.

CESM simulations
We use CESM2 scientifically validated component set FHIST with historical CAM6 using 1 degree finite volume dycore with prescribed sea surface temperature and sea ice concentrations to simulate the BC concentration and DRF from 1980 to 2010, with the first year (1980) as the model spin-up time.The major components of the CESM2 include Community Atmospheric Model version 6 (CAM6) with Modal Aerosol Mode with 4 modes (MAM4) and the Community Land Model version 5.0 (CLM5).

X Liu et al
The CESM has 32 layers in the vertical, and the horizontal spatial resolution is 0.9 • × 1.25 • latitude by longitude.The model runs at a time step of 30 min.
In MAM4, aerosols are divided into four modes according to their states and particle sizes.Newlyemitted BC from wildfire emission is categorized as the primary carbon mode and is gradually transferred to accumulation mode as its size grows.The resulting BC concentrations come from both the primary carbon mode bin and the accumulation mode bin.Considering that the air quality near the surface is more important to the ecosystems and human health, we mainly focus on the BC concentration at the lowest atmospheric level above the ground surface in the CESM2.The BC DRF is calculated by the difference between the net radiation at the top of the atmosphere (short-and long-wave radiation) with and without the inclusion of BC.Modifying the radiatively active aerosol species definition for different modes allows the model to calculate the BC DRF diagnostically with the physics packages that use the rapid radiative transfer model for general circulation models.The CESM2 outputs monthly mean BC concentration and instantaneous BC DRF every 73 time steps to reach a good balance of computational cost and data size (Conley et al 2013).
The uncertainties of modeling results mainly come from BC emissions, dry and wet depositions, and aging processes (Liu et al 2021).BC emission from biomass burning in CESM is treated as emission from the land surface by default.It should be noted that the increasing trend of plume height observed in some wildfire events might affect the release height of wildfires and alter BC emission strength (Wilmot et al 2022).However, such release height is not uniform and hence is not straightforward to be taken into consideration in a global BC emission inventory.It is also worth mentioning that the input BC inventory in CESM2 is from wildfire sources rather than agricultural sectors, so the confidence factor of the emission process employed during the calculation of the integrated confidence factor could be different.The PKU wildfire emission inventory was constructed based on the GFED3 and reanalysis of the tropospheric chemical composition (RETRO) subject to biomass burning reconstruction dataset over the past 40 years (Wang et al 2014).As uncertainty exists in the magnitude of biomass burning emission (Randerson et al 2012, Huang et al 2013, van der Werf et al 2017), the confidence of this factor should be taken into account.As van der Werf et al suggested, a 1σ of about 50% of uncertainty is assessed for GFED3, where small fires were excluded (2017).However, the confidence factor with a 95% confidence limit cannot be quantified due to large knowledge gaps in many other critical input data and parameters.We, therefore, still use the literature-recommended integrated confidence factor of 5.5 for the model outputs, which integrates the uncertainties from BC emissions, dry depositions, and aging processes.The confidence factor of the emission process is estimated to be 4 (Liu et al 2021).Since this modeling investigation focuses on the BC burden from wildfire emissions, we did not attempt to verify modeling results against field measurements of BC concentrations associated with all emission sources.The SPEI considers the temperature factor and introduces the influence of surface evaporation changes, which is more sensitive to the drought reaction caused by global temperature rise.Positive and negative values of SPEI suggest wetness and dryness, respectively.The selection of these two indices is also based on the consideration that they are both monthly-averaged gridded data (0.5 • × 0.5 • latitude by longitude) covering the entire globe with an extended period, which provides the long-term climate change trend, in line with our model outputs.

Mann-Kendall trend test and Theil-Sen slope
To detect the monotonic trends in time series, Mann-Kendall trend test (Mann 1945, Kendall 1949), a nonparametric test which has been widely applied in scientific research (Jolly et al 2015, Bari et al 2016, Sarkar et al 2018), is used in this study.The original Mann-Kendall test did not consider the seasonality and serial correlation of data (Hussain and Mahmud 2019).So, it is used here only to analyze the significance of the trends of different annual time series.
The Theil-Sen trend estimate, calculated by the median of slopes between all data pairs, is a robust estimate of the magnitude of linear trend (Theil 1949, Sen 1968, Hussain and Mahmud 2019).In the following context, the Theil-Sen slope with the significance of Mann-Kendall test higher than 0.95 is used as a measure of temporal trend.

Multidimensional ensemble empirical mode decomposition (MEEMD)
MEEMD based on the ensemble empirical mode decomposition (EEMD), has been widely used in climate research (Prasad et al 2007, Ji et al 2014, Ali et al 2019) to elucidate temporal trends for multidimensional grid datasets (Wu et al 2009, 2016Rehman and Mandic 2010).Each time series at a specific spatial grid cell, expressed by x(t), can be decomposed to a residual term R n (t) and a sum of several oscillatory components ∑ C j (t) (j = 1, 2, …, n), given by The residual term is a time-varying curve containing one extremum at most after removing intrinsic interannual or interdecadal variability, whose dimension is the same as the times series.Residuals from all grid cells enable us to analyze the large-scale temporal trend of a given meteorological variant.
In the EEMD calculation of both indices, the amplitude of the noise added to the data is 0.05, and the trial number is 400.The maximum number of intrinsic mode functions (IMFs) is taken as 5. Detailed information about MEEMD can be found in Ji et al (2014).

Wildfire BC emission and land cover change
Wildfire emission is a major component of biomass burning with natural interannual variations.The intensity of wildfire-induced BC emission density is associated with the severity of the wildfire.We examine the distribution of global wildfire-induced BC emissions to identify the key areas prone to wildfires.The distribution of mean BC wildfire emissions from 1981 to 2010 is shown in figure 1(a).Wildfire BC emissions in sub-Saharan Africa, northern Australia, mid-South America, Southeast Asia, and the northern boreal forest zone are significantly higher than other regions worldwide.The most significant increasing trend in wildfire BC emission, estimated by the Theil-Sen slope, occurs in the Amazon Basin (figure 1(b)), while wildfire BC emissions in Southeast Asia, Australia, southeast US, and northern boreal forest zone decreased notably during the three decades.Global wildfire BC emission between 1981 and 2010, as displayed in figure S1, ranged from 1.59 to 2.44 Tg yr −1 , which was close to global fire BC emission reported by GFED (Bond et al 2013, Dong et al 2019).It fluctuated annually without showing a significant increasing or decreasing trend from 1981 to 2010, which is also the case for wildfire BC emission in sub-Saharan Africa, where no significant temporal trend of wildfire BC emission can be seen in most places (figure 1(b)).On average, BC wildfire emissions in sub-Saharan Africa accounted for 62.4% of the global total during the study period, indicating its significant impact on non-anthropogenic BC emission.
Wildfire acts, to some extent, as a disturbance to land cover change (Ehrlich et al 1997) By comparing figures 1 and 2, we cannot discern much similarity in spatial distribution pattern between wildfire BC emission and the frequency of FSG cover change.For example, in Africa, where the land cover change in the FSG is insignificant, the wildfire-induced BC emission is the highest across the globe.As a result, the impact of land cover change on long-term wildfire-induced BC emission in Africa is negligible.On the contrary, in the Mediterranean region, where wildfire BC emission is relatively low, FSG transition is much more significant, suggesting the intensity of wildfire emission is not strongly associated with the frequency of land cover transition.This discrepancy can be attributed to different plant types, climate, topography, and management (Cochrane 2009).Also, it should be emphasized that fire is not the sole reason that causes land cover change.Anthropogenic activities such as deforestation and afforestation can also change FSG cover.

BC concentration from wildfire emission
The distribution of near-surface wildfire BC concentration, hereafter referred to as BC s , averaged from 1981 to 2010, is shown in figure 2(a).As seen, most wildfire BC emissions occurred in the areas to the north of 45 • S. Likewise, BC s are relatively high in sub-Saharan Africa, northern Australia, Amazon Basin, and Southeast Asia.The annual average BC s exceeded 100 ng kg −1 in sub-Saharan African regions, with significantly higher BC s levels extending from the latitude of 10 • N to 10 • S in Africa, suggesting heavy BC pollution from forest fires.Both wildfire BC emission and BC s show prominent increasing trend in Amazon Basin (figures 1(b) and 2(b)).Even though the northern boreal forest zone is identified as a crucial region where lightning-caused large wildfires account for more than 85% of the total area burned in this region (Macias Fauria and Johnson 2008), we can observe a decreasing trend of wildfire BC emissions and BC s from 1981 to 2010.

Response of BC concentration to climate indices
As aforementioned, the severity of warming and drought are crucial climate factors causing wildfires, so we examined the associations between BC s and climate indices by quantifying the correlations between gridded BC s and near-SAT T and SPEI (method).To distinguish the impact of climate change on wildfire BC concentration, MEEMD (method) is used to decompose the temporal trend of these two indices from 1901 to 2018.Before decomposition, we calculated the annual average of T and SPEI to filter out seasonal variations.The residuals of T and SPEI after the MEEMD are defined as R T and R SPEI , respectively.It has been suggested that the occurrence of wildfires depend on weather and climate change has a potential impact on it (Wotton et al 2017), so we analyzed the intra-annual and interannual relationship between wildfire BC emission and these two climate indices at each grid cell, where monthly data were used in the intra-annual regression model, and annual data were used in the interannual regression model.The R 2 , estimated by the ratio of the regression sum of squares to the total sum of squares, was used to evaluate the performance of regression models.For the regression model with annual data as input, the area-weighted average R 2 is 0.10.For the regression model with monthly data as input, the estimated average R 2 is 0.06.The result suggests that the two climate indices can better explain annual variations of wildfire BC.Moreover, as the main focus of this study is to analyze the impact of interannual climate change on wildfires, the effects of seasonal climate variability are not further discussed in the following context.
We then explored the interannual correlations between R T /R SPEI and BC s at different time lags.In general, the area-weighted correlation coefficients across the globe yield better results when the time lag is 0 (table S1).Even though all correlation coefficients are not very significant, it can be inferred that a rise in T, along with severer drought (decreasing SPEI), tends to raise the level of BC s .Figure 3 shows the correlation coefficients between R T /R SPEI and contemporary wildfire-induced BC s .Significantly positive correlations can be observed in most areas worldwide, especially South America (except for Bolivia), where both wildfire BC emissions (figure 1(a)) and T (figure S4) enhanced.The positive correlation between R T and BC concentration was also very prominent in regions like North Africa, Western and Central Asia, and Greenland, where local wildfire emissions were low.In contrast, the correlation coefficients between R T and BC s in sub-Saharan Africa and Australia, where wildfire emissions were relatively high, were relatively insignificant.The correlation between R SPEI and BC s varies in different regions, as shown in figure 3(b).Though the distribution of correlation coefficients is inhomogeneous with small-scale features, we can identify a prominent negative correlation between R SPEI and BC s in North Africa and Western Asia, where both SPEI and the difference in R SPEI are negative (figure S5).
The occurrence of wildfire is modulated by different meteorological factors, which may strengthen or counterbalance each other.For example, T in the northern boreal forest zone is inclined, which favors wildfires, but the concurrent increase in SPEI tends to counteract it.To quantify the effects of both R T and R SPEI , we established nine regression models with BC s as the dependent variables at each land grid.Ray et al have applied the same method to identify the contribution of temperature and precipitation variation to the crop yield variation (2015).Except for the R T and R SPEI terms, we also considered the R T and R SPEI squares and the interaction terms of R T and R SPEI as described in table 1, which can help us explain interactions among these variables from a linear regression perspective.To avoid overfitting of the model, we did not consider more combinations of R T and R SPEI .The regression period was set from 1981 to 2010 with annual data as input.The number of land grids with both available R T and R SPEI is 16 967.To determine whether the models are better than the null model with a random climate, F tests are conducted for each model at P = 0.10 level.Also, valid R 2 of the established models must fall into the range between [0, 1].If a model meets these two constraints, it can be regarded as a significant model.Global area-weighted R 2 and the number of land grid cells with significant models are presented in table 1.As seen, the number of grid cells with a statistically significant fit is smallest for Model 9 which has the greatest complexity of variables.Comparing Model 4 with Models 1, 2 and 3, the implement of the interaction term significantly enhanced the R 2 without large loss in valid grid number, suggesting the important composite impact of T and SPEI on BC s .
Then, regression analysis was conducted at each grid cell to find the best-fit model for different regions.The best-fit model defined here is referred to the significant regression model with the largest R 2 , and the number of regression model for selection at each land grid ranges from 0 to 9. Selected best-fit models for different regions and their R 2 are presented in figure 4. The best-fit models at 53.28% of land grid cells are better than the null model.Among all significant best-fit models, Model 8 accounts for the largest percentage (24.20%) of the selected models.Model 9 and Model 7 account for 22.61% and 20.29% of the best-fit models, respectively.Model 5 and Model 6 were not adopted.Of these models, Model 8 was mainly applied in Africa, South America, southern North America, South Asia, Southeast Asia, and Australia.Model 9 was primarily used in the northern boreal forest zone.The R 2 for the best-fit model is significantly larger in South America, high-latitude Asian regions, the southeastern US, and Australia than in the other regions of the world.The global area-weighted R 2 for the bestfit models is 0.33, which manifests that about one third of the variations in wildfire-induced BC concentration can be explained by the variations in near-SAT, precipitation, and evapotranspiration.If only linear regression models (Model 1, 2, and 3) are used for selection, the number of land grids with best-fit models accounts for only 36.05% of all land grids.As shown in figure S6(a), among all linear regression models, Model 3 yields a better prediction of BC concentrations at most land grids, thereby suggesting the importance of both R T and R SPEI used as indicators of BC s .The global area-weighted R 2 for the linear best-fit model is 0.17, about half of the best-fit R 2 for all nine regression models, indicating that further improvement in BC s prediction attributable to wildfire-induced emissions could be achieved by implementing more quadratic and interaction terms.

BC DRF caused by wildfire emission
The mean wildfire-induced BC DRF averaged over the study period is shown in figure 5(a).Similar to the spatial distribution of wildfire-related BC emissions (figure 1(a)) and BC concentrations (figure 2(a)), the largest wildfire-induced BC DRF can be observed in sub-Saharan Africa, with an average value greater than +800 mW m −2 (figure 5(a)).The line graphs in figure 5(a) shows that the amplitude of wildfire-induced BC DRF varies significantly among different regions and countries.While the regional wildfire-related BC DRF in sub-Saharan Africa exceeded +1000 mW m −2 in 1998, the global wildfire-induced BC DRF was less than +300 mW m −2 simultaneously.Also, the 30 year averaging wildfire BC DRF in sub-Saharan Africa (+860.87mW m −2 ) transcended the BC DRF in other regions and countries.Without the contribution of sub-Saharan Africa, the global wildfireinduced BC DRF could drop by about one-fifth.This share was calculated by the relative change in global BC DRF with (+179.19mW m −2 ) and without the inclusion of sub-Saharan Africa (+146.65 mW m −2 ).It should be emphasized that BC DRF in the downwind areas of sub-Saharan Africa could be also affected by emissions from the continent, so the estimated contribution from wildfire emission in sub-Saharan Africa to global BC DRF reduction was likely underestimated.The result suggests that wildfires in sub-Saharan Africa account for a large portion of the global wildfire-induced BC DRF.Also, we identify that wildfire-induced BC DRF in Brazil, Southeast Asia, Australia, Canada, and Russia are significantly higher than the global average value.In contrast, relatively lower BC DRF can be seen in Europe, India, and the US.Despite the significant fluctuation in global area-weighted wildfire-induced BC DRF, we observe an increasing trend during the two decades with the p-value of Mann-Kendall trend test less than 0.05, suggesting an enhanced climate effect caused by wildfire BC emissions.
The Theil-Sen slope of BC DRF between 2010 and 1981, as displayed in figure 5(b), is also most evident in the Amazon Basin, where a prominent positive trend can be identified.Annual BC DRF in Brazil, illustrated by the line graph in figure 5(a), also shows a significant increasing trend.In tropical Africa and the Arctic, wildfire-induced BC DRF inclined.However, there is no evident change in wildfireinduced BC DRF in vast regions like the northern boreal forest zone, the US, and Australia, suggesting a relatively small impact of wildfire-induced BC emissions in these regions on global climate change.

Discussion
Given the significant decreasing BC emissions from fossil fuel combustion and other anthropogenic sources subject to worldwide emission control strategies, BC emissions from wildfire biomass burning and other natural sources are expected to account for a larger portion in total BC emissions.The present study assessed global BC concentrations' temporal and spatial evolution and climate forcing sourced in wildfires from 1981 to 2010.We first tried to link wildfires with the land cover change of FGS.However, insignificant land cover change could be found in sub-Saharan Africa where wildfire-induced BC emission was the most prominent worldwide.We then simulated BC concentrations and radiative forcing during this period across the globe.The modeling results provided a long-period (30 years) gridded concentration dataset that enabled us to establish statistically significant relationships between wildfire-emitted BC concentrations and climate variables contributing to wildfire biomass burning.Using multivariate regression models, we identified that about 33% variation in wildfire BC concentration was attributable to the changes in near-SAT, precipitation, and evapotranspiration in 53% of the land grid cells around the world, implying that global climate change could drive BC emissions, environmental contamination, and radiative forcing from wildfire biomass burning.
Although past warming seems not to drive the enhancement of wildfires significantly from a global perspective, increasing wildfire frequency and intensity were found to respond to a warming climate in some climate and ecologically sensitive regions, such as the Arctic (Luo et al 2020).Global wildfireassociated BC DRF shows a significant increasing trend during the study period.Our result in the present study revealed that sub-Saharan Africa is a critical area with significant wildfire-induced BC emission, concentration, and DRF.The social and ecological function of tropical rainforests in Africa is irreplaceable.However, the proneness of wildfire makes it hard to preserve the tropical rain forests where the social economy is among the less developed.
Additionally, the significant emission caused by wildfires in tropical Africa can threaten our efforts to mitigate air pollution and global climate change.Collective efforts from other countries may be needed to guide and support native African people in fire management for the good of the world (Moreira et al 2012, Moritz et al 2014).Concerns have been also raised about the higher risks of escaping fires resulted from deforestation and other agricultural practices such as crop residue burning under global climate change (Brando et al 2020).The increasing deforestation rate across the Amazon Basin and other regions in the globe in recent years (Mataveli et al 2021) could potentially play an increasingly important role in growing BC emissions from wildfire biomass burning.These should be taken into consideration in the future investigation.
Here we select two climate indices, the global nearsurface air temperature(SAT) (T) from the Climatic Research Unit (Mitchell and Jones 2005, Harris et al 2020) and the standardized precipitation evapotranspiration index (SPEI) (Vicente-Serrano et al 2010, Svoboda and Fuchs 2017), to quantify the responses of wildfire-induced BC concentrations to climate.
The differences in R T /R SPEI between 2010 and 1981 can then be considered indicators of climate change.The spatial distribution of mean T and R T from 1981 to 2010 are shown in figures S4(a) and (b).Both T and R T to the south of 30 • N and the Southern Hemisphere are mostly positive, with large values of T and R T centralized in tropical and subtropical regions and negative values of T and R T across the Arctic.Likewise, the Theil-Sen slope of T and the difference in R T between 2010 and 1981 show a similar pattern.The most prominent positive trend of T occurred in the Arctic, high latitudes of North America, and the regions around 30 • N, featured by the downdraft of the Hadley Cell.Large-scale positive differences of R T can be identified in Central Asia and Canada, in line with areas with prominent wildfires.Also, the decreasing T measured by both the Theil-Sen slope of T and the difference of R T are noteworthy in Bolivia.Nevertheless, we can identify non-significant Theil-Sen slope of T and a positive difference in R T in mid-Russia at the same time.Figures S5(a) and (b) illustrate the mean SPEI and R SPEI from 1981 to 2010.Compared to zonally distributed T and R T , SPEI, and R SPEI seem relatively spatially inhomogeneous.Distinct differences in their spatial distribution can be identified in North Africa.The significant negative values of SPEI in Africa and Western Asia indicate high severity of drought in these regions.The positive SPEI in the northern boreal forest zone suggest relatively higher wetness which is not conducive to the occurrence of wildfires.The Theil-Sen slope of SPEI and the difference

Figure 2 .
Figure 2. Mean wildfire-induced near-surface BC concentrations (ng kg −1 ) from 1981 to 2010 and Theil-Sen slope of BC concentration from 1981 to 2010.(a) Mean BC concentrations averaged from 1981 to 2010; (b) Theil-Sen slope of near-surface BC concentration from 1981 to 2010.

Figure 3 .
Figure 3. Correlation coefficients between R T and RSPEI and wildfire-induced BC concentrations from 1981 to 2010.(a) Correlation coefficients between R T and BC concentrations from 1981 to 2010; (b) correlations between RSPEI and BC concentrations from 1981 to 2010.

Figure 4 .
Figure 4. Selected best-fit models and their R 2 .(a) Best-fit models.The subscript 1-9 of R 2 below the color bar of figure (a) represents the serial number of the selected model at a given model grid cell; (b) R 2 , estimated by the ratio of the regression sum of squares to the total sum of squares.Area-weighted R 2 across the globe with the best-fit significant model is 0.33.

Figure 5 .
Figure 5. Mean wildfire-related BC DRF (mW m −2 ) averaged from 1981 to 2010 and its Theil-Sen slope between 1981 and 2010.(a) Mean wildfire-induced BC DRF across the globe and temporal variations of annual wildfire-induced BC DRF in selected regions and countries.The colors of solid lines in line charts are based on the color bar of figure (a); (b) Theil-Sen slope of wildfire-induced BC DRF between 2010 and 1981.
. Though wildfire does not necessarily lead to a change in land cover, we can identify the hotspots of land cover change in savannas and forests.We use HIstoric Land Dynamics Assessment + (HILDA +) land use and land cover data to estimate the annual land cover change in those areas with forests, grasslands, and shrublands cover (FGS) (Winkler et al 2021).The frequency of land cover changes at one spatial grid from 1981 to 2010, estimated by the times of annual transition of FGS into other land types, is shown in figure S2.The map of land cover in 1981 is shown in figure S3.The most evident and frequent FSG cover changes occurred in the Mediterranean, where shrub encroachment and forest fires have been widely reported (Pinol et al 1998, Cristofanelli et al 2009, Athanasopoulou et al 2014, Bekar and Tavs ¸anoǧlu 2017, Bowman et al 2017, Winkler et al 2021).Land cover change in India, Australia, northern boreal forest zone, forest zone in Africa, and central South America is also very prominent.
Stevenson S, Swain D L, Singh D, Kalashnikov D A and Huang X 2022 Climate change increases risk of extreme rainfall following wildfire in the western United States Sci.Adv. 8 1-12 Urbanski S P 2013 Combustion efficiency and emission factors for wildfire-season fires in mixed conifer forests of the northern Rocky Mountains, US Atmos.Chem.Phys.13 7241-62 Urbanski S P, Hao W M and Nordgren B 2011 The wildland fire emission inventory: Western United States emission estimates and an evaluation of uncertainty Atmos.Chem.Phys.11 12973-3000 van der Werf G R et al 2017 Global fire emissions estimates during 1997-2016 Earth Syst.Sci.Data 9 697-720 Veira A, Lasslop G and Kloster S 2016 Wildfires in a warmer climate: emission fluxes, emission heights, and black carbon concentrations in 2090-2099 J. Geophys.Res.121 3195-223 Vicente-Serrano S M, Beguería S and López-Moreno J I 2010 A multiscalar drought index sensitive to global warming: the standardized precipitation evapotranspiration index J.Clim.23 1696-718 Wang R et al 2014 Trend in global black carbon emissions from 1960 to 2007 Environ.Sci.Technol.48 6780-7 Ward D S, Kloster S, Mahowald N M, Rogers B M, Randerson J T and Hess P G 2012 The changing radiative forcing of fires: global model estimates for past, present and future Atmos.Chem.Phys.12 10857-86 Westerling A L and Bryant B P 2007 Climate change and wildfire in California Clim.Change 87 231-49 Westerling A L, Bryant B P, Preisler H K, Holmes T P, Hidalgo H G, Das T and Shrestha S R 2011 Climate change and growth scenarios for California wildfire Clim.Change 109 445-63 Whitman E, Parisien M A, Thompson D K and Flannigan M D 2019 Short-interval wildfire and drought overwhelm boreal forest resilience Sci.Rep. 9 1-12 Williams A P et al 2019 Observed impacts of anthropogenic climate change on wildfire in California Earth's Future 7 892-910 Wilmot T Y, Mallia D V, Hallar A G and Lin J C 2022 Wildfire plumes in the Western US are reaching greater heights and injecting more aerosols aloft as wildfire activity intensifies Sci.Rep. 12 1-14 Winkler K, Fuchs R, Rounsevell M and Herold M 2021 Global land use changes are four times greater than previously estimated Nat. Commun.12 1-10 Wotton B M, Flannigan M D and Marshall G A 2017 Potential climate change impacts on fire intensity and key wildfire suppression thresholds in Canada Environ.Res.Lett.12 095003 Wu Z, Feng J, Qiao F and Tan Z M 2016 Fast multidimensional ensemble empirical mode decomposition for the analysis of big spatio-temporal datasets Phil.Trans.R. Soc.A 374 20150197 Wu Z, Huang N E and Chen X 2009 The multi-dimensional ensemble empirical mode decomposition method Adv.Adapt.Data Anal. 1 339-72 Wulder M A, White J C, Alvarez F, Han T, Rogan J and Hawkes B 2009 Characterizing boreal forest wildfire with multi-temporal Landsat and LIDAR data Remote Sens. Environ.113 1540-55 Yokelson R J, Karl T, Artaxo P, Blake D R, Christian T J, Griffith D W T, Guenther A and Hao W M 2007 The tropical forest and fire emissions experiment: overview and airborne fire emission factor measurements Atmos.Chem.Phys.7 5175-96 Yu P et al 2019 Black carbon lofts wildfire smoke high into the stratosphere to form a persistent plume Science 365 587-90