# Quantifying Spatiotemporal Dynamics of Solar Radiation over the Northeast China Based on ACO-BPNN Model and Intensity Analysis

Authors: Xiangqian Li; Zhijun Tong; Enliang Guo; Xiaolong Luo
Journal: Advances in Meteorology
 (2017)

---

## Abstract

Reliable information on the spatiotemporal dynamics of solar radiation plays a crucial role in studies relating to global climate change. In this study, a new backpropagation neural network (BPNN) model optimized with an Ant Colony Optimization (ACO) algorithm was developed to generate the ACO-BPNN model, which had demonstrated superior performance for simulating solar radiation compared to traditional BPNN modelling, for Northeast China. On this basis, we applied an intensity analysis to investigate the spatiotemporal variation of solar radiation from 1982 to 2010 over the study region at three levels: interval, category, and conversion. Research findings revealed that (1) the solar radiation resource in the study region increased from the 1980s to the 2000s and the average annual rate of variation from the 1980s to the 1990s was lower than that from the 1990s to the 2000s and (2) the gains and losses of solar radiation at each level were in different conditions. The poor, normal, and comparatively abundant levels were transferred to higher levels, whereas the abundant level was transferred to lower levels. We believe our findings contribute to implementing ad hoc energy management strategies to optimize the use of solar radiation resources and provide scientific suggestions for policy planning.


---

## Body

## 1. Introduction

Solar radiation is one of the most renewable and sustainable energy sources for the Earth’s ecosystem, affecting nearly all physical, chemical, and biological processes on the planet [1–3]. Its spatiotemporal variation on the Earth’s surface directly determines the formation and evolution of climate [4, 5]. Therefore, complete and accurate information on the spatiotemporal dynamics of solar radiation has been the focus of global climate change studies during the past several decades [6–8].However, solar radiation observation sites are relatively sparse compared to conventional weather stations that measure conventional meteorological parameters (i.e., temperature and precipitation); this may be attributed to the equipment and maintenance costs [9]. This imposes a significant issue and serious hindrance to forming a better understanding of the variations in solar radiation, especially in developing countries. To that end, high precision models for simulating solar radiation at low costs are imperative, and utilizing easily accessible meteorological variables and available tools is an effective way to achieve this goal. According to the literature, we can summarize four main types of models for estimating solar radiation: (1) radiative transfer models, (2) climatic empirical models, (3) remote sensing-based models, and (4) artificial intelligence models. Radiative transfer models take into account the effects of physical factors (e.g., scattering of atmospheric molecules and aerosols and absorption of carbon dioxide and ozone absorption) and meteorological factors (e.g., cloud amount and sunshine hours) on solar radiation, which make for strong physical attribute characterization architecture. Spectral radiation models and broadband radiation models are representative theoretical models [10–20]. Nevertheless, the structures of these classical models are fairly complex. In addition, some input parameters (e.g., ozone thickness and atmospheric precipitable water) are difficult to acquire at large scales. These disadvantages confine the application of theoretical models. The climate empirical models mainly utilize conventional weather observation data such as sunshine percentage, cloud amount, and temperature to establish empirical formulas, which have been extensively used in simulating solar radiation due to the convenient calculations and relatively high precision. Angstrom-Prescott and Hargreaves-Samani models are common climate empirical models [21–23]; however, these are difficult to apply universally as their empirical coefficients usually change with time and place. Artificial intelligence models, in particular the widely used back propagation (BP) neural network (BPNN), have made significant contributions to the study of solar radiation simulations [24–28]. BPNN has the ability of self-organization, self-adaptation, and self-learning, which avoids the complex process of analyzing data and constructing model as it summarizes the experience and reduces error from human interpretation [29]. Additionally, its ability to adaptively process information can be used to address the complex input-output nonlinear mapping problem, especially regarding fuzzy information problems that require simultaneous consideration of many factors and imprecise conditions [30]. The causes of solar radiation variations are quite complicated as there are many meteorological elements that contribute to its value [31–33]. In general, using the BPNN model to simulate solar radiation is a better method to overcome the existing problems. Nevertheless, the traditional BPNN model has its own limitations, such as a slower convergence rate and easily defaulting into a local optimal solution. Consequently, we propose ACO-BPNN, a modified BPNN model optimized by an ant colony algorithm, to overcome those disadvantages.To date, almost all studies on spatiotemporal solar radiation variation adopt climatological statistical diagnosis methods such as the empirical orthogonal function, wavelet analysis, and the Mann-Kendall test. These were designed to detect spatial patterns, periodic regularity, and abrupt point, respectively [34–36]. Although these methods could quantitatively analyze the spatiotemporal characteristics of solar radiation to some extent, they could not reveal the dynamic process of changes in solar radiation and reflect dynamic conversion regularity between different levels of solar radiation over the study region. Hence, it is necessary to introduce a new method to improve analysis for the spatiotemporal dynamics of regional solar radiation. In recent years, the intensity analysis method has been applied to land use and land cover change research, which can fully apply transfer matrix to analyze the intensity of land changes at three levels during different time periods: time interval, categories, and conversion [37]. There is no doubt that we will have a clearer understanding of the spatiotemporal dynamics of solar radiation if the precisely simulated characteristics of the ACO-BPNN model are combined with the dynamic analysis ability of intensity analysis.Northeast China (NEC), known as a region very sensitive to global climate change [38–40], has suffered enormous economic losses due to frequent and increasing natural disasters (i.e., floods, droughts, and soil erosion) [41–43]. Restoring vegetation and carrying out sustainable agriculture have become necessary for ecological remediation and the sustainable development of the economy in this region. Solar radiation is necessary for the growth of vegetation; therefore, it is critical to quantitatively analyze its spatiotemporal variation to implement proper and reasonable vegetation restoration plans and agricultural ecological patterns. However, the majority of previous studies regarding estimating solar radiation over NEC adopted climate empirical models and solar radiative models, and dynamic changes in solar radiation of this region were limited to traditional trend analyses that did not reflect the continuous, dynamic changes [44–46]. Therefore, the main objectives of this study were to (1) estimate solar radiation over NEC using the ACO-BPNN model and (2) apply the intensity analysis method to characterize the dynamic variations in solar radiation of this region.


## 2. Materials and Methods

### 2.1. Study Area Description

The NEC region, spanning approximately 1.24 million km2 (15.32% of the total area of China), and located between 38°42′–53°35′N and 115°25′–135°09′E (Figure 1), was selected as the study area of interest. It starts from the ChangBai Mountains in the east, reaches the Greater Khingan Mountains in the west, borders the Bohai Sea in the south, and extends to the Lesser Khingan Mountains in the north. These mountain ranges surround the central plains, including the Liao River Plain and the Songnen Plain. The Sanjiang Plain is located in the northeast corner of NEC, and the central plains constitute the Northeast Plain that is a significant marketable grain base. The NEC region lies in the east coast of Eurasia and thus is sensitive to the activity and intensity of a temperate monsoon climate. The annual precipitation varies from 400 mm in the northwest to 1,000 mm in the southeast and is concentrated mostly in summer (50–70%), with approximately equal amounts occurring in spring and autumn and the least amount occurring in the winter. The annual mean temperature ranges from −5.5°C in the northwest to 11.5°C in the southeast.Figure 1
Map of the study area, Northeast China.


### 2.2. Data Collection and Preparation

For a better understanding of the spatiotemporal dynamics of solar radiation, data on meteorological observations and solar radiation of NEC were collected from the China National Meteorological Centre (CNMC), including monthly sum values of air temperature(at), wind speed (ws), vapor pressure (vp), relative humidity (rh), cloud amount (ca), sunshine duration (sd), and daily temperature range (dtr). Data from 1982 to 2010 were initially included from 131 meteorological stations distributed over the entire NEC region. Solar radiation data from 18 radiation stations also adopted monthly sum values (MJ/m2). The longitude, latitude, and elevation above sea level of all stations were also recorded; as these stations were established at different years, some stations are missing data and therefore the records are of different lengths. To ensure consistency of the data series, meteorological data sets from a total of 114 stations, including 18 stations with radiation observations and 96 stations without them, were finally adopted for analysis, covering a period of 29 years from January 1982 to December 2010.The data sequence was reviewed to ensure that there were no missing data, and three radiation stations (SuoLun, TongLiao, and ChangChun) were chosen as the verification points. Therefore, these stations were not involved in the establishment of the ACO-BPNN model.


### 2.3. Models

#### 2.3.1. Overview of the BPNN Model

The BPNN is named after the back propagation (BP) algorithm, which is a typical multilayer network including the input layer, hidden layer, and the output layer. Hecht Nielsen proved that the 3-layer feedforward network with one hidden layer could approximate any multivariate function [47]. Therefore, BPNNs usually adopt 3-layer structures, and higher fitting precision can be obtained by adding the number of neurons in each layer. Moreover, the model adopts full connection between layers, and there is no connection between the same layer neurons (Figure 2).Figure 2
Backpropagation neural network with three layers.The BP algorithm aims to minimize the total error of the network by feeding back the results of the learning to the hidden layer neurons to adjust the weights. The specific algorithm is as follows:(1)netj=∑iwijxi+Θj,yj=fnetj,where xi is the output value of the previous node i; wij is the connection weights; Θj is the threshold of unit j, generally between -1 and 1; and f usually adopts a Sigmoid-type activation function:(2)fa=11+e-a.For a training sample, the error between actual output and expected output is defined as follows:(3)E=12∑k=1ltk-ok2,where l is the number of output layer units; tk is the expected output of unit k of output layer; ok is the actual output of unit k of output layer.The time required to process all training samples is referred to as a cycle. When training a multilayer neural network, the BP algorithm often can make the errorE less than the specified threshold set by researchers after a number of cycles. At this point, the network achieves convergence and ends the iterative training process.


#### 2.3.2. Incorporating ACO into the BPNN

Due to simple principles and strong operability, the BP algorithm has been used extensively in a large number of research fields. However, the BP algorithm also revealed some disadvantages with its popularization. In essence, the BP algorithm is a gradient descent algorithm that will significantly reduce the speed of model convergence when the objective function approaches the optimal value. If the initial parameter values are not properly set, the algorithm easily falls into a local optimum. By contrast, the Ant Colony Optimization (ACO) algorithm simulates the mode of an ant colony in finding the optimal foraging path, which is characterized by an intelligent universal search optimization ability and a more robust performance [48]. Therefore, in this study, an ACO algorithm was chosen to offer several optimized groups of initial weight values for the BPNN, which is then further refined to the most appropriate weights. The BPNN not only maintains the nonlinear mapping ability of the network but also possesses the universal optimization ability.To apply the ACO algorithm to the network, the definition domains of all parameters consisting of weights and thresholds were separated into a series of discrete points. Consequently, each point was an optional value for the corresponding parameter. Once an ant was determined, it can only select a value for each parameter from numerous optional points and records its counter mark concurrently. In other words, each parameter required a corresponding pheromone table (Table1).Table 1
Pheromone table for each parameter.
w
i
Tag
1
2
⋯
n
Dividing calibration
a
1
a
2
⋯
a
n
Pheromone amount
τ
1
τ
2
⋯
τ
nIn Table1, wi is the ith parameter to be optimized, ai is the dividing calibrated value that can be seen as a point, τi represents the pheromone amount of point ai, and n is the number of optional points. So the domain was divided into n-1 parts.When thekth ant reached the parameter wi, it selected the point value according to the probability that can be calculated by the following formula:(4)Pki=τi∑1≤j≤mτj.After selecting values for all the parameters, thekth ant accomplishes a traversal. Then it returns to the original nest and updates the pheromone table concurrently through the following equation:(5)τi+1=1-ρτi+Δτi,where ρ is the evaporation coefficient of pheromone, which ranges from 0 to 1. Δτi refers to a pheromone increase, which can be obtained based on the following equation:(6)Δτi=QEr,where Q denotes the pheromone quantity, which is a constant. Er is the error between the actual output and expected output.When all ants finished selecting values, a better combination of parameters for neural network was obtained. However, it was determined that the output precision of the neural network was relatively low when these parameters were directly applied, which may have been attributed to the fact that the definition domain was split into a series of discrete points. Therefore, a further searching algorithm, BP, was required to improve accuracy and facilitate the neural network in achieving rapid convergence. However, the BP algorithm provided random values for initializing the neural network parameters, which increased the risk of being trapped into the local optimum. We took full advantage of the ACO to provide a more suitable group of initial weight values for the BPNN. As a result, the speed of convergence and output precision was advanced. The framework of the ACO-BPNN model is displayed in Figure3.Figure 3
The flow chart of the improved BPNN: ACO-BPNN model.


#### 2.3.3. Error Analysis Methods

In this study, the mean absolute percentage error (MAPE), mean absolute biased error (MABE), and root mean square error (RMSE) were employed to analyze the error between measured and simulated radiation values. MAPE, MABE, and RMSE most often reflected the comprehensive condition, long term fitting performance, and sensitivity and short-term fitting information of the model. The formulas to calculate them are as follows:(7)MAPE=1N∑i=1NSm-SrSr,MABE=1N∑i=1NSm-Sr,RMSE=∑i=1NSm-Sr2N,where Sm represents the simulated value of the monthly sum of solar radiation, Sr represents the actual value of monthly solar radiation, and N refers to the number of samples.


### 2.4. Spatial Interpolation Method

Spatial interpolation was performed to provide a monthly spatial distribution of the dynamic analysis for the solar radiation across the NEC region from 1982 to 2010. Here the inverse distance weighted (IDW) method was applied, as it yielded higher interpolation precision than other general methods when distributing the solar radiation spatial pattern in China [49]. The background theory and specific calculations steps of the IDW are detailed in the study by Goovaerts [50].


### 2.5. Intensity Analysis

The intensity analysis method developed by Aldwaik and Pontius [37] fully utilized a transition matrix to quantitatively analyze the intensity of land use and land cover changes at three distinct levels (time interval, land class hierarchy, and conversion) during different time periods [37]. Similarly, the solar radiation within the different levels in the study region was also usually dynamically changing during the different periods. The intensity analysis method was introduced to counter these dynamic changes and to address three problems on spatiotemporal variation in solar radiation. (1) Is the total annual change of solar radiation occurring quickly or slowly within a certain time period? Formulas (8) and (9) were used to answer this problem. (2) Based on the previous answer, does the solar radiation at different levels vary gently or actively? Formulas (10) were used to answer this problem. (3) Based on the two former answers, which conversion form dominated the process of mutual transition between different levels of solar radiation? Formulas (11)–(14) were used to answer this problem. By addressing these three problems, the overall tendency of solar radiation, the dynamic change regulations of solar radiation at different levels, and the mutual conversion regularity between different levels across NEC were revealed. The specific calculations of intensity analysis are shown in formulas (8)–(14):(8)U=∑t=1T-1∑i=1JCtij-Ctjj/∑j=1J∑i=1JCtijYT-Y1×100%,(9)St=∑j=1J∑i=1JCtij-Ctjj/∑j=1J∑i=1JCtijYt+1-Yt×100%,where U and St refer to the uniform line value for the whole time interval and the annual variation intensity for a certain time interval, respectively. J and T represent the solar radiation level and the number of points in time, respectively. The subscript indexes i and j henceforth stand for the solar radiation levels at the initial and final times for a certain time interval, respectively, and subscript index t henceforth indicates a certain point in time ranging from 1 to T-1. Yt henceforth denotes the year of t, and Ctij henceforth means the number of pixels transferred from the Yt year at level i to the point in time Yt+1 at level j.(10)Gtj=∑i=1JCtij-Ctjj/Yt+1-Yt∑i=1JCtij×100%,Lti=∑j=1JCtij-Ctii/Yt+1-Yt∑j=1JCtij×100%.In the above formulas,Gtj and Lti are the annual total increasing and decreasing intensity level j and level i during the time interval [Yt,Yt+1], respectively.(11)Wtn=∑i=1JCtin-Ctnn/Yt+1-Yt∑j=1J∑i=1JCtij-Ctnj×100%.Here,Wtn refers to the uniform intensity of conversion from non-n levels at time point Yt to level n during the time interval [Yt,Yt+1], and subscript index n refers to the solar radiation level converted from other levels.(12)Rtin=Ctin/Yt+1-Yt∑j=1JCtij×100%.In this formula,Rtin refers to the annual conversion intensity from level i to level n(i≠n) during the time interval [Yt,Yt+1].(13)Vtm=∑j=1JCtmj-Ctmm/Yt+1-Yt∑j=1J∑j=1JCtij-Ctim×100%.In formula (13), Vtm refers to the uniform intensity of conversion from level m at time point Yt+1 to non-m levels during the time interval [Yt,Yt+1], and subscript index m refers to the solar radiation level converting to other levels.(14)Qtmj=Ctmj/Yt+1-Yt∑j=1JCtij×100%.In the above,Qtmj refers to the annual conversion intensity from level m to level j(m≠j) during the time interval [Yt,Yt+1].In summation, the intensity analysis method was performed based on the distribution conditions of solar radiation at different levels at each point in time. Therefore, the classification criteria at all points in time must be uniform. In this study we used previous research to base the division of solar radiation resources across NEC [51, 52] and take the national current standard and specifications into consideration simultaneously. Finally, the solar radiation resources from 1982 to 2010 in NEC were divided into the four levels shown in Table 2, and the 1980s, 1990s, and 2000s were selected as the time points.Table 2
Solar radiation classification in NEC.
Levels
Name
Classification threshold (MJ⋅m−2⋅a−1)
I
Abundant
G
∗ ≥ 5200
II
Comparatively abundant
5000 ≤G < 5200
III
Normal
4800 ≤G < 5000
IV
Poor
G < 4800
N
o
t
e
∗. G represents the annual total solar radiation.


## 2.1. Study Area Description

The NEC region, spanning approximately 1.24 million km2 (15.32% of the total area of China), and located between 38°42′–53°35′N and 115°25′–135°09′E (Figure 1), was selected as the study area of interest. It starts from the ChangBai Mountains in the east, reaches the Greater Khingan Mountains in the west, borders the Bohai Sea in the south, and extends to the Lesser Khingan Mountains in the north. These mountain ranges surround the central plains, including the Liao River Plain and the Songnen Plain. The Sanjiang Plain is located in the northeast corner of NEC, and the central plains constitute the Northeast Plain that is a significant marketable grain base. The NEC region lies in the east coast of Eurasia and thus is sensitive to the activity and intensity of a temperate monsoon climate. The annual precipitation varies from 400 mm in the northwest to 1,000 mm in the southeast and is concentrated mostly in summer (50–70%), with approximately equal amounts occurring in spring and autumn and the least amount occurring in the winter. The annual mean temperature ranges from −5.5°C in the northwest to 11.5°C in the southeast.Figure 1
Map of the study area, Northeast China.


## 2.2. Data Collection and Preparation

For a better understanding of the spatiotemporal dynamics of solar radiation, data on meteorological observations and solar radiation of NEC were collected from the China National Meteorological Centre (CNMC), including monthly sum values of air temperature(at), wind speed (ws), vapor pressure (vp), relative humidity (rh), cloud amount (ca), sunshine duration (sd), and daily temperature range (dtr). Data from 1982 to 2010 were initially included from 131 meteorological stations distributed over the entire NEC region. Solar radiation data from 18 radiation stations also adopted monthly sum values (MJ/m2). The longitude, latitude, and elevation above sea level of all stations were also recorded; as these stations were established at different years, some stations are missing data and therefore the records are of different lengths. To ensure consistency of the data series, meteorological data sets from a total of 114 stations, including 18 stations with radiation observations and 96 stations without them, were finally adopted for analysis, covering a period of 29 years from January 1982 to December 2010.The data sequence was reviewed to ensure that there were no missing data, and three radiation stations (SuoLun, TongLiao, and ChangChun) were chosen as the verification points. Therefore, these stations were not involved in the establishment of the ACO-BPNN model.


## 2.3. Models

### 2.3.1. Overview of the BPNN Model

The BPNN is named after the back propagation (BP) algorithm, which is a typical multilayer network including the input layer, hidden layer, and the output layer. Hecht Nielsen proved that the 3-layer feedforward network with one hidden layer could approximate any multivariate function [47]. Therefore, BPNNs usually adopt 3-layer structures, and higher fitting precision can be obtained by adding the number of neurons in each layer. Moreover, the model adopts full connection between layers, and there is no connection between the same layer neurons (Figure 2).Figure 2
Backpropagation neural network with three layers.The BP algorithm aims to minimize the total error of the network by feeding back the results of the learning to the hidden layer neurons to adjust the weights. The specific algorithm is as follows:(1)netj=∑iwijxi+Θj,yj=fnetj,where xi is the output value of the previous node i; wij is the connection weights; Θj is the threshold of unit j, generally between -1 and 1; and f usually adopts a Sigmoid-type activation function:(2)fa=11+e-a.For a training sample, the error between actual output and expected output is defined as follows:(3)E=12∑k=1ltk-ok2,where l is the number of output layer units; tk is the expected output of unit k of output layer; ok is the actual output of unit k of output layer.The time required to process all training samples is referred to as a cycle. When training a multilayer neural network, the BP algorithm often can make the errorE less than the specified threshold set by researchers after a number of cycles. At this point, the network achieves convergence and ends the iterative training process.


### 2.3.2. Incorporating ACO into the BPNN

Due to simple principles and strong operability, the BP algorithm has been used extensively in a large number of research fields. However, the BP algorithm also revealed some disadvantages with its popularization. In essence, the BP algorithm is a gradient descent algorithm that will significantly reduce the speed of model convergence when the objective function approaches the optimal value. If the initial parameter values are not properly set, the algorithm easily falls into a local optimum. By contrast, the Ant Colony Optimization (ACO) algorithm simulates the mode of an ant colony in finding the optimal foraging path, which is characterized by an intelligent universal search optimization ability and a more robust performance [48]. Therefore, in this study, an ACO algorithm was chosen to offer several optimized groups of initial weight values for the BPNN, which is then further refined to the most appropriate weights. The BPNN not only maintains the nonlinear mapping ability of the network but also possesses the universal optimization ability.To apply the ACO algorithm to the network, the definition domains of all parameters consisting of weights and thresholds were separated into a series of discrete points. Consequently, each point was an optional value for the corresponding parameter. Once an ant was determined, it can only select a value for each parameter from numerous optional points and records its counter mark concurrently. In other words, each parameter required a corresponding pheromone table (Table1).Table 1
Pheromone table for each parameter.
w
i
Tag
1
2
⋯
n
Dividing calibration
a
1
a
2
⋯
a
n
Pheromone amount
τ
1
τ
2
⋯
τ
nIn Table1, wi is the ith parameter to be optimized, ai is the dividing calibrated value that can be seen as a point, τi represents the pheromone amount of point ai, and n is the number of optional points. So the domain was divided into n-1 parts.When thekth ant reached the parameter wi, it selected the point value according to the probability that can be calculated by the following formula:(4)Pki=τi∑1≤j≤mτj.After selecting values for all the parameters, thekth ant accomplishes a traversal. Then it returns to the original nest and updates the pheromone table concurrently through the following equation:(5)τi+1=1-ρτi+Δτi,where ρ is the evaporation coefficient of pheromone, which ranges from 0 to 1. Δτi refers to a pheromone increase, which can be obtained based on the following equation:(6)Δτi=QEr,where Q denotes the pheromone quantity, which is a constant. Er is the error between the actual output and expected output.When all ants finished selecting values, a better combination of parameters for neural network was obtained. However, it was determined that the output precision of the neural network was relatively low when these parameters were directly applied, which may have been attributed to the fact that the definition domain was split into a series of discrete points. Therefore, a further searching algorithm, BP, was required to improve accuracy and facilitate the neural network in achieving rapid convergence. However, the BP algorithm provided random values for initializing the neural network parameters, which increased the risk of being trapped into the local optimum. We took full advantage of the ACO to provide a more suitable group of initial weight values for the BPNN. As a result, the speed of convergence and output precision was advanced. The framework of the ACO-BPNN model is displayed in Figure3.Figure 3
The flow chart of the improved BPNN: ACO-BPNN model.


### 2.3.3. Error Analysis Methods

In this study, the mean absolute percentage error (MAPE), mean absolute biased error (MABE), and root mean square error (RMSE) were employed to analyze the error between measured and simulated radiation values. MAPE, MABE, and RMSE most often reflected the comprehensive condition, long term fitting performance, and sensitivity and short-term fitting information of the model. The formulas to calculate them are as follows:(7)MAPE=1N∑i=1NSm-SrSr,MABE=1N∑i=1NSm-Sr,RMSE=∑i=1NSm-Sr2N,where Sm represents the simulated value of the monthly sum of solar radiation, Sr represents the actual value of monthly solar radiation, and N refers to the number of samples.


## 2.3.1. Overview of the BPNN Model

The BPNN is named after the back propagation (BP) algorithm, which is a typical multilayer network including the input layer, hidden layer, and the output layer. Hecht Nielsen proved that the 3-layer feedforward network with one hidden layer could approximate any multivariate function [47]. Therefore, BPNNs usually adopt 3-layer structures, and higher fitting precision can be obtained by adding the number of neurons in each layer. Moreover, the model adopts full connection between layers, and there is no connection between the same layer neurons (Figure 2).Figure 2
Backpropagation neural network with three layers.The BP algorithm aims to minimize the total error of the network by feeding back the results of the learning to the hidden layer neurons to adjust the weights. The specific algorithm is as follows:(1)netj=∑iwijxi+Θj,yj=fnetj,where xi is the output value of the previous node i; wij is the connection weights; Θj is the threshold of unit j, generally between -1 and 1; and f usually adopts a Sigmoid-type activation function:(2)fa=11+e-a.For a training sample, the error between actual output and expected output is defined as follows:(3)E=12∑k=1ltk-ok2,where l is the number of output layer units; tk is the expected output of unit k of output layer; ok is the actual output of unit k of output layer.The time required to process all training samples is referred to as a cycle. When training a multilayer neural network, the BP algorithm often can make the errorE less than the specified threshold set by researchers after a number of cycles. At this point, the network achieves convergence and ends the iterative training process.


## 2.3.2. Incorporating ACO into the BPNN

Due to simple principles and strong operability, the BP algorithm has been used extensively in a large number of research fields. However, the BP algorithm also revealed some disadvantages with its popularization. In essence, the BP algorithm is a gradient descent algorithm that will significantly reduce the speed of model convergence when the objective function approaches the optimal value. If the initial parameter values are not properly set, the algorithm easily falls into a local optimum. By contrast, the Ant Colony Optimization (ACO) algorithm simulates the mode of an ant colony in finding the optimal foraging path, which is characterized by an intelligent universal search optimization ability and a more robust performance [48]. Therefore, in this study, an ACO algorithm was chosen to offer several optimized groups of initial weight values for the BPNN, which is then further refined to the most appropriate weights. The BPNN not only maintains the nonlinear mapping ability of the network but also possesses the universal optimization ability.To apply the ACO algorithm to the network, the definition domains of all parameters consisting of weights and thresholds were separated into a series of discrete points. Consequently, each point was an optional value for the corresponding parameter. Once an ant was determined, it can only select a value for each parameter from numerous optional points and records its counter mark concurrently. In other words, each parameter required a corresponding pheromone table (Table1).Table 1
Pheromone table for each parameter.
w
i
Tag
1
2
⋯
n
Dividing calibration
a
1
a
2
⋯
a
n
Pheromone amount
τ
1
τ
2
⋯
τ
nIn Table1, wi is the ith parameter to be optimized, ai is the dividing calibrated value that can be seen as a point, τi represents the pheromone amount of point ai, and n is the number of optional points. So the domain was divided into n-1 parts.When thekth ant reached the parameter wi, it selected the point value according to the probability that can be calculated by the following formula:(4)Pki=τi∑1≤j≤mτj.After selecting values for all the parameters, thekth ant accomplishes a traversal. Then it returns to the original nest and updates the pheromone table concurrently through the following equation:(5)τi+1=1-ρτi+Δτi,where ρ is the evaporation coefficient of pheromone, which ranges from 0 to 1. Δτi refers to a pheromone increase, which can be obtained based on the following equation:(6)Δτi=QEr,where Q denotes the pheromone quantity, which is a constant. Er is the error between the actual output and expected output.When all ants finished selecting values, a better combination of parameters for neural network was obtained. However, it was determined that the output precision of the neural network was relatively low when these parameters were directly applied, which may have been attributed to the fact that the definition domain was split into a series of discrete points. Therefore, a further searching algorithm, BP, was required to improve accuracy and facilitate the neural network in achieving rapid convergence. However, the BP algorithm provided random values for initializing the neural network parameters, which increased the risk of being trapped into the local optimum. We took full advantage of the ACO to provide a more suitable group of initial weight values for the BPNN. As a result, the speed of convergence and output precision was advanced. The framework of the ACO-BPNN model is displayed in Figure3.Figure 3
The flow chart of the improved BPNN: ACO-BPNN model.


## 2.3.3. Error Analysis Methods

In this study, the mean absolute percentage error (MAPE), mean absolute biased error (MABE), and root mean square error (RMSE) were employed to analyze the error between measured and simulated radiation values. MAPE, MABE, and RMSE most often reflected the comprehensive condition, long term fitting performance, and sensitivity and short-term fitting information of the model. The formulas to calculate them are as follows:(7)MAPE=1N∑i=1NSm-SrSr,MABE=1N∑i=1NSm-Sr,RMSE=∑i=1NSm-Sr2N,where Sm represents the simulated value of the monthly sum of solar radiation, Sr represents the actual value of monthly solar radiation, and N refers to the number of samples.


## 2.4. Spatial Interpolation Method

Spatial interpolation was performed to provide a monthly spatial distribution of the dynamic analysis for the solar radiation across the NEC region from 1982 to 2010. Here the inverse distance weighted (IDW) method was applied, as it yielded higher interpolation precision than other general methods when distributing the solar radiation spatial pattern in China [49]. The background theory and specific calculations steps of the IDW are detailed in the study by Goovaerts [50].


## 2.5. Intensity Analysis

The intensity analysis method developed by Aldwaik and Pontius [37] fully utilized a transition matrix to quantitatively analyze the intensity of land use and land cover changes at three distinct levels (time interval, land class hierarchy, and conversion) during different time periods [37]. Similarly, the solar radiation within the different levels in the study region was also usually dynamically changing during the different periods. The intensity analysis method was introduced to counter these dynamic changes and to address three problems on spatiotemporal variation in solar radiation. (1) Is the total annual change of solar radiation occurring quickly or slowly within a certain time period? Formulas (8) and (9) were used to answer this problem. (2) Based on the previous answer, does the solar radiation at different levels vary gently or actively? Formulas (10) were used to answer this problem. (3) Based on the two former answers, which conversion form dominated the process of mutual transition between different levels of solar radiation? Formulas (11)–(14) were used to answer this problem. By addressing these three problems, the overall tendency of solar radiation, the dynamic change regulations of solar radiation at different levels, and the mutual conversion regularity between different levels across NEC were revealed. The specific calculations of intensity analysis are shown in formulas (8)–(14):(8)U=∑t=1T-1∑i=1JCtij-Ctjj/∑j=1J∑i=1JCtijYT-Y1×100%,(9)St=∑j=1J∑i=1JCtij-Ctjj/∑j=1J∑i=1JCtijYt+1-Yt×100%,where U and St refer to the uniform line value for the whole time interval and the annual variation intensity for a certain time interval, respectively. J and T represent the solar radiation level and the number of points in time, respectively. The subscript indexes i and j henceforth stand for the solar radiation levels at the initial and final times for a certain time interval, respectively, and subscript index t henceforth indicates a certain point in time ranging from 1 to T-1. Yt henceforth denotes the year of t, and Ctij henceforth means the number of pixels transferred from the Yt year at level i to the point in time Yt+1 at level j.(10)Gtj=∑i=1JCtij-Ctjj/Yt+1-Yt∑i=1JCtij×100%,Lti=∑j=1JCtij-Ctii/Yt+1-Yt∑j=1JCtij×100%.In the above formulas,Gtj and Lti are the annual total increasing and decreasing intensity level j and level i during the time interval [Yt,Yt+1], respectively.(11)Wtn=∑i=1JCtin-Ctnn/Yt+1-Yt∑j=1J∑i=1JCtij-Ctnj×100%.Here,Wtn refers to the uniform intensity of conversion from non-n levels at time point Yt to level n during the time interval [Yt,Yt+1], and subscript index n refers to the solar radiation level converted from other levels.(12)Rtin=Ctin/Yt+1-Yt∑j=1JCtij×100%.In this formula,Rtin refers to the annual conversion intensity from level i to level n(i≠n) during the time interval [Yt,Yt+1].(13)Vtm=∑j=1JCtmj-Ctmm/Yt+1-Yt∑j=1J∑j=1JCtij-Ctim×100%.In formula (13), Vtm refers to the uniform intensity of conversion from level m at time point Yt+1 to non-m levels during the time interval [Yt,Yt+1], and subscript index m refers to the solar radiation level converting to other levels.(14)Qtmj=Ctmj/Yt+1-Yt∑j=1JCtij×100%.In the above,Qtmj refers to the annual conversion intensity from level m to level j(m≠j) during the time interval [Yt,Yt+1].In summation, the intensity analysis method was performed based on the distribution conditions of solar radiation at different levels at each point in time. Therefore, the classification criteria at all points in time must be uniform. In this study we used previous research to base the division of solar radiation resources across NEC [51, 52] and take the national current standard and specifications into consideration simultaneously. Finally, the solar radiation resources from 1982 to 2010 in NEC were divided into the four levels shown in Table 2, and the 1980s, 1990s, and 2000s were selected as the time points.Table 2
Solar radiation classification in NEC.
Levels
Name
Classification threshold (MJ⋅m−2⋅a−1)
I
Abundant
G
∗ ≥ 5200
II
Comparatively abundant
5000 ≤G < 5200
III
Normal
4800 ≤G < 5000
IV
Poor
G < 4800
N
o
t
e
∗. G represents the annual total solar radiation.


## 3. Results

### 3.1. Comparative and Validating Analysis of ACO-BPNN Model

A total of 5220 input/output pairs were employed to train the BPNN and ACO-BPNN to fully establish the relationship between the input vectorX=[at,ws,vp,rh,ca,sd,dtr]T and the output vector Y=[solar  radiation]T at 15 stations from 1982 to 2010. Following this, there was a performance comparison between the trained BPNN and ACO-BPNN using iteration times, network errors (E), and correlation coefficients (R).As shown in Figure4(a), an optimal solution was obtained after only five iterations through the ACO-BPNN model. However, the BPNN iterated 15 times to reach the convergence state, and the minimum network error was larger than the one obtained by ACO-BPNN. This indicated that the process of identifying an optimum solution by BPNN fell into the local optimal solution when modelling solar radiation. In addition to the rate of convergence, the MAPE, MABE, and RMSE of the ACO-BPNN were all less than the BPNN, and a higher correlation coefficient of observed and simulated monthly solar radiation was achieved by the ACO-BPNN (0.9808) than by BPNN (0.9759) (Table 3). Therefore, it was concluded that the modified ACO-BPNN performed superiorly to the traditional BPNN.Table 3
Statistic indicators of the BPNN and ACO-BPNN models.
Model
MAPE/%
MABE/MJ⋅m−2
RMSE/MJ⋅m−2
R
BPNN
6.49
23.68
35.60
0.9759
ACO-BPNN
6.45
23.03
31.45
0.9808Figure 4
Performance comparison between BPNN and ACO-BPNN using network error, iteration times (a), and correlation coefficient ((b) and (c)).
(a)
(b)
(c)Although the optimized ACO-BPNN accurately characterized the complicated functional relationship between solar radiation and the independent variables, further testing was required before applying it to estimate solar radiation at other meteorological stations without radiation observations. To achieve this, the modified ACO-BPNN needed validation at other radiation stations where observed data did not contribute to its training process.As shown in Figure5, a total of three verification stations are reserved for testing, including SuoLun, TongLiao, and ChangChun; radiation observations are collected at these locations from 2000 to 2010. The contrast result reveals a strong fit between the observed and simulated values.Figure 5
Comparison between observed and simulated values at SuoLun, TongLiao, and ChangChun from 2000 to 2010.Table4 lists the statistical error metrics and the model fit of the stations, where the MAPE was between 4.54% and 7.15%, the MABE ranged from 14.7073 MJ·m−2 to 27.1159 MJ·m−2, the RMSE varied from 21.6252 MJ·m−2 to 38.5571 MJ·m−2, and all the correlation coefficients of the modelled and observed values were greater than 0.97. This showed that the monthly solar radiation simulated by the ACO-BPNN achieved high accuracy.Table 4
Statistical indicators of the ACO-BPNN model at the verification stations.
Stations
MAPE/%
MABE/MJ⋅m−2
RMSE/MJ⋅m−2
R
SuoLun
4.54
17.2776
22.7772
0.9919
TongLiao
7.15
27.1159
38.5571
0.9729
ChangChun
4.50
14.7073
21.6252
0.9911Based on comparative and validating analysis of the ACO-BPNN model, it was determined that the ACO-BPNN model had superior performance in estimating solar radiation compared to traditional models and could be applied to other 96 weather stations lacking observations.


### 3.2. Spatial Patterns of Solar Radiation during Different Periods

Based on monthly solar radiation of 114 stations including 96 simulated weather stations and 18 observed radiation stations for the period 1982–2010, we obtained spatial patterns of solar radiation over NEC in the 1980s, 1990s, and 2000s by utilizing inverse distance weighting (IDW) interpolation method (Figure6). This not only merged solar radiation information at each year into the corresponding decade but also satisfied the precondition of following intensity analysis.Figure 6
Spatial patterns of annual average solar radiation over Northeast China from 1982 to 2010.On this basis, we classified solar radiation for different periods according to the established division standard (Table2). As shown in Figure 7, higher levels of solar radiation resources appeared in the south and west of the study region, while lower ones occurred in the north and east. Overall, the solar radiation resources show a steady decreasing trend from the southwest to the northeast across the study region.Figure 7
Spatial distribution of annual average solar radiation levels during the major temporal periods in Northeast China.To acquire further detailed variations of solar radiation resources over NEC during the three periods, statistics on the areas of solar radiation resources at different levels during various time periods were calculated. The results in Table5 highlight that the poor level area accounted for the largest proportion, while the area ratio of the abundant level was the smallest, and the area of normal and comparatively abundant levels accounted for the difference between the two from 1980s to 2000s. Compared to the 1980s, in the other decades the proportion of the normal level area was greater than the comparative level. The ratio of poor level area was more than 40% for all three periods, indicating that although the solar radiation resources presented an upward trend in recent years, the amount of resources was still relatively general overall.Table 5
Area statistics of vegetation coverage levels in each period.
Periods
Area & ratio
I
II
III
IV
1982–1990
Area (km2)
18923
275287
222010
725953
Ratio (%)
1.52%
22.16%
17.87%
58.44%
1991–2000
Area (km2)
65289
256145
270977
649762
Ratio (%)
5.26%
20.62%
21.81%
52.31%
2001–2010
Area (km2)
155302
256224
327766
502881
Ratio (%)
12.5%
20.63%
26.39%
40.48%
Variable area (km2)
1980s–1990s
46366
−19142
48967
−76191
1990s–2000s
90013
79
56789
−146881


### 3.3. Intensity Analysis of Spatiotemporal Dynamics of Solar Radiation

We analyzed the spatial pattern of solar radiation at each level, and the corresponding statistics results from 1982–2010 are reported in this section. The intensity analysis was introduced to analyze the dynamic changes of the solar radiation from three perspectives: time interval, category, and conversion. To achieve this, the ArcGIS software platform was used to obtain the solar radiation transition matrix at two time stages (the 1980s–1990s and 1990s–2000s), which was a prerequisite for intensity analysis. The transition matrix is displayed in Table6.Table 6
Transition matrix of different levels of solar radiation over two study periods (pixel).
Terminal point
I
II
III
IV
Initiation
Gross loss
Initial point
I
8931
9992
0
0
18923
9992
60524
2202
2563
0
65289
4765
II
56201
205098
13988
0
275287
70189
94778
145262
16105
0
256145
110883
III
157
39718
165295
16840
222010
56715
0
105535
160976
4466
270977
110001
IV
0
1337
91694
632922
725953
93031
0
3225
148122
498415
649762
151347
Gross gain
56358
51047
105682
16840
94778
110962
166790
4466
Termination
65289
256145
270977
649762
155302
256224
327766
502881The bold font refers to the variation of different levels of solar radiation from the 1980s to the 1990s, and the italicized font refers to the different solar radiation variation levels from the 1990s to the 2000s. On this basis, there is a further comprehensive analysis at the three levels consisting of time interval, category, and conversion for the table in the following sections.


#### 3.3.1. Results of Time Interval Analysis

Figure8 displays the analysis results for time intervals of solar radiation resources from the 1980s to the 2000s. The bars to the left of the central axis represent the overall intensity of solar radiation variation at each time interval. The bars right of the axis refer to time intensity calculated by formula (9). The dashed line on the right side was the uniform line acquired using formula (8). If the right-hand bar surpassed the uniform line, it was concluded that the variation of solar radiation resources levels at that time interval was comparatively fast. Otherwise, it meant that the variation was relatively slow.Figure 8
Time intensity analysis for the 1980s–1990s and 1990s–2000s over Northeast China.As shown in Figure8, there were notably more interval variations in solar radiation during the 1990s–2000s than in the 1980s–1990s. In addition, the rate of change in solar radiation from the 1990s to the 2000s was faster than that in the previous decade, as it passed the uniform line value (2.44%) for the whole period. Overall, regardless of the interval variations or changing rates, the solar radiation of the NEC region presented a tendency for growth.The results of the time interval analysis provided general variations in the regularity of solar radiation resources during different periods, and the variation laws of solar radiation resources at different levels are outlined in the next section.


#### 3.3.2. Results of Category Analysis

Figure9 gives the results of the category analysis, which was used to identify whether the variation of solar radiation resources at different levels over two time intervals was active or gentle. Each pair of horizontal bars denotes a level of solar radiation resources. The bars on the right refer to the annual intensity gains and losses at each time period, which were obtained using formulas (10). The bars on the left refer to variation area (represented by pixel numbers) at each time period and were derived using the numerator of formulas (10). If the right-hand bar surpassed the uniform line (the dashed line) which was the output result of formula (9), it was confirmed that the variation of the solar radiation resources at that time stage was comparatively active. Otherwise, it meant that the variation was relatively stable.Figure 9
Category intensity analysis for the 1980s–1990s (a1, a2) and 1990s–2000s (b1, b2); (a1)/(b1) denote gains and (a2)/(b2) are losses; A, CA, N, and P represent the classifications abundant, comparatively abundant, normal, and poor.The increasing and decreasing solar radiation conditions at each level from the 1980s to the 1990s are displayed in Figures9(a1) and 9(a2). The results revealed that the increasing area of the abundant level was larger than its decreasing area, and the growing area of the normal level was slightly larger than its diminishing area. Conversely, the increasing area of the poor level was less than its decreasing area, and this was likewise for the comparatively abundant level. Moreover, in terms of the state of intensity variation, both the increases and decreases of the abundant level solar radiation were active. Unlike the abundant level, the comparatively abundant and normal level of solar radiation presented relatively active increases and decreases. Contrary to the former three levels, both the increases and decreases of the poor level solar radiation presented were minimal. Overall, the solar radiation resources in the 1990s were richer than in the 1980s.Figures9(b1) and 9(b2) show the result of the category analysis for the 1990s–2000s. The analysis revealed that the increasing area of abundant level was far larger than its decreasing area, and the growing area of comparatively abundant and normal level almost equalled their diminishing area. Unlike the former three levels, the increasing area of poor level was nearly zero, which was far less than its decreasing area. For the variation in intensity, the abundant level gains were fairly active, whereas its losses were relatively minimal; both the gains and losses of the comparatively abundant and normal levels were relatively active. In contrast, the poor level exhibited both gentle gains and losses. On the whole, the solar radiation resources of the 2000s were richer than the 1990s.The category analysis not only confirmed the theory that the general intensity of solar radiation resources was increasing at an accelerated rate but also answered the question regarding the amount of activity variation between levels. However, the results regarding the mutual conversion relationships between solar radiation resources at each level were not involved and instead are presented in Section3.3.3.


#### 3.3.3. Results of Conversion Analysis

By comparing the conversion intensities of solar radiation resources from one level to other during each time interval, we identified which type of conversion was stronger. To achieve this goal, we summarized the dominant conversions for each time period according to the results of the category analyses obtained from formulas (11)–(14). As shown in Table 7, the leading form of solar radiation resources at each level during the 1980s–1990s was generally consistent with that from the 1990s–2000s. The poor, normal, and comparatively abundant levels converted to normal, comparatively abundant, and abundant levels, respectively. The abundant level predominantly converted into comparatively abundant and normal levels. In general, the poor, normal, and normal level resources mainly converted into higher levels, but the abundant level mainly converted to a lower level.Table 7
Conversion of dominant solar radiation in Northeast China in the 1980s–1990s and the 1990s–2000s.
Convert from
Convert into
1980s–1990s
1990s–2000s
Poor
Normal
Normal
Normal
Comparatively abundant
Comparatively abundant
Comparatively abundant
Abundant
Abundant
Abundant
Comparatively abundant
Comparatively abundant & normalThe conversion analysis further verified the result of the interval and category analysis and uncovered the conversion relationships between the different levels.


## 3.1. Comparative and Validating Analysis of ACO-BPNN Model

A total of 5220 input/output pairs were employed to train the BPNN and ACO-BPNN to fully establish the relationship between the input vectorX=[at,ws,vp,rh,ca,sd,dtr]T and the output vector Y=[solar  radiation]T at 15 stations from 1982 to 2010. Following this, there was a performance comparison between the trained BPNN and ACO-BPNN using iteration times, network errors (E), and correlation coefficients (R).As shown in Figure4(a), an optimal solution was obtained after only five iterations through the ACO-BPNN model. However, the BPNN iterated 15 times to reach the convergence state, and the minimum network error was larger than the one obtained by ACO-BPNN. This indicated that the process of identifying an optimum solution by BPNN fell into the local optimal solution when modelling solar radiation. In addition to the rate of convergence, the MAPE, MABE, and RMSE of the ACO-BPNN were all less than the BPNN, and a higher correlation coefficient of observed and simulated monthly solar radiation was achieved by the ACO-BPNN (0.9808) than by BPNN (0.9759) (Table 3). Therefore, it was concluded that the modified ACO-BPNN performed superiorly to the traditional BPNN.Table 3
Statistic indicators of the BPNN and ACO-BPNN models.
Model
MAPE/%
MABE/MJ⋅m−2
RMSE/MJ⋅m−2
R
BPNN
6.49
23.68
35.60
0.9759
ACO-BPNN
6.45
23.03
31.45
0.9808Figure 4
Performance comparison between BPNN and ACO-BPNN using network error, iteration times (a), and correlation coefficient ((b) and (c)).
(a)
(b)
(c)Although the optimized ACO-BPNN accurately characterized the complicated functional relationship between solar radiation and the independent variables, further testing was required before applying it to estimate solar radiation at other meteorological stations without radiation observations. To achieve this, the modified ACO-BPNN needed validation at other radiation stations where observed data did not contribute to its training process.As shown in Figure5, a total of three verification stations are reserved for testing, including SuoLun, TongLiao, and ChangChun; radiation observations are collected at these locations from 2000 to 2010. The contrast result reveals a strong fit between the observed and simulated values.Figure 5
Comparison between observed and simulated values at SuoLun, TongLiao, and ChangChun from 2000 to 2010.Table4 lists the statistical error metrics and the model fit of the stations, where the MAPE was between 4.54% and 7.15%, the MABE ranged from 14.7073 MJ·m−2 to 27.1159 MJ·m−2, the RMSE varied from 21.6252 MJ·m−2 to 38.5571 MJ·m−2, and all the correlation coefficients of the modelled and observed values were greater than 0.97. This showed that the monthly solar radiation simulated by the ACO-BPNN achieved high accuracy.Table 4
Statistical indicators of the ACO-BPNN model at the verification stations.
Stations
MAPE/%
MABE/MJ⋅m−2
RMSE/MJ⋅m−2
R
SuoLun
4.54
17.2776
22.7772
0.9919
TongLiao
7.15
27.1159
38.5571
0.9729
ChangChun
4.50
14.7073
21.6252
0.9911Based on comparative and validating analysis of the ACO-BPNN model, it was determined that the ACO-BPNN model had superior performance in estimating solar radiation compared to traditional models and could be applied to other 96 weather stations lacking observations.


## 3.2. Spatial Patterns of Solar Radiation during Different Periods

Based on monthly solar radiation of 114 stations including 96 simulated weather stations and 18 observed radiation stations for the period 1982–2010, we obtained spatial patterns of solar radiation over NEC in the 1980s, 1990s, and 2000s by utilizing inverse distance weighting (IDW) interpolation method (Figure6). This not only merged solar radiation information at each year into the corresponding decade but also satisfied the precondition of following intensity analysis.Figure 6
Spatial patterns of annual average solar radiation over Northeast China from 1982 to 2010.On this basis, we classified solar radiation for different periods according to the established division standard (Table2). As shown in Figure 7, higher levels of solar radiation resources appeared in the south and west of the study region, while lower ones occurred in the north and east. Overall, the solar radiation resources show a steady decreasing trend from the southwest to the northeast across the study region.Figure 7
Spatial distribution of annual average solar radiation levels during the major temporal periods in Northeast China.To acquire further detailed variations of solar radiation resources over NEC during the three periods, statistics on the areas of solar radiation resources at different levels during various time periods were calculated. The results in Table5 highlight that the poor level area accounted for the largest proportion, while the area ratio of the abundant level was the smallest, and the area of normal and comparatively abundant levels accounted for the difference between the two from 1980s to 2000s. Compared to the 1980s, in the other decades the proportion of the normal level area was greater than the comparative level. The ratio of poor level area was more than 40% for all three periods, indicating that although the solar radiation resources presented an upward trend in recent years, the amount of resources was still relatively general overall.Table 5
Area statistics of vegetation coverage levels in each period.
Periods
Area & ratio
I
II
III
IV
1982–1990
Area (km2)
18923
275287
222010
725953
Ratio (%)
1.52%
22.16%
17.87%
58.44%
1991–2000
Area (km2)
65289
256145
270977
649762
Ratio (%)
5.26%
20.62%
21.81%
52.31%
2001–2010
Area (km2)
155302
256224
327766
502881
Ratio (%)
12.5%
20.63%
26.39%
40.48%
Variable area (km2)
1980s–1990s
46366
−19142
48967
−76191
1990s–2000s
90013
79
56789
−146881


## 3.3. Intensity Analysis of Spatiotemporal Dynamics of Solar Radiation

We analyzed the spatial pattern of solar radiation at each level, and the corresponding statistics results from 1982–2010 are reported in this section. The intensity analysis was introduced to analyze the dynamic changes of the solar radiation from three perspectives: time interval, category, and conversion. To achieve this, the ArcGIS software platform was used to obtain the solar radiation transition matrix at two time stages (the 1980s–1990s and 1990s–2000s), which was a prerequisite for intensity analysis. The transition matrix is displayed in Table6.Table 6
Transition matrix of different levels of solar radiation over two study periods (pixel).
Terminal point
I
II
III
IV
Initiation
Gross loss
Initial point
I
8931
9992
0
0
18923
9992
60524
2202
2563
0
65289
4765
II
56201
205098
13988
0
275287
70189
94778
145262
16105
0
256145
110883
III
157
39718
165295
16840
222010
56715
0
105535
160976
4466
270977
110001
IV
0
1337
91694
632922
725953
93031
0
3225
148122
498415
649762
151347
Gross gain
56358
51047
105682
16840
94778
110962
166790
4466
Termination
65289
256145
270977
649762
155302
256224
327766
502881The bold font refers to the variation of different levels of solar radiation from the 1980s to the 1990s, and the italicized font refers to the different solar radiation variation levels from the 1990s to the 2000s. On this basis, there is a further comprehensive analysis at the three levels consisting of time interval, category, and conversion for the table in the following sections.


### 3.3.1. Results of Time Interval Analysis

Figure8 displays the analysis results for time intervals of solar radiation resources from the 1980s to the 2000s. The bars to the left of the central axis represent the overall intensity of solar radiation variation at each time interval. The bars right of the axis refer to time intensity calculated by formula (9). The dashed line on the right side was the uniform line acquired using formula (8). If the right-hand bar surpassed the uniform line, it was concluded that the variation of solar radiation resources levels at that time interval was comparatively fast. Otherwise, it meant that the variation was relatively slow.Figure 8
Time intensity analysis for the 1980s–1990s and 1990s–2000s over Northeast China.As shown in Figure8, there were notably more interval variations in solar radiation during the 1990s–2000s than in the 1980s–1990s. In addition, the rate of change in solar radiation from the 1990s to the 2000s was faster than that in the previous decade, as it passed the uniform line value (2.44%) for the whole period. Overall, regardless of the interval variations or changing rates, the solar radiation of the NEC region presented a tendency for growth.The results of the time interval analysis provided general variations in the regularity of solar radiation resources during different periods, and the variation laws of solar radiation resources at different levels are outlined in the next section.


### 3.3.2. Results of Category Analysis

Figure9 gives the results of the category analysis, which was used to identify whether the variation of solar radiation resources at different levels over two time intervals was active or gentle. Each pair of horizontal bars denotes a level of solar radiation resources. The bars on the right refer to the annual intensity gains and losses at each time period, which were obtained using formulas (10). The bars on the left refer to variation area (represented by pixel numbers) at each time period and were derived using the numerator of formulas (10). If the right-hand bar surpassed the uniform line (the dashed line) which was the output result of formula (9), it was confirmed that the variation of the solar radiation resources at that time stage was comparatively active. Otherwise, it meant that the variation was relatively stable.Figure 9
Category intensity analysis for the 1980s–1990s (a1, a2) and 1990s–2000s (b1, b2); (a1)/(b1) denote gains and (a2)/(b2) are losses; A, CA, N, and P represent the classifications abundant, comparatively abundant, normal, and poor.The increasing and decreasing solar radiation conditions at each level from the 1980s to the 1990s are displayed in Figures9(a1) and 9(a2). The results revealed that the increasing area of the abundant level was larger than its decreasing area, and the growing area of the normal level was slightly larger than its diminishing area. Conversely, the increasing area of the poor level was less than its decreasing area, and this was likewise for the comparatively abundant level. Moreover, in terms of the state of intensity variation, both the increases and decreases of the abundant level solar radiation were active. Unlike the abundant level, the comparatively abundant and normal level of solar radiation presented relatively active increases and decreases. Contrary to the former three levels, both the increases and decreases of the poor level solar radiation presented were minimal. Overall, the solar radiation resources in the 1990s were richer than in the 1980s.Figures9(b1) and 9(b2) show the result of the category analysis for the 1990s–2000s. The analysis revealed that the increasing area of abundant level was far larger than its decreasing area, and the growing area of comparatively abundant and normal level almost equalled their diminishing area. Unlike the former three levels, the increasing area of poor level was nearly zero, which was far less than its decreasing area. For the variation in intensity, the abundant level gains were fairly active, whereas its losses were relatively minimal; both the gains and losses of the comparatively abundant and normal levels were relatively active. In contrast, the poor level exhibited both gentle gains and losses. On the whole, the solar radiation resources of the 2000s were richer than the 1990s.The category analysis not only confirmed the theory that the general intensity of solar radiation resources was increasing at an accelerated rate but also answered the question regarding the amount of activity variation between levels. However, the results regarding the mutual conversion relationships between solar radiation resources at each level were not involved and instead are presented in Section3.3.3.


### 3.3.3. Results of Conversion Analysis

By comparing the conversion intensities of solar radiation resources from one level to other during each time interval, we identified which type of conversion was stronger. To achieve this goal, we summarized the dominant conversions for each time period according to the results of the category analyses obtained from formulas (11)–(14). As shown in Table 7, the leading form of solar radiation resources at each level during the 1980s–1990s was generally consistent with that from the 1990s–2000s. The poor, normal, and comparatively abundant levels converted to normal, comparatively abundant, and abundant levels, respectively. The abundant level predominantly converted into comparatively abundant and normal levels. In general, the poor, normal, and normal level resources mainly converted into higher levels, but the abundant level mainly converted to a lower level.Table 7
Conversion of dominant solar radiation in Northeast China in the 1980s–1990s and the 1990s–2000s.
Convert from
Convert into
1980s–1990s
1990s–2000s
Poor
Normal
Normal
Normal
Comparatively abundant
Comparatively abundant
Comparatively abundant
Abundant
Abundant
Abundant
Comparatively abundant
Comparatively abundant & normalThe conversion analysis further verified the result of the interval and category analysis and uncovered the conversion relationships between the different levels.


## 3.3.1. Results of Time Interval Analysis

Figure8 displays the analysis results for time intervals of solar radiation resources from the 1980s to the 2000s. The bars to the left of the central axis represent the overall intensity of solar radiation variation at each time interval. The bars right of the axis refer to time intensity calculated by formula (9). The dashed line on the right side was the uniform line acquired using formula (8). If the right-hand bar surpassed the uniform line, it was concluded that the variation of solar radiation resources levels at that time interval was comparatively fast. Otherwise, it meant that the variation was relatively slow.Figure 8
Time intensity analysis for the 1980s–1990s and 1990s–2000s over Northeast China.As shown in Figure8, there were notably more interval variations in solar radiation during the 1990s–2000s than in the 1980s–1990s. In addition, the rate of change in solar radiation from the 1990s to the 2000s was faster than that in the previous decade, as it passed the uniform line value (2.44%) for the whole period. Overall, regardless of the interval variations or changing rates, the solar radiation of the NEC region presented a tendency for growth.The results of the time interval analysis provided general variations in the regularity of solar radiation resources during different periods, and the variation laws of solar radiation resources at different levels are outlined in the next section.


## 3.3.2. Results of Category Analysis

Figure9 gives the results of the category analysis, which was used to identify whether the variation of solar radiation resources at different levels over two time intervals was active or gentle. Each pair of horizontal bars denotes a level of solar radiation resources. The bars on the right refer to the annual intensity gains and losses at each time period, which were obtained using formulas (10). The bars on the left refer to variation area (represented by pixel numbers) at each time period and were derived using the numerator of formulas (10). If the right-hand bar surpassed the uniform line (the dashed line) which was the output result of formula (9), it was confirmed that the variation of the solar radiation resources at that time stage was comparatively active. Otherwise, it meant that the variation was relatively stable.Figure 9
Category intensity analysis for the 1980s–1990s (a1, a2) and 1990s–2000s (b1, b2); (a1)/(b1) denote gains and (a2)/(b2) are losses; A, CA, N, and P represent the classifications abundant, comparatively abundant, normal, and poor.The increasing and decreasing solar radiation conditions at each level from the 1980s to the 1990s are displayed in Figures9(a1) and 9(a2). The results revealed that the increasing area of the abundant level was larger than its decreasing area, and the growing area of the normal level was slightly larger than its diminishing area. Conversely, the increasing area of the poor level was less than its decreasing area, and this was likewise for the comparatively abundant level. Moreover, in terms of the state of intensity variation, both the increases and decreases of the abundant level solar radiation were active. Unlike the abundant level, the comparatively abundant and normal level of solar radiation presented relatively active increases and decreases. Contrary to the former three levels, both the increases and decreases of the poor level solar radiation presented were minimal. Overall, the solar radiation resources in the 1990s were richer than in the 1980s.Figures9(b1) and 9(b2) show the result of the category analysis for the 1990s–2000s. The analysis revealed that the increasing area of abundant level was far larger than its decreasing area, and the growing area of comparatively abundant and normal level almost equalled their diminishing area. Unlike the former three levels, the increasing area of poor level was nearly zero, which was far less than its decreasing area. For the variation in intensity, the abundant level gains were fairly active, whereas its losses were relatively minimal; both the gains and losses of the comparatively abundant and normal levels were relatively active. In contrast, the poor level exhibited both gentle gains and losses. On the whole, the solar radiation resources of the 2000s were richer than the 1990s.The category analysis not only confirmed the theory that the general intensity of solar radiation resources was increasing at an accelerated rate but also answered the question regarding the amount of activity variation between levels. However, the results regarding the mutual conversion relationships between solar radiation resources at each level were not involved and instead are presented in Section3.3.3.


## 3.3.3. Results of Conversion Analysis

By comparing the conversion intensities of solar radiation resources from one level to other during each time interval, we identified which type of conversion was stronger. To achieve this goal, we summarized the dominant conversions for each time period according to the results of the category analyses obtained from formulas (11)–(14). As shown in Table 7, the leading form of solar radiation resources at each level during the 1980s–1990s was generally consistent with that from the 1990s–2000s. The poor, normal, and comparatively abundant levels converted to normal, comparatively abundant, and abundant levels, respectively. The abundant level predominantly converted into comparatively abundant and normal levels. In general, the poor, normal, and normal level resources mainly converted into higher levels, but the abundant level mainly converted to a lower level.Table 7
Conversion of dominant solar radiation in Northeast China in the 1980s–1990s and the 1990s–2000s.
Convert from
Convert into
1980s–1990s
1990s–2000s
Poor
Normal
Normal
Normal
Comparatively abundant
Comparatively abundant
Comparatively abundant
Abundant
Abundant
Abundant
Comparatively abundant
Comparatively abundant & normalThe conversion analysis further verified the result of the interval and category analysis and uncovered the conversion relationships between the different levels.


## 4. Discussion

### 4.1. Estimation of Solar Radiation

The reasons for variations in solar radiation are complicated, as the latter has complex coupling relationships with many meteorological factors. Therefore, the BPNN model can be applied to estimate solar radiation, as it is suitable for addressing fuzzy problems that require many factors and uncertainties to be considered. However, the traditional BPNN model also has some disadvantages; for example, if the initial values of set parameters are not appropriate, the optimization process easily generates default local minimums. To solve this problem and simulate more accurate monthly solar radiation values, we utilized an ACO algorithm to modify the traditional BPNN and constructed a new ACO-BPNN model. Our idea was to make full use of the global optimization ability of the ACO algorithm to provide several optimized groups of initial weighted values for the BPNN, to further refine the most appropriate weights for the BPNN, and finally to obtain an optimal global solution.Comparative analysis revealed that when considering convergence speed and the network error, the improved ACO-BPNN had a superior performance to the BPNN model. Through further contrast and verification of the observed and simulated values, it was revealed that the modified ACO-BPNN had excellent generalization in the study area and could be used to estimate the solar radiation of other nonradiation data stations.The ACO-BPNN developed in this study is a promising method for radiation prediction, especially for the remote regions without radiation observation, but there are possibilities for improvement in future studies. Namely, as the site data used in this study were discrete, solar radiation with continuous spatial distribution from remotely sensed data would be a desirable addition for the ACO-BPNN model.


### 4.2. Dynamic Intensity of Solar Radiation

In this study, the intensity analysis method was introduced into the dynamic study of solar radiation. Our research results revealed that, on the whole, the solar radiation over NEC had an obviously upward trend during the study period of 1980–2010, which was consistent with results of previous studies [44–46]. This estimated increasing trend was also in accordance with the findings of Wild et al. (2005), who deemed that the surface of the Earth has shifted from “dark” to “bright” over the past decade. Specifically, the solar radiation reaching the Earth’s surface presented an overall increasing trend from the 1980s onwards, which indicated that the variation trend of solar radiation might have global characteristics [53].The intensity analysis method was applied to deeply excavate the spatiotemporal variation of solar radiation across NEC from category and conversion aspects. Compared with previous trend, cycle, and principal component analysis methods based on station data, this study focused on the area of solar radiation at each level, which made a more detailed spatiotemporal analysis possible by considering both the variation and conversion of different levels of solar radiation. Nevertheless, there were still some shortcomings in the intensity analysis method that need to be further developed in future studies; one such example is the lack of locational information in the analysis process.


## 4.1. Estimation of Solar Radiation

The reasons for variations in solar radiation are complicated, as the latter has complex coupling relationships with many meteorological factors. Therefore, the BPNN model can be applied to estimate solar radiation, as it is suitable for addressing fuzzy problems that require many factors and uncertainties to be considered. However, the traditional BPNN model also has some disadvantages; for example, if the initial values of set parameters are not appropriate, the optimization process easily generates default local minimums. To solve this problem and simulate more accurate monthly solar radiation values, we utilized an ACO algorithm to modify the traditional BPNN and constructed a new ACO-BPNN model. Our idea was to make full use of the global optimization ability of the ACO algorithm to provide several optimized groups of initial weighted values for the BPNN, to further refine the most appropriate weights for the BPNN, and finally to obtain an optimal global solution.Comparative analysis revealed that when considering convergence speed and the network error, the improved ACO-BPNN had a superior performance to the BPNN model. Through further contrast and verification of the observed and simulated values, it was revealed that the modified ACO-BPNN had excellent generalization in the study area and could be used to estimate the solar radiation of other nonradiation data stations.The ACO-BPNN developed in this study is a promising method for radiation prediction, especially for the remote regions without radiation observation, but there are possibilities for improvement in future studies. Namely, as the site data used in this study were discrete, solar radiation with continuous spatial distribution from remotely sensed data would be a desirable addition for the ACO-BPNN model.


## 4.2. Dynamic Intensity of Solar Radiation

In this study, the intensity analysis method was introduced into the dynamic study of solar radiation. Our research results revealed that, on the whole, the solar radiation over NEC had an obviously upward trend during the study period of 1980–2010, which was consistent with results of previous studies [44–46]. This estimated increasing trend was also in accordance with the findings of Wild et al. (2005), who deemed that the surface of the Earth has shifted from “dark” to “bright” over the past decade. Specifically, the solar radiation reaching the Earth’s surface presented an overall increasing trend from the 1980s onwards, which indicated that the variation trend of solar radiation might have global characteristics [53].The intensity analysis method was applied to deeply excavate the spatiotemporal variation of solar radiation across NEC from category and conversion aspects. Compared with previous trend, cycle, and principal component analysis methods based on station data, this study focused on the area of solar radiation at each level, which made a more detailed spatiotemporal analysis possible by considering both the variation and conversion of different levels of solar radiation. Nevertheless, there were still some shortcomings in the intensity analysis method that need to be further developed in future studies; one such example is the lack of locational information in the analysis process.


## 5. Conclusions

In this study, we first utilized the ACO-BPNN model modified by the ACO algorithm to simulate monthly solar radiation and then used IDW interpolation to acquire monthly spatial patterns based on the radiation data. This included simulated radiation of 96 stations and observed radiation of 18 stations from 1982 to 2010. Based on the results, the spatiotemporal variation of solar radiation resources over NEC was analyzed from the aspects of time interval, category, and conversion by intensity analysis. The main conclusions can be summarized as follows:(1) The ACO-BPNN model outperformed the traditional BPNN, whether based on the network error and convergence rate or fit between simulation and observation, in all aspects. The contrast verification between simulated and observed solar radiation at the SuoLun, TongLiao, and ChangChun stations showed that the time series of simulated solar radiation were highly consistent with observed ones, which indicated that the ACO-BPNN we developed can be applied to model solar radiation in regions lacking radiation data.(2) The solar radiation of the study region had a downward trend from the southwest to the northeast during the 1980s–2000s. In the 1980s, the areas of varying levels of solar radiation from large to small were as follows: the poor level, the comparatively abundant level, the normal level, and the abundant level. However, during the 1990s and 2000s, the area of solar radiation at each level was in the following order: poor > normal > comparatively abundant > abundant. Overall, the solar radiation resources in NEC continued to increase during the study period of 1980s–2000s.(3) The annual average rate of variation of solar radiation during the 1980s–1990s was lower than that during the 1990s–2000s. From the 1980s–1990s, solar radiation was actively increasing and decreasing in abundance, and the abundant and normal levels exhibited relatively active increases and decreases whereas the variation level experienced comparatively gentle fluctuations. During the 1990s–2000s, the abundant level had active gains, but the losses were subtle. The increases and decreases of the comparatively abundant and normal levels were active. In contrast, the poor level showed a gentle increase and decrease. The poor, normal, and comparatively abundant levels were converted to higher levels, whereas the abundant level transformed to lower levels.



---

*Source: 1042603-2017-08-24.xml*