Spatiotemporal Optimization for Short-Term Solar Forecasting Based on Satellite Imagery

: Solar forecasting is essential for optimizing the integration of solar photovoltaic energy into a power grid. This study presents solar forecasting models based on satellite imagery. The cloud motion vector (CMV) model is the most popular satellite-image-based solar forecasting model. However, it assumes constant cloud states, and its accuracy is, thus, inﬂuenced by changes in local weather characteristics. To overcome this limitation, satellite images are used to provide spatial data for a new spatiotemporal optimized model for solar forecasting. Four satellite-image-based solar forecasting models (a persistence model, CMV, and two proposed models that use clear-sky index change) are evaluated. The error distributions of the models and their spatial characteristics over the test area are analyzed. All models exhibited different performances according to the forecast horizon and location. Spatiotemporal optimization of the best model is then conducted using best-model maps, and our results show that the skill score of the optimized model is 21% better than the previous CMV model. It is, thus, considered to be appropriate for use in short-term forecasting over large areas. The results of this study are expected to promote the use of spatial data in solar forecasting models, which could improve their accuracy and provide various insights for the planning and operation of photovoltaic plants.


Introduction
The growing importance of solar forecasting has been emphasized in many published studies and in review articles comparing different forecasting techniques. However, review papers [1][2][3] have categorized previous studies according to the method, forecast horizon, and other characteristics. In addition, text mining of literature [4] and those focused on machine learning methods [5]. However, it can be identified that these reviews did not propose one optimal method [6]. It is difficult to evaluate the superiority of models as their performances depend on evaluating their metrics, forecast horizons, and the region studied. Although several objective comparison methods have been proposed [6,7], their performances depend on spatiotemporal features, which need to be thoroughly considered.
Although ensemble or machine learning models have been shown to perform well, models suited to specific forecasting conditions are required. As it is difficult to develop one model that provides optimal results for every geographical region and forecast horizon. Thus, hybrid or ensemble models containing several different models have been introduced, and these have been shown to exhibit better performances than individual models [8]. Many high-performance specific models have been selected and combined in this respect [9]. Machine learning techniques have yielded remarkable performances when using such models as input. However, to improve the results of these models, it is essential to gain physical interpretation and domain knowledge of the specific models. able performance was observed for a CMV model using spatial averaging as a smoothing method [11,14], but the smoothing was performed using a simple filter application. Spatiotemporal correlation [27] analyses have been conducted, but analyses were only conducted to enhance CMV extraction in a small-scale area. Some simple comparisons between different regions have been made, such as between different climate zones within the same area [8]; however, most studies have used site-specific or overall-averaging results. As the available solar energy incident on Earth depends on various weather conditions, it must be correlated with the terrain (land or sea, and coastal or mountain regions). Therefore, conducting a spatial analysis using distribution patterns and topography could contribute to error analyses and spatial optimization. In this respect, satellite images provide uninterrupted information about an entire region, compared to specific site-measurement data, and are, hence, suitable for spatial analysis studies.
The main aim of this study was to analyze spatial characteristics relevant to the development of solar-forecasting models with the objective of improving forecasting performance. The following steps were applied: (1) The characteristics of satellite-imagebased forecasting models containing CMV were compared; (2) a spatial analysis of the performance of these models was conducted; and (3) spatiotemporal optimization was conducted to maximize the forecasting performance. A new spatiotemporal forecasting model was then developed as a hybrid of CMV and additional models considering their spatial accuracy. The structure of this model is shown in Figure 1, where it is also compared to the structure of most previously studied models. The new models in red, such as delta persistence (DEL), are described in detail in the methodology chapter.
To date, there have been insufficient evaluations of the spatial accuracy and distribution and their correlations with other factors in forecasting studies. Although CMV models use spatial information and improvements have been made, only ambient spatial information is included, and other factors or distributions are not considered. A remarkable performance was observed for a CMV model using spatial averaging as a smoothing method [11,14], but the smoothing was performed using a simple filter application. Spatiotemporal correlation [27] analyses have been conducted, but analyses were only conducted to enhance CMV extraction in a small-scale area. Some simple comparisons between different regions have been made, such as between different climate zones within the same area [8]; however, most studies have used site-specific or overall-averaging results. As the available solar energy incident on Earth depends on various weather conditions, it must be correlated with the terrain (land or sea, and coastal or mountain regions). Therefore, conducting a spatial analysis using distribution patterns and topography could contribute to error analyses and spatial optimization. In this respect, satellite images provide uninterrupted information about an entire region, compared to specific site-measurement data, and are, hence, suitable for spatial analysis studies.
The main aim of this study was to analyze spatial characteristics relevant to the development of solar-forecasting models with the objective of improving forecasting performance. The following steps were applied: (1) The characteristics of satellite-image-based forecasting models containing CMV were compared; (2) a spatial analysis of the performance of these models was conducted; and (3) spatiotemporal optimization was conducted to maximize the forecasting performance. A new spatiotemporal forecasting model was then developed as a hybrid of CMV and additional models considering their spatial accuracy. The structure of this model is shown in Figure 1, where it is also compared to the structure of most previously studied models. The new models in red, such as delta persistence (DEL), are described in detail in the methodology chapter.

Study Area
The study area was the South Korean peninsula, which has complex topography. The location of the study area and associated basic spatial information are shown in Figure 2, and the area corresponds to a region spanning 636 × 459 km. The digital elevation model (DEM), also shown in Figure 2, was resampled from Shuttle Radar Topography Mission (SRTM) data. The study area experiences complex meteorological phenomena as it is surrounded by sea and mountainous regions (a long mountain range in the east and Jeju Island, which is a large island in the south, are visible on the map in Figure 2). The area lies within a temperate climate zone, experiences four

Study Area
The study area was the South Korean peninsula, which has complex topography. The location of the study area and associated basic spatial information are shown in Figure 2, and the area corresponds to a region spanning 636 × 459 km. The digital elevation model (DEM), also shown in Figure 2, was resampled from Shuttle Radar Topography Mission (SRTM) data. The study area experiences complex meteorological phenomena as it is surrounded by sea and mountainous regions (a long mountain range in the east and Jeju Island, which is a large island in the south, are visible on the map in Figure 2). The area lies within a temperate climate zone, experiences four seasons, and is affected by the East Asian monsoon. The mean global horizontal irradiance (GHI) of South Korea is approximately 3.5 kWh/m 2 /day. The spatial distribution of GHI and the clear-sky index (kt) are also shown; these were estimated from satellite images, and the process used is described in detail later. Only data from 2018 were used, and it is acknowledged that annual volatility may result in different distributions. Nevertheless, differences relating to the land and sea areas are obvious, regardless of such inter-year volatility. Although all variables (mean GHI, mean kt, and standard deviation (STD) of kt) are related to the solar potential, they show slightly different distributions.

Satellite Images
The Communication, Ocean and Meteorological Satellite (COMS), a geostationary satellite operated by the Korea Meteorological Administration, was used to obtain the satellite images for this study. The satellite images had a spatial resolution of 1 km and a temporal resolution of 15 min. Only processed daytime images (08:30 to 17:00) for the entire period of 2018 were used, and no raw images were employed. The kt and clear-sky GHI images were extracted from the COMS images using the University of Arizona Solar Irradiance Based on Satellite-Korea Institute of Energy Research (UASIBS/KIER) model [35], and subsequently used for solar forecasting, and the detailed algorithm and accuracy results can be found in the associated reference. Figure 1 shows that the kt images were used in forecasting, and clear-sky GHI images were used to calculate GHI from kt.

Satellite-Image-Based Solar Forecasting Models
Four models were used for satellite-image-based solar forecasting. The persistence (PER) and CMV models have been reported in previous studies and were used as the reference and standard models, respectively. The delta persistence (DEL) model and a hybrid model of the CMV and DEL (CMV+DEL) were initially proposed in this study to overcome the limitations of the existing CMV models. Simple diagrams of all the models are shown in Figure 3. The models were analyzed using 100 × 100 km sample images over a period of 5 h from processed satellite images, and the forecasting models were implemented using Python code and open-source libraries, such as OpenPIV, OpenCV, and TensorFlow. The spatial distribution of GHI and the clear-sky index (k t ) are also shown; these were estimated from satellite images, and the process used is described in detail later. Only data from 2018 were used, and it is acknowledged that annual volatility may result in different distributions. Nevertheless, differences relating to the land and sea areas are obvious, regardless of such inter-year volatility. Although all variables (mean GHI, mean k t , and standard deviation (STD) of k t ) are related to the solar potential, they show slightly different distributions.

Satellite Images
The Communication, Ocean and Meteorological Satellite (COMS), a geostationary satellite operated by the Korea Meteorological Administration, was used to obtain the satellite images for this study. The satellite images had a spatial resolution of 1 km and a temporal resolution of 15 min. Only processed daytime images (08:30 to 17:00) for the entire period of 2018 were used, and no raw images were employed. The k t and clearsky GHI images were extracted from the COMS images using the University of Arizona Solar Irradiance Based on Satellite-Korea Institute of Energy Research (UASIBS/KIER) model [35], and subsequently used for solar forecasting, and the detailed algorithm and accuracy results can be found in the associated reference. Figure 1 shows that the k t images were used in forecasting, and clear-sky GHI images were used to calculate GHI from k t .

Satellite-Image-Based Solar Forecasting Models
Four models were used for satellite-image-based solar forecasting. The persistence (PER) and CMV models have been reported in previous studies and were used as the reference and standard models, respectively. The delta persistence (DEL) model and a hybrid model of the CMV and DEL (CMV+DEL) were initially proposed in this study to overcome the limitations of the existing CMV models. Simple diagrams of all the models are shown in Figure 3. The models were analyzed using 100 × 100 km sample images over a period of 5 h from processed satellite images, and the forecasting models were implemented using Python code and open-source libraries, such as OpenPIV, OpenCV, and TensorFlow.

Persistence Model
The PER model is commonly used as a reference model [3]. It has several variations, but we used the "clear-sky persistence" model [6], which uses the last-measured kt for all periods. The PER model is described by Equation (1), where x,y denotes the location of kt, and tn denotes the nth time point of kt. The calculated kt images at t0 are applied to images at t1 to tn with different clear-sky GHI values. This model assumes that cloud covers remains unchanged, and its expected accuracy is, therefore, low. Compared to previous site-specific studies, the spatial analysis of the PER model can show the reference accuracy in each location. In addition, the distribution pattern and distinct locations can be identified.

Persistence Model
The PER model is commonly used as a reference model [3]. It has several variations, but we used the "clear-sky persistence" model [6], which uses the last-measured k t for all periods. The PER model is described by Equation (1), where x,y denotes the location of k t , and t n denotes the nth time point of k t . The calculated k t images at t 0 are applied to images at t 1 to t n with different clear-sky GHI values. This model assumes that cloud covers remains unchanged, and its expected accuracy is, therefore, low. Compared to previous site-specific studies, the spatial analysis of the PER model can show the reference accuracy in each location. In addition, the distribution pattern and distinct locations can be identified.

Cloud Motion Vector Models
CMVs are based on the extraction of movement information from sequential images, and they have been used in many applications. Block matching is the most popular CMV model used for solar forecasting [11,14,15], and it identifies the vectors by estimating the displacement of similar areas in the before and after images. It is commonly used in particle image velocimetry (PIV) [36] for determining fluid dynamics. PIV uses a crosscorrelation for evaluation, whereas the CMV model applied in Germany [14] uses the mean squared error (MSE). The MSE is simple and easy to implement, but it requires extremely high computational power; therefore, a low spatial resolution of 43 km was used in the previous study. In contrast, the cross-correlation method (CCM), which calculates cross-correlation using Fourier transforms, is more appropriate for high-resolution images. However, it has limitations in similar pixel groups [37] and shows a lower performance than the OF method [20].
The OF method extracts movement using an image processing technique and is being actively researched for computer vision. Although block matching can also be categorized as OF, we use these two terms separately here. OF, thus, refers to differential methods, such as the Horn-Schunck method [38] or deep-learning-based models, which show smoother CMVs. Many associated models have been developed recently. Some studies have used OF methods, such as the Lucas-Kanade (LK) method [19,20] and Classic++ [21] for solar forecasting. However, as these models were optimized using different images and have many hyperparameters, they need to be objectively compared. Therefore, we compared several CMV extraction methods using processed satellite images: CCM, Farneback [39], Dual TVL1 [40], Deep Flow [41], Brox [42], and PWC-Net [43]. All of these methods are OF methods, except for CCM and PWC-Net, which are deep learning models. The optimized parameters for the Dual TVL1 were taken from similar previous studies [44].
The extracted CMV is applied to k t images using the warp method for image processing. When the time of the first image is t 0 , the CMV is estimated using images at t 0 and t 1 . First, forecasting is conducted using the image at t 1 , followed by forecasting using the forecasted image. Data from multiple locations can be considered for actual CMV processing. However, at a particular area at the boundary of the image, no data are recorded, due to no coming CMV data, and the PER model has been applied to such areas. Equation (2) describes the CMV model, where u and v denote the two vectors of CMV.

Delta Persistence (DEL) Model
A new type of persistence model, named delta persistence, was developed, which assumes that the change in k t is constant (see Equation (3)). The concept is similar to ramp persistence [28], but ramp persistence uses a summation term, whereas the DEL model uses a product term. The product term for the change in k t is referred to as the delta (∆) term in this study. As the rate of change is constant, the magnitude of the change continuously decreases over time. This model can consider the change in cloud thickness, but it cannot consider cloud movements. As k t is used as a denominator, the rate of change can be extremely large when k t is close to 0. Thus, the rate of change is set to be less than 1 to form an inverse index with increasing k t ,

CMV + DEL Model
The hybrid model of CMV and DEL model uses ∆ terms estimated from images from the CMV model. As shown in Equation (4), the difference between this model and the DEL is the association of k t_CMV with some terms, This model uses the CMV to accurately forecast future k t images, and the ∆ term is then used to consider the formation and disappearance of clouds. The ∆ terms are estimated using actual images at t 1 and forecasted images using the CMV. Theoretically, the model should provide a higher performance than the CMV and DEL separately, but its overall accuracy was not significantly improved because of errors in the extracted CMV and ∆ terms. A detailed error analysis was, thus, conducted to verify the accumulation of errors.

Spatial Averaging
Spatial averaging was applied as a smoothing method. This is not strictly a part of the forecasting model but is rather a post-processing technique for increasing the performance. Previous studies have shown the importance of smoothing in PER and CMV models [14]. During smoothing, the average value of nearby locations was used to reduce the spatial uncertainty. Rectangular filters with different sizes according to the forecast horizon [14] were used in this study. The size of the filter was selected based on the resulting accuracy, and was specifically 10, 20, 30, and 50 km for 0.25 h, 0.5 h, 0.75 h, and >1 h, respectively. Equation (5) is the simple equation used for smoothing, where FOR denotes all forecasting models, s denotes smoothing, and N denotes the number of nearby pixels during smoothing.
Smoothing was also applied to the ∆ term in the DEL and hybrid models, as this term also has spatial uncertainty. The equation for smoothing is similar to Equation (5), but it uses the ∆ term instead of the k t_FOR term. The size of the filter was 150 km for both models.

Spatial Accuracy Analysis
Three accuracy metrics were used in the analysis: The root mean squared error (RMSE), mean bias error (MBE), and the skill score (SS). The relative value, which divides the metric into a mean value was used for normalization. The SS was used to effectively visualize and compare the performances of the models. The simplified equation in Equation (6) [3] was used to calculate SS, where RMSE is the RMSE of the model being analyzed, and RMSE p is the RMSE of the persistence model. Although a combination of climatology and persistence models is recommended for SS calculations [6], we used the persistence model as a reference, as it is considered to be more accurate than a climatology model for short-term forecasting. The spatial accuracy of the models was visualized and analyzed by temporally averaging the metric in each location. However, the influence of weather conditions was greater than the difference between each model, therefore, the SS effectively showed the difference between the accuracy of the models and the reference model. A similar large distribution in the forecast accuracy was identified in the actual results. Accordingly, the RMSE and MBE maps of the PER model were analyzed, and then the SS and MBE maps of the other models were analyzed. To quantify the similarity of the distributions or interconnections, the correlations between several variables and metrics were considered in this analysis.

Spatiotemporal Optimization
Previous studies have considered temporal patterns according to the models and proposed ensemble models. The terms "intra-hour", "intra-day", and "a day ahead" have been introduced to represent these temporal patterns. However, as the majority of satellite-image-based forecasting models used in previous studies were CMV, temporal optimization has not been previously applied. It was applied here to compare the four models in the current study, and the additional consideration of spatial characteristics was expected to improve the model performance. Both temporal and spatial accuracy patterns were analyzed and used in the optimization process, and the temporal patterns were also analyzed by averaging the spatial values.
A new spatiotemporal optimization method was developed by mapping the best models according to the location. The best-model map was calculated by comparing the RMSE map for each model. The temporal best-model maps were estimated in each forecast horizon and used for the actual optimized (OPT) model. When calculating the best-model maps, the training and test data were separately used for validation, and twofold cross-validation was conducted to estimate the results in all periods. As the data were a time-series, they were separated daily. A physical interpretation of the best-model map was made considering the models' characteristics and other spatial information. Only the land area of South Korea was used for model comparison considering the practical usage.

Performance of Forecasting Models
The six CMV extraction models in the methodology chapter were compared. The differences between the models for extracting CMV values were found to be insignificant, probably because they all assumed that the shape or brightness of an object was constant. The maximum RMSE difference between the OF models was up to 3% in k t image forecasting, and the RMSE was 5% lower in the CCM model compared to the best OF models. This was evidenced by the fact that the difference between the CMV persistence forecasting and CMV updating for each image was marginal, and only a 2% difference in the RMSE was observed for k t image forecasting. Therefore, one of the most popular OF models, Dual TVL1, was selected as the CMV in this study, as it performed slightly better than the other models. Figure 4 shows an example of the results of the forecasting models using sample images. This example shows the evolution of complex weather phenomena in the study area with 15-min intervals. The low-k t regions, shown in blue, represent the clouds. The extracted CMV and ∆ terms are also identified in this figure. In general, the extracted CMV or DEL was applied to the image at t 1 ; however, in this figure, it was applied to the image at t 0 to analyze the performance of the model. The additional processing for the area with no data in the CMV model was not complemented in the images. Smoothing was applied to the ∆ terms for the DEL model, but it was not applied to the forecasted images in this figure. and DEL models were observed for the hybrid model. Although the hybrid model combines the advantages of both models, it also contains accumulated errors from both models, and this could make it difficult to accurately predict clouds in some cases. In the dataset used in this study, the clouds rarely moved or changed while maintaining their shape.  The forecast images clearly show the limitations and characteristics of the forecasting models for these complex weather conditions. The CMV model considers the movement of clouds, but the forecast shape is different from the ground truth (GT) because of the errors propagated by the CMV model and the rapid changes in cloud conditions. For the DEL model, the shape of the clouds differs from the GT, but the use of the ∆ term enables the overall trend of clearing clouds to be visualized. The characteristics of both the CMV and DEL models were observed for the hybrid model. Although the hybrid model combines the advantages of both models, it also contains accumulated errors from both models, and this could make it difficult to accurately predict clouds in some cases. In the dataset used in this study, the clouds rarely moved or changed while maintaining their shape.

Spatial Analysis
The spatial distributions of the RMSE and MBE of the PER model are shown in Figure 5. The metrics were calculated over 4 h using 16 images. There were no significant changes in the spatial distributions of the error as a function of the forecast horizon for the PER model.

Spatial Analysis
The spatial distributions of the RMSE and MBE of the PER model are shown in Figure 5. The metrics were calculated over 4 h using 16 images. There were no significant changes in the spatial distributions of the error as a function of the forecast horizon for the PER model.   The RMSE map is seen to follow the geographical map well, in which the land and sea regions were clearly different. The maximum difference in the RMSE was approximately 26%, which means that the RMSE of the area with the lowest accuracy is more than twice that of the area with the highest accuracy area. The color bar in the figure was normalized according to the quantile; thus, the actual minimum and maximum value may differ from the color bar. The mountain on Jeju Island shows the highest RMSE (>50%), but the area in the peninsula excluding Jeju has RMSE values of up to 40%. In contrast, the lowest RMSE of 25% is found in the sea near the peninsula. In the case of the MBE, values up to 5% were observed, and there was a tendency to overestimate the values in the inland highlands and underestimate the values in the remaining areas. The MBE map suggests that there was considerable cloud formation in the high-altitude areas, and less cloud cover in the low-altitude areas.
The SS and MBE values of the six forecasting models are shown in Figures 6 and 7, respectively; the differences between each model are clearly identifiable in these figures.
Smoothing increased the overall accuracy, as shown by the SS map of the smoothed model. However, this improvement was small, or even negative, for the eastern mountain range and Jeju Island, and this trend was confirmed for all the other models. The simple PER model was found to be better than the other models for forecasting complex weather at high altitude and at the coast. The performance of the CMVs model was improved by smoothing in many regions, but negligible or negative improvements were observed at the edge of the image. This was expected because there were no incoming clouds at the boundary, but the western region was more affected by this due to the influence of the westerly wind in the study area. The DELs and hybrid model generally underestimated the MBE values. This performance was improved by modifying the delta terms of the DELs model, although this effect was not significant for the RMSE values. The RMSE map is seen to follow the geographical map well, in which the land and sea regions were clearly different. The maximum difference in the RMSE was approximately 26%, which means that the RMSE of the area with the lowest accuracy is more than twice that of the area with the highest accuracy area. The color bar in the figure was normalized according to the quantile; thus, the actual minimum and maximum value may differ from the color bar. The mountain on Jeju Island shows the highest RMSE (>50%), but the area in the peninsula excluding Jeju has RMSE values of up to 40%. In contrast, the lowest RMSE of 25% is found in the sea near the peninsula. In the case of the MBE, values up to 5% were observed, and there was a tendency to overestimate the values in the inland highlands and underestimate the values in the remaining areas. The MBE map suggests that there was considerable cloud formation in the high-altitude areas, and less cloud cover in the low-altitude areas.
The SS and MBE values of the six forecasting models are shown in Figures 6 and 7, respectively; the differences between each model are clearly identifiable in these figures.
Smoothing increased the overall accuracy, as shown by the SS map of the smoothed model. However, this improvement was small, or even negative, for the eastern mountain range and Jeju Island, and this trend was confirmed for all the other models. The simple PER model was found to be better than the other models for forecasting complex weather at high altitude and at the coast. The performance of the CMVs model was improved by smoothing in many regions, but negligible or negative improvements were observed at the edge of the image. This was expected because there were no incoming clouds at the boundary, but the western region was more affected by this due to the influence of the westerly wind in the study area. The DELs and hybrid model generally underestimated the MBE values. This performance was improved by modifying the delta terms of the DELs model, although this effect was not significant for the RMSE values.  The correlations between several variables and metrics were calculated, as shown in Figure 8, where part Figure 8a shows data for the entire study region (land and seas) and Figure 8b shows data for the land region only. The correlations were made using spatial data, without considering the time-series variability. There was no significant correlation between the DEM and other variables, unlike the visualized pattern. The mean GHI showed a correlation of 60% with the DEM. The influence of the terrain was correlated with the sea; therefore, a high correlation was not seen for the DEM alone (a meaningful correlation is expected when complex phenomena that are correlated with the sea or altitude are considered). The RMSE of the PER model had a high correlation with mean kt and was close to 0.9 in all regions. Other similar variables showed low correlations, and these results were similar for the other models. Although the correlation was lower when only the land area was considered, the correlation was still over 0.7. This result can be used to estimate the reference forecasting accuracy. In terms of the SS of the forecasting models, PERs and DELs, the CMVs and hybrid models showed different distributions over the study areas because of the boundary effect of the CMV. This tendency was weaker inland; however, the difference between the CMVs and DELs models was still clear with respect to the low correlation value. Other models also showed similar but different SS distributions for the inland regions. The correlations between several variables and metrics were calculated, as shown in Figure 8, where part Figure 8a shows data for the entire study region (land and seas) and Figure 8b shows data for the land region only. The correlations were made using spatial data, without considering the time-series variability. There was no significant correlation between the DEM and other variables, unlike the visualized pattern. The mean GHI showed a correlation of 60% with the DEM. The influence of the terrain was correlated with the sea; therefore, a high correlation was not seen for the DEM alone (a meaningful correlation is expected when complex phenomena that are correlated with the sea or altitude are considered). The RMSE of the PER model had a high correlation with mean k t and was close to 0.9 in all regions. Other similar variables showed low correlations, and these results were similar for the other models. Although the correlation was lower when only the land area was considered, the correlation was still over 0.7. This result can be used to estimate the reference forecasting accuracy. In terms of the SS of the forecasting models, PERs and DELs, the CMVs and hybrid models showed different distributions over the study areas because of the boundary effect of the CMV. This tendency was weaker inland; however, the difference between the CMVs and DELs models was still clear with respect to the low correlation value. Other models also showed similar but different SS distributions for the inland regions. Energies 2021, 14, x FOR PEER REVIEW 13 of 21

Optimized Model
The effects of the forecast horizon on the performance results of the models are shown in Figures 9 and 10, and the mean accuracy of each forecasting model as a function of the forecast horizon is shown in Figure 9. As the mean value of only the inland area of South Korea is considered, the influence of the CMV boundary effect is insignificant. The RMSE continuously increases with the increasing forecast horizon because of the accumulation of image forecasting errors. However, the RMSE is expected to converge at approximately

Optimized Model
The effects of the forecast horizon on the performance results of the models are shown in Figures 9 and 10, and the mean accuracy of each forecasting model as a function of the forecast horizon is shown in Figure 9. As the mean value of only the inland area of South Korea is considered, the influence of the CMV boundary effect is insignificant. The RMSE continuously increases with the increasing forecast horizon because of the accumulation of image forecasting errors. However, the RMSE is expected to converge at approximately 50%, which is the RMSE of the model using the constant annual average k t . As the effect of the forecast horizon on the RMSE is larger than the difference between any of the models, the SS in Figure 10 is more effective in the analysis. The SS curves clearly show the differences in accuracy compared to the reference PER model. 50%, which is the RMSE of the model using the constant annual average kt. As the effect of the forecast horizon on the RMSE is larger than the difference between any of the models, the SS in Figure 10 is more effective in the analysis. The SS curves clearly show the differences in accuracy compared to the reference PER model.   The performances of the PERs and CMVs models improved after smoothing (~10% increase in SS). A comparison between the PERs, CMVs, and DELs models showed that the CMVs model provided a higher accuracy than other models for most of the forecast horizons. The performance of the PER model was lower than all the other models. The hybrid model showed the mean accuracy was worse than that of CMVs. In the case of the temporally optimized model, which uses the best model for each forecast horizon, the Energies 2021, 14, 2216 14 of 18 accuracy of that model will be the highest line in these models. This study suggests that additional factors, such as spatial optimization, can be used to further improve the model, as seen in the OPT model. Each forecasting model exhibited a different spatiotemporal behavior and performance for specific regions. Therefore, a best-model map produced using the RMSE values of the PER, CMVs, DELs, and hybrid models for a 4-h forecast horizon is shown in Figure 11. Although the PER model provided low accuracy, it was the best model for solar forecasting in some highland regions near the sea. The DELs model was selected as the best model in areas near those of the PER model and for most of the boundary area; the DELs model was also uniquely selected for the area near the southwest coast. The CMVs model was the best model for most of the central area, while the hybrid model was selected for some inland areas. In particular, the inland lowlands near the sea in the west were selected for the CMVs model. As the best forecasting model depends on the chosen forecast horizon, best-model maps for various forecast horizons (1, 2, 3, and 4 h) are also shown in Figure 12. With an increasing forecast horizon, the CMVs and hybrid models were initially the best models in many areas, but the DELs model later has expanded its area. In addition, there was an increase in the total area in which the hybrid model was the best model compared to the CMVs model when only the land area was considered. Therefore, the influence of the delta terms increased for longer forecast horizons. A larger boundary effect was also identified when the clouds moved. The PER and DELs models continuously outperformed in some areas in the eastern mountain range and Jeju Island. As the best forecasting model depends on the chosen forecast horizon, best-model maps for various forecast horizons (1, 2, 3, and 4 h) are also shown in Figure 12. With an increasing forecast horizon, the CMVs and hybrid models were initially the best models in many areas, but the DELs model later has expanded its area. In addition, there was an increase in the total area in which the hybrid model was the best model compared to the CMVs model when only the land area was considered. Therefore, the influence of the delta terms increased for longer forecast horizons. A larger boundary effect was also identified when the clouds moved. The PER and DELs models continuously outperformed in some areas in the eastern mountain range and Jeju Island. Spatiotemporal optimization was conducted using the temporal best-model maps. The best-model maps were estimated for each forecast horizon with 15 min intervals. The OPT model had an RMSE of 27% and SS of 18%, and the latter was 21% higher than that of the best previous model (CMVs). The RMSE of all models is shown in Figure 13. Only the accuracy of the inland region was considered. The CMVs, DELs, and hybrid models provided similar accuracies, but with respect to the temporal pattern, the OPT model outperformed the highest performance models for all forecast horizons, as shown in Figure  10, particularly for forecast horizons longer than 1 h. This result can be attributed to the consideration of the ∆ terms in areas with highly dynamic cloud activity. Spatiotemporal optimization was conducted using the temporal best-model maps. The best-model maps were estimated for each forecast horizon with 15 min intervals. The OPT model had an RMSE of 27% and SS of 18%, and the latter was 21% higher than that of the best previous model (CMVs). The RMSE of all models is shown in Figure 13. Only the accuracy of the inland region was considered. The CMVs, DELs, and hybrid models provided similar accuracies, but with respect to the temporal pattern, the OPT model outperformed the highest performance models for all forecast horizons, as shown in Figure 10, particularly for forecast horizons longer than 1 h. This result can be attributed to the consideration of the ∆ terms in areas with highly dynamic cloud activity.

Conclusions
This study compared satellite-image-based solar forecasting models and their spatial characteristics. The error distributions of the CMVs and DELs models were significantly different, and the simple PER model showed the highest accuracy for some highland areas near the sea. The best model was determined with respect to the forecast horizon and location. Overall, the proposed OPT model, which includes spatiotemporal optimization, outperformed all other models, particularly for forecast horizons of more than 1 h. These results indicate that it is necessary to consider spatial data when conducting solar forecasting. However, as it is difficult to consider spatial information in site-specific studies, the findings of this study are expected to be useful for complementing previous studies.
In addition to improving the forecasting performance, other significant information was derived from conducting spatial analyses. The accuracy of the forecasting models was related to the elevation of terrain and proximity to the sea, and it was highly correlated

Conclusions
This study compared satellite-image-based solar forecasting models and their spatial characteristics. The error distributions of the CMVs and DELs models were significantly different, and the simple PER model showed the highest accuracy for some highland areas near the sea. The best model was determined with respect to the forecast horizon and location. Overall, the proposed OPT model, which includes spatiotemporal optimization, outperformed all other models, particularly for forecast horizons of more than 1 h. These results indicate that it is necessary to consider spatial data when conducting solar forecast-ing. However, as it is difficult to consider spatial information in site-specific studies, the findings of this study are expected to be useful for complementing previous studies.
In addition to improving the forecasting performance, other significant information was derived from conducting spatial analyses. The accuracy of the forecasting models was related to the elevation of terrain and proximity to the sea, and it was highly correlated with the mean k t . This effect should, thus, be considered when calculating and mapping solar resources. Accuracy maps could be used to identify areas that need to be analyzed by models supplemented for complex weather, or areas that are advantageous for forecasting. These maps could also be used to select sites for future photovoltaic power plants, or for designing future forecasting models. In this respect, spatial data obtained from satellite images are highly valuable. The availability of such data and spatial analysis methods are expected to make continuous contributions to solar forecasting research and to the wider photovoltaic industry.