Unraveling the Impact of Land Cover Changes on Climate Using Machine Learning and Explainable Artificial Intelligence

Abstract

A general issue in climate science is the handling of big data and running complex and computationally heavy simulations. In this paper, we explore the potential of using machine learning (ML) to spare computational time and optimize data usage. The paper analyzes the effects of changes in land cover (LC), such as deforestation or urbanization, on local climate. Along with green house gas emission, LC changes are known to be important causes of climate change. ML methods were trained to learn the relation between LC changes and temperature changes. The results showed that random forest (RF) outperformed other ML methods, and especially linear regression models representing current practice in the literature. Explainable artificial intelligence (XAI) was further used to interpret the RF method and analyze the impact of different LC changes on temperature. The results mainly agree with the climate science literature, but also reveal new and interesting findings, demonstrating that ML methods in combination with XAI can be useful in analyzing the climate effects of LC changes. All parts of the analysis pipeline are explained including data pre-processing, feature extraction, ML training, performance evaluation, and XAI.

1. Introduction

One of the peculiar features of climate science is the accumulation of an enormous amount of data [1,2]. The estimated size of climate data exceeds ten petabytes and continues to grow exponentially [3]. Furthermore, the number of different and diverse data sources is also increasing. Initially, information is collected by thousands of ground-based weather instruments all over the world, such as weather stations, as well as by a large number of satellites that perform measurements from kilometers above the ground. These data need to be processed and transformed to formats that are comparable with each other. Climate science multimedia systems exists, but are still underinvestigated compared to other areas such as social media or medicine [4,5].

Motivated by the big data challenges and the need for multimedia systems within climate science, we will in this paper address the important, yet less studied, problem of analyzing the climate effects of land cover (LC) changes, such as deforestation or urbanization. We will go through all the steps of the multimedia pipeline to handle and gain knowledge from the large amounts of climate data to address this problem.

Nowadays, climate changes and global warming are indisputable facts [6,7,8,9,10]. Global surface temperature has been methodically collected since 1850. According to these records, the last 30 years exceed any previous decade in temperature. Furthermore, in some regions, the temperature has been reconstructed over the last 1400 years, and the period between 1983 and 2012 was the warmest 30-year period during this time [6,7,10,11].

Climate changes and temperature growth have a huge impact on natural and anthropogenic systems on all continents and oceans: melting of snow and ice, sea level rise, decrease in fresh water volume and quality, changes in precipitation patterns, behavior alterations of marine organisms and animals, negative effects on agriculture, and many other effects [9].

Anthropogenic factors, such as CO${}_2$ emission, are considered as the main cause of global warming in the second half of the 20th century [9]. LC changes are also known to influence the regional climate because they alter biophysical mechanisms, such as evapotranspiration, albedo, and surface roughness [12,13,14]. LC type transformation has various causes. On the one hand, it can be caused by natural factors such as floods, sea level rise, or wildfires. Climate change is also causing LC changes, such as favoring tree expansion in mountain areas, or greening, or forest degradation [15,16,17]. Anthropogenic factors, such as deforestation or expansion of agricultural land, also have a significant and often dominant impact on an LC transition. LC plays a significant role in energy and water exchange between atmosphere and the Earth’s surface. The terrestrial areas not only produce the greenhouse gases (such as CO${}_2$), but also absorb them [7]. Therefore, sustainable land management is an important tool for climate change mitigation. The Intergovernmental Panel on Climate Change (IPCC) [7] states in a recent report that the development of appropriate policies can considerably contribute to the climate change adaptation and affect the rate of temperature rise. Some of mechanisms that have already been implemented, confirm the efficiency of this approach [7]. Good examples of these measures are sustainable food production and forest management, food waste reduction, and avoidance and prevention of deforestation and land degradation. Even more political actions can be adopted. Nevertheless, in order to develop efficient policies, it is important to understand how different changes in LC affect local and global climate [18]. Researchers pay special attention to the importance of long-term monitoring of various types of LC transformations and their relation to climate changes [7,9,14,19].

Nowadays, simulations based on numerical climate models are the largest source of climate information [3]. They allow researchers to model a climate response to some specific changes in a climate system. To perform an experiment, one should run a climate model with different input variables a few times, and then compare the results to understand the impact of these input parameters. For instance, a climate model can show what kind of changes occurs in a climate system if the input data differs only in LC. Nevertheless, a result obtained from the climate model can have complex non-linear patterns that are difficult to identify, but can help to get new insight and spare costly simulations. For example, machine learning (ML) is of special interest to researchers as a powerful tool for such kind of tasks, as well as for other problems within climate science. However, while ML is widely used in different scientific fields, it still has a limited application within climate science [20].

There are many ways to study the impact of LC transformations on temperature. Some works on LC changes and their impact on climate are based on directly observed data [21,22], while others use mainly simulations and climate models [12,23,24]. In general, scientists agree that the climate system is very complex and depends on many factors. The impact of LC change can vary on global and regional scale. Moreover, the same transformations can lead to different consequences depending on the region where it happened. Obviously, different LC changes have a unique effect on climate [25]. Most publications are focused on some individual LC changes, for example, deforestation [24,26] or urbanization effects [27,28]. However, this question is rarely studied in a broad perspective, taking into account all types of LC transitions simultaneously. Recently, Huang et al. analyzed the regional impact of cumulative LC changes on European climate [14]. The key point was to take into account all types of LC simultaneously and further to distinguish the individual impact of different LC changes in regional climate.

In this paper, we follow the same approach as Huang et al. but given the aforementioned complexities in how LC changes affect temperature, we will explore the potential of using well working ML methods, such as support vector regression (SVR), random forest (RF), multiple linear regression (MLR), and least absolute shrinkage and selection operator regression (LASSO) to learn these complex relations [29]. The method that learns the relations best, will further be used to study the effects of LC changes on temperature, using a new suggested framework based on explainable artificial intelligence methods (XAI) [30].

2. Data

The dataset used in this paper is the same dataset as used in [14]. We provide a brief description of the most important properties of the data, and refer to [14] for further details. The dataset consists of two parts: (1) land cover data described in Section 2.1 and (2) temperature data described in Section 2.2. The dataset approximately covers Europe from about 22° W to 45° E longitude and from 27° N to 72° N latitude [31]. The data have a resolution of 467 cells in the south–north direction and 479 cells in the west-east direction. By excluding grid cells over water, the dataset consists of a total of 121,849 grid cells.

2.1. Land Cover Dataset

The European Space Agency (ESA) has produced detailed global LC maps for the period from 1992 to 2015 as a part of the Climate Change Initiative (CCI) [32]. These maps have a spatial resolution of 0.002778 degree (around 300 m at the equator) at latitude and longitude directions, and they contain 37 LC classes following the United Nations LC Classification System (UNLCCS) [33]. To obtain the dataset used in this paper, the 37 UNLCCS LC classes were transformed to the IGBP-MODIS land cover classification system following the cross-walking table given by Huang et al. [14]. The IGBP-MODIS system consists of 21 categories that are described in

Table 1. IGBP-MODIS classification system.

Land Cover Categories
Evergreen Needleleaf Forest	Cropland
Evergreen Broadleaf Forest	Urban and Built-Up
Deciduous Needleleaf Forest	Cropland/Natural Vegetation Mosaic
Deciduous Broadleaf Forest	Snow and Ice
Mixed Forest	Barren or Sparsely Vegetated
Closed Shrublands	Water
Open Shrublands	Wooden Tundra
Woody Savannas	Mixed Tundra
Savannas	Barren Tundra
Grassland	Lake
Permanent Wetland

. The spatial resolution of the re-classified data was further aggregated to a 0.11 degree (approximately 12 km at the equator) to agree with the climate model simulations described below. Each cell of the aggregated LC dataset contains information about the portion of each of the 21 LC classes.

In the period from 1992 to 2015, some categories of LC underwent more substantial changes than others. Figure 1 shows the most prominent LC changes in the dataset, such as the expansion of urban and built-up cover and changes in Evergreen Needleleaf forest. Different colors represent the proportion of a certain LC in each cell on the grid.

2.2. Simulated Temperature Data

The Weather Research and Forecasting (WRF) model version 3.9.1 was used to make simulations based on the input data that include the LC data for 1992 and 2015 and the settings of the international Coordinated Regional Climate Downscaling Experiment (CORDEX) initiative (EURO-CORDEX) [14,34]. The WRF model is one of the most accurate models for region climate simulations and has been validated and widely used in Europe [34,35]. Average climate data for 24 years (from 1 January 1992 to 31 December 2015) were produced from the ERA-Interim data and used as initial and lateral boundary conditions [36]. The result of the WRF model simulations we focus on is the simulated 2-m air temperature in degrees Celsius for each day of the year [14]. The simulated temperatures are a result of two runs of the WRF model. In these two independent runs, the boundary conditions are the same, but the LC dataset is different: one is with the input data of LC in 1992 and the other is with the input data of LC in 2015. Therefore, these simulations illustrate how the temperature would change if only LC changed. We refer to Huang et al. [14] for a detailed introduction of the setting of the simulations.

3. Machine Learning and Explainable Artificial Intelligence

Let $x_{i,j},i=1,\dots ,N,j=1\dots ,p$ represent N observations of p features, and let $y_i,i=1,\dots ,N$ represent some associated response. The aim of ML is to learn a function that is able to predict the response from these features. In this paper we consider four well-known models, namely MLR, LASSO, SVM, and RF.

MLR assumes a linear association between the features and the response

$$y_i={\beta }_0+{\beta }_1{x}_{i,1}+\cdots {\beta }_p{x}_{i,p}+{\epsilon }_i,i=1,\dots ,N$$

(1)

where ${\epsilon }_i$ represent zero mean Gaussian distributed error terms. The parameter estimates are usually found by minimizing the least squares error

$${\widehat{\beta }}_0,\dots ,{\widehat{\beta }}_p=\underset{{\beta }_0,\dots ,{\beta }_p}{arg\; min}\sum _{i=1}^N{\left(y_i-{\beta }_0-\sum _{j=1}^p{\beta }_j{x}_{i,j}\right)}^2$$

(2)

Given many features, a potential challenge with linear regression is that the model not only fits the signal in the data, but also the noise, usually resulting in poor prediction performance on held-out data. The problem is referred to as over-fitting. Regularization is a popular technique to address this issue. For example the LASSO model adds the sum of the absolute value of the parameter estimates as a penalty term to the optimization [37]

$${\widehat{\beta }}_0,\dots ,{\widehat{\beta }}_p=\underset{{\beta }_0,\dots ,{\beta }_p}{arg\; min}\sum _{i=1}^N{\left(y_i-{\beta }_0-\sum _{j=1}^p{\beta }_j{x}_{i,j}\right)}^2+\lambda \sum _{j=1}^p\left|{\beta }_j\right|$$

(3)

A positive property of the LASSO, is that the resulting model often will be sparse in the sense that most of the parameter estimates are set to zero, making model interpretation easier. A higher value of the regularization parameter $\lambda$ results in a more sparse solution, and less chance of over-fitting. In this paper, we adjusted the value of $\lambda$ to optimize prediction performance on held-out data.

We also consider two other very popular ML methods, namely the SVM and the RF. When the response is continuous (as it is in this work), SVM is often referred to as support vector regression (SVR). The idea behind SVR is to find the regression plane such that as many of the observations are within a (support) region around the regression plane as possible. The width of the support region is also part of the optimizing procedure.

The RF model consists of an ensemble of decision trees and, thus, is called random forest. A decision tree is a flowchart-like structure in which each internal node represents a decision based on a single feature or linear combination of a subset of features. The classification or prediction decision is based on a series of such individual decisions. RF is based on using different bootstrapping techniques to train multiple decision trees. RF makes decisions based on all the trees, for example through the average output from the trees or the majority output. In this paper, we based the decisions on the average outputs.

High dimensional data or a complex model, can make model interpretation difficult. Regression models can to some extent be interpreted by studying the size of the regression parameters ${\beta }_0,\dots ,{\beta }_p$, and represent the core of statistical inference. However, other models, such as the RF, are far more difficult to interpret. Recently the field of XAI has received a lot of interest trying to provide explanations for such opaque models. The core idea of XAI techniques is quite simple, and based on analyzing how changes in the input features affect the model output, but more sophisticated methods have also been developed [30]. In this paper, we will resort to a quite simple XAI approach based on analyzing how changes in a single feature will affect the output. The analysis will be explained in further detail in Section 4.2.

4. Experiments

In this section, we will describe the experiments with the aim of measuring the effects of LC changes on temperature. The section is organized as follows. In Section 4.1, we describe how to extract features from the dataset, and in Section 4.2, we describe our XAI-based method to analyze the effects of LC changes on temperature.

4.1. Feature Extraction

Our approach is based on using the difference in LC as features

$$x_{i,j}=L{C}_{i,j}^{2015}-L{C}_{i,j}^{1992},i=1,\dots ,N,j=1,\dots ,21$$

(4)

where $L{C}_{i,j}^{1992}$ and $L{C}_{i,j}^{2015}$ refer to the portion of LC type j in grid cell i in 1992 and 2015, respectively. Since the $x_{i,j}$’s represent differences between two portions, it follows that

$$\sum _{j=1}^{21}{x}_{i,j}=0\forall i$$

(5)

To be able to study the effects of LC changes on temperature we define the differences in average temperatures

$${\overline{y}}_{i,D}=\frac{1}{\left|D\right|}\sum _{d\in D}{T}_{i,d}^{2015}-T_{i,d}^{1992},i=1,\dots ,N$$

(6)

where N refers to the number of grid cells and $T_{i,d}^{1992}$ and $T_{i,d}^{2015}$ to the temperatures from the simulations described in Section 2.2 in grid cell i at day d. D refers to some part of the whole year, and $\left|D\right|$ the number of days in this period. In the computer experiments, five periods were used, namely winter (December, January, February), spring (March, April, May), summer (June, July, August), and autumn (September, October, November), and the whole year.

We will predict ${\overline{y}}_{i,D}$ using LC changes in the same geographic location in line with the recent literature [14,25]. For a given period D, the dataset used in the experiments therefore were as follows

$$\mathrm{Output}:\left[\begin{array}{c}{\overline{y}}_{1,D}\\ {\overline{y}}_{2,D}\\ ⋮\\ {\overline{y}}_{N,D}\end{array}\right],\mathrm{Input}:\left[\begin{array}{cccc}{x}_{1,1}& x_{1,2}& \dots & x_{1,21}\\ x_{2,1}& x_{2,2}& \dots & x_{2,21}\\ & & ⋮& \\ x_{N,1}& x_{N,2}& \dots & x_{N,21}\end{array}\right]$$

(7)

4.2. Analyzing Effects of LC Changes on Temperature Using XAI

In this section, we explain our approach to analyze the effects of LC changes on temperature. We suggest to use a XAI technique which is based on inserting different LC changes into to the trained models to study the resulting effects on temperature changes. We point out three important considerations when using this approach:

• To analyze the effect of some LC changes, we must check that the given change is frequently present in the training dataset to ensure that the ML method has learned the relation between this LC change and the temperature change well. If a LC change is not present in the training dataset, we cannot trust the model prediction related to this change. We will, therefore, only show the effects for the most frequent LC changes in the dataset.
• In statistical inference, multicollinearity is a well-known issue, resulting in a range of models that have about the same agreement with the data, but may represent different inferential conclusions. The XAI technique used in this paper has to tackle the same potential challenge. However, except from the "grand" collinearity in Equation (5), we did not observe any strong multicollinarity in the data, and we can, therefore, reliably use the suggested XAI approach.
• It is well-known that the association between LC changes and temperature is complex and associated with noise. It is, therefore, important to quantify the uncertainty of the ML prediction. We quantify uncertainty by using standard $(1-\alpha )·100\%$ model output prediction intervals

  $$\widehat{y}\pm z_{\alpha /2}\widehat{\sigma }$$

  (8)

  where $\widehat{y}$ is the model prediction, $z_{\alpha /2}$ the $\alpha /2$ quantile of the standard normal distribution and $\widehat{\sigma }$ the estimated standard deviation prediction error [38]. The standard deviation was estimated by the prediction error on unseen test examples in a ten fold CV experiment over the data samples in Equation (7).

5. Results

In this section, we summarize the results from the experiments described above. In Section 5.1, we compare the prediction performance of the different ML methods introduced in Section 3, and in Section 5.2, we show the effects of LC changes on temperature.

5.1. Evaluation of ML Methods to Predict Temperature Changes from LC Changes

In this section, we represent the performance of MLR, LASSO, RF, and SVR to predict temperature changes ${\overline{y}}_{i,D}$ from LC changes $x_{i,j}$. In addition, we evaluated the performance of a baseline predicting the temperature without using the LC features, i.e., the baseline predicted using the average temperature in the training data. The whole geographic area was divided into sectors, sizes 25 × 25, 50 × 50, and 75 × 75 cells. The methods were evaluated on a spatial cross validation (CV) procedure, where the methods were trained on data from all sectors except one. The remaining sector was used as a test set. This approach evaluates how well the methods can learn from some parts of the geographic area, and predict on others. We also evaluated the algorithms where the cells used for training and testing were randomly selected over the whole geographic area. The results from the sector and the randomization experiments were consistent, and only the results for the sector approach are shown. The prediction performance were measured using root mean squared error (RMSE) between the temperature differences ${\overline{y}}_{i,D}$ and model predictions. The results are shown in

Table 2. RMSE for the different methods. The values in boldface refers to the method that performed the best.

Sector Size	75 × 75	50 × 50	25 × 25
Baseline	0.1727	0.1730	0.1638
MLR	0.1713	0.1726	0.1618
LASSO	0.1727	0.1729	0.1638
RF	0.1642	0.1631	0.1511
SVM	0.1718	0.1745	0.1634

.

We see that RF outperforms the other ML methods, especially the linear regression models that represents the current practice in literature [14,25]. Using the 5 × 2 CV test [39], we verified that the RF performed significantly better than all the other algorithms with $p-\mathrm{values}<{10}^{-14}$. We further see that the other algorithms performed about equally well to the baseline and thus had trouble taking advantage of the information in the LC features. The performance of the methods are reduced with the larger sectors, and is as expected since the difference between the properties of the training and test areas are greater and the number of data samples used for training are reduced.

5.2. Effects of LC Changes on Temperature

The analyses are based on the RF model which had the best performance in our experiments. Since the RF model outperformed the linear regression models, that represent current practice [14,25], it is reasonable to assume that our analysis will be more robust and potentially reveal new insights from a climate science perspective. We will discuss this further in Section 6. The analyses are based on the suggested approach in Section 4.2.

Table 3. Temperature changes in the whole of Europe depending on LC transformation. Blue and red cells represent LC changes resulting in statistically significant cooling and warming, respectively. Intense blue and red cells show cells where LC changes resulted in at least a 0.5 °C temperature change with 95% certainty.

LC Change	Event Rate	Year, °C	Winter (DJF), °C	Spring (MAM), °C	Summer (JJA), °C	Autumn (SON), °C
Cropland to Urban and Built-Up	29,694	[0.0074, 0.6232]	[−0.3763, 0.3403]	[−0.743, 0.6816]	[0.131, 1.3258]	[−0.5429, 0.2748]
Cropland/Natural Vegetation Mosaic to Open Shrublands	21,539	[−0.2683, 0.3475]	[−0.5078, 0.2088]	[−0.9844, 0.4402]	[0.1809, 1.3757]	[−0.4492, 0.3685]
Evergreen Needleleaf Forest to Open Shrublands	16,561	[−0.3856, 0.2302]	[0.2147, 0.9313]	[−0.8066, 0.6179]	[−0.3429, 0.8519]	[0.2222, 1.0399]
Cropland to Open Shrublands	13,891	[−0.4164, 0.1994]	[−0.5177, 0.1988]	[−0.6879, 0.7367]	[0.1593, 1.3542]	[−0.4127, 0.405]
Cropland/Natural Vegetation Mosaic to Evergreen Needleleaf Forest	12,671	[−0.3734, 0.2425]	[−0.4416, 0.275]	[−0.8785, 0.546]	[0.0558, 1.2506]	[−0.5779, 0.2398]
Barren or Sparsely Vegetated to Urban and Built-Up	12,350	[−0.1896, 0.4262]	[−0.3071, 0.4095]	[−1.0188, 0.4057]	[−1.2128, −0.018]	[−0.6348, 0.1829]
Deciduous Broadleaf Forest to Urban and Built-Up	10,773	[0.0911, 0.7069]	[−0.433, 0.2836]	[−0.1756, 1.2489]	[0.0826, 1.2774]	[−0.2724, 0.5453]
Evergreen Needleleaf Forest to Urban and Built-Up	10,629	[0.0343, 0.6501]	[−0.1379, 0.5787]	[−0.444, 0.9805]	[−0.2144, 0.9804]	[−0.4064, 0.4113]
Deciduous Broadleaf Forest to Cropland	10,587	[−0.943, −0.3271]	[−0.3909, 0.3257]	[−1.3147, 0.1098]	[−1.7565, −0.5617]	[−0.8687, −0.051]
Grassland to Urban and Built-Up	10,448	[0.0068, 0.6227]	[−0.3413, 0.3753]	[−0.6201, 0.8045]	[−0.2194, 0.9754]	[−0.4165, 0.4012]
Barren or Sparsely Vegetated to Open Shrublands	10,131	[0.1269, 0.7427]	[−0.1964, 0.5202]	[1.1181, 2.5427]	[−0.7157, 0.4791]	[−0.5438, 0.2739]
Evergreen Needleleaf Forest to Cropland	9983	[−0.9488, −0.3329]	[−0.1522, 0.5643]	[−1.2604, 0.1642]	[−1.7591, −0.5643]	[−0.5507, 0.267]
Permanent Wetland to Open Shrublands	9935	[−0.1733, 0.4425]	[−0.4717, 0.2449]	[0.2164, 1.6409]	[−0.3666, 0.8282]	[−0.7281, 0.0897]
Mixed Forest to Open Shrublands	9320	[−0.2837, 0.3322]	[−0.617, 0.0996]	[0.1774, 1.602]	[−0.359, 0.8359]	[−0.627, 0.1907]
Cropland to Mixed Forest	8943	[−0.2976, 0.3182]	[−0.3925, 0.3241]	[−0.7132, 0.7113]	[0.4207, 1.6156]	[−0.6425, 0.1752]

shows the temperature change for the whole of Europe, if a grid cell completely changed from one LC to another. The table shows the 15 most frequent LC changes in the dataset ensuring reliability in our results (recall consideration 1. in Section 4.2).

Table A1. Temperature changes in northern Europe depending on LC transformation. Blue and red cells represent LC changes resulting in statistically significant cooling and warming, respectively. Intense red cells show cells where LC changes resulted in at least a 0.5 °C temperature change with 95% certainty.

LC Change	Event Rate	Year, °C	Winter (DJF), °C	Spring (MAM), °C	Summer (JJA), °C	Autumn (SON), °C
Permanent Wetland to Open Shrublands	8444	[−0.1177, 0.4338]	[−0.5883, 0.4403]	[0.3553, 1.6446]	[−0.5858, 0.9682]	[−0.7983, 0.1586]
Evergreen Needleleaf Forest to Open Shrublands	7728	[−0.1976, 0.3539]	[−0.5871, 0.4416]	[0.6168, 1.9061]	[−0.5377, 1.0164]	[−0.7309, 0.226]
Permanent Wetland to Evergreen Needleleaf Forest	6685	[−0.0996, 0.452]	[−0.5532, 0.4754]	[0.0165, 1.3058]	[0.1051, 1.6591]	[−0.6916, 0.2653]
Barren or Sparsely Vegetated to Open Shrublands	5860	[0.1398, 0.6913]	[−0.313, 0.7156]	[1.1442, 2.4334]	[−1.0866, 0.4674]	[−0.6445, 0.3124]
Permanent Wetland to Mixed Forest	5450	[−0.1356, 0.4159]	[−0.7163, 0.3123]	[0.4486, 1.7379]	[0.1063, 1.6604]	[−0.5375, 0.4194]
Evergreen Needleleaf Forest to Permanent Wetland	5191	[−0.1806, 0.3709]	[−0.775, 0.2537]	[−0.5402, 0.7491]	[0.1807, 1.7347]	[−0.661, 0.2959]
Mixed Forest to Open Shrublands	4935	[−0.2339, 0.3176]	[−0.4794, 0.5492]	[0.5399, 1.8292]	[−0.5467, 1.0074]	[−0.7691, 0.1877]
Permanent Wetland to Deciduous Broadleaf Forest	4648	[−0.166, 0.3855]	[−0.4265, 0.6021]	[0.174, 1.4633]	[−0.1017, 1.4523]	[−0.7823, 0.1746]
Deciduous Broadleaf Forest to Open Shrublands	4396	[−0.1111, 0.4405]	[−0.6223, 0.4063]	[0.5051, 1.7944]	[−0.7121, 0.842]	[−0.5298, 0.4271]
Cropland/Natural Vegetation Mosaic to Open Shrublands	4293	[−0.1426, 0.409]	[−0.4483, 0.5803]	[0.4013, 1.6906]	[−0.7162, 0.8378]	[−0.6488, 0.3081]

–

Table A3. Temperature changes in southern Europe depending on LC transformation. Blue and red cells represent LC changes resulting in statistically significant cooling and warming, respectively. Intense blue and red cells show cells where LC changes resulted in at least a 0.5 °C temperature change with 95% certainty.

LC Change	Event Rate	Year, °C	Winter (DJF), °C	Spring (MAM), °C	Summer (JJA), °C	Autumn (SON), °C
Barren or Sparsely Vegetated to Urban and Built-Up	5233	[−0.1294, 0.3696]	[−0.1071, 0.219]	[−1.0459, 0.0215]	[−0.0131, 0.7063]	[−0.1367, 0.4365]
Cropland to Urban and Built-Up	4473	[−0.1651, 0.3339]	[−0.0778, 0.2483]	[−0.8228, 0.2446]	[0.4718, 1.1912]	[−0.2072, 0.3661]
Barren or Sparsely Vegetated to Cropland	3455	[−0.6361, −0.1371]	[−0.1009, 0.2253]	[−0.8968, 0.1706]	[−0.9814, −0.262]	[−0.6312, −0.0579]
Cropland/Natural Vegetation Mosaic to Urban and Built-Up	3417	[−0.2611, 0.238]	[−0.1155, 0.2107]	[−0.9243, 0.1432]	[0.4579, 1.1773]	[−0.2176, 0.3557]
Cropland/Natural Vegetation Mosaic to Open Shrublands	3071	[−0.2362, 0.2628]	[−0.2812, 0.045]	[−0.7906, 0.2769]	[0.2316, 0.951]	[−0.405, 0.1683]
Cropland/Natural Vegetation Mosaic to Cropland	2029	[−0.3848, 0.1142]	[−0.1325, 0.1937]	[−1.1184, −0.051]	[−1.5166, −0.7972]	[−0.3396, 0.2337]
Barren or Sparsely Vegetated to Open Shrublands	1926	[−0.2398, 0.2593]	[−0.2134, 0.1128]	[−0.7118, 0.3556]	[0.0753, 0.7947]	[−0.3321, 0.2412]
Cropland/Natural Vegetation Mosaic to Evergreen Needleleaf Forest	1424	[−0.2675, 0.2315]	[−1.1594, −0.8333]	[−1.0964, −0.029]	[−0.2992, 0.4202]	[−0.525, 0.0483]
Cropland to Open Shrublands	1422	[−0.313, 0.1861]	[−0.2122, 0.1139]	[−1.1058, −0.0384]	[0.2722, 0.9916]	[−0.0905, 0.4828]
Cropland/Natural Vegetation Mosaic to Deciduous Broadleaf Forest	1412	[−0.4676, 0.0315]	[−0.1385, 0.1877]	[−1.1793, −0.1119]	[−0.4408, 0.2786]	[−0.4877, 0.0856]

in Appendix A, show results for the northern, central, and southern regions of Europe. The values in the parentheses in the tables, show 95% prediction intervals. Blue and red cells show LC changes resulting in statistically significant cooling and warming, respectively. Intense blue and red cells show cells where LC changes resulted in at least a 0.5 °C temperature change with 95% certainty. The column ’Event rate’ shows the number of cells where the LC change was observed in the dataset.

6. Discussion

Despite the substantial amount of uncertainty in the predictions, the results reveal several statistically significant temperature changes. They also show that the most frequent LC changes result in mainly warming in northern and central Europe and primarily cooling in the southern Europe. The most frequent LC changes are largely different for the different parts of Europe, which makes sense since the different parts of Europe mainly consist of different types of vegetation. However, for the LC changes that are frequent in more than one part of Europe, we observe a consistency in temperature change. For example cropland to urban built-up result in significant warming in all three parts of Europe and for the whole of Europe. There is also a consistency between seasons of the year in the sense that a LC change either results in warming or cooling for every season, and interestingly this observation was not detected by Huang et al. [14] with the regression based approach (there are no rows with both red and blue cells). For example, for the whole of Europe, deciduous broadleaf forest to cropland results in statistically significant cooling for both summer and autumn and no statistically significant warming (or cooling) for the other seasons.

To further verify the validity of our suggested approach, we now analyze how consistent our results are with other studies based on statistical approaches and climate model simulations. Many studies revealed a strong correlation between temperature increase and growth in shrub species [6,40,41,42,43]. Some of these researchers discussed the positive feedback loop when LC transitions affect climate, while temperature changes also influence LC transformation [40,43,44]. Firstly, a warming increases a spreading of shrublands. Then, LC transition to shrublands influences the energy exchange, increasing the absorption of solar radiation due to lower surface albedo. This, in turn, results in a temperature rise. However, it can be complicated to distinguish what is the main driver in this feedback loop. In this paper, we analyze only the impact of LC on temperature change, ignoring the effect of a warming on LC. We observed that transition to open shrublands alone leads to a temperature increase in northern and southern Europe.

Some works demonstrate that shrubland increase in Arctic can lead to an annual temperature increase [41,42,45], which is consistent with our own findings. However, most articles only consider the growth of shrubs and do not pay attention to the initial cover. Therefore, our approach can help in understanding how prominent is the effect of LC transformation to shrubs depending on the initial LC. For instance, the replacement of barren or sparsely vegetated cover to shrublands causes a more significant warming than a temperature rise associated with transition from permanent wetland to open shrublands.

Urbanization and its impact on temperature is another subject which draws the interest of climate scientists. In general, researchers conclude that the transition to urban and built-up covers causes a warming [7,14,46,47]. Indeed, we also observed that most of the LC changes to urban and built-up covers results in a temperature growth during the whole year, as well as seasonally.

Deforestation and its contribution to a temperature increase, is an important research subject that has been explored by many authors [14,48,49]. In this paper, we also observed a similar trend. Most LC changes associated with deforestation observed in our work lead to a significant temperature increase.

Afforestation is considered as a possible solution to the problem of the warming effect of deforestation because of its contribution to cooling [7,48,49]. In this paper, we detected such a trend in southern Europe where the shift from cropland or natural vegetation mosaic to Evergreen Needleleaf or deciduous broadleaf forest results in a significant cooling. However, in central Europe, we could not identify a clear pattern in temperature change associated with afforestation. Moreover, the transition from permanent wetland to any kind of forest contributes to a warming in northern Europe. This is consistent with the results of Li et al. where a transition of any LC to forest leads to a cooling in tropical regions but to warming in high latitudes [49].

Summarizing, we can conclude that our predictions of the LC-change-impact on temperature are consistent with the main trends described by the IPCC [6,7] and other studies. Our analyses also revealed new insights which supports the assumption that the ML techniques can be a useful tool in climate science, and it is possible to develop a model that can make a meaningful prediction. In addition, our approach allows us to extract more complex patterns and gain a more clear understanding of the effect of different LC transitions. This demonstrates that the ML techniques can help to figure out the effect of LC changes on surface temperature which opens up for a myriad of future work to explore and exploit this further.

7. Conclusions

In this paper, we have presented a framework based on ML and XAI to analyze the effects of LC changes on temperature. The results show that the RF model documented better prediction performance that linear regression based models, that is the current practice in the literature [14,25]. Our framework based on RF is able to find several statistically significant relations that align with other research. Our analyses also revealed new insight from a climate science perspective. For example the consistency between seasons.

We train models that predict temperature changes using LC change at the same geographic location as features. However, it is expected that temperature changes can also be affected by LC changes at other geographic locations. An interesting direction for future research is, therefore, to develop models to predict temperature using also LC changes from other geographic locations as features. This will, however, complicate the XAI analyses since temperature changes in the model now depend on LC change from multiple geographic locations. Another interesting direction is to analyze the effects of telecoupling, how LC changes in one place affect the climate in other locations.