Global Spatial and Temporal Variation of the Combined Effect of Aerosol and Water Vapour on Solar Radiation

: This study aims to calculate the combined and individual effects of the optical thickness of aerosols (AOT) and precipitable water vapour (PWV) on the solar radiation reaching the Earth’s surface at a global scale and to analyse its spatial and temporal variation. For that purpose, a novel but validated methodology is applied to CERES SYN1deg products for the period 2000–2019. Spatial distributions of AOT and PWV effects, both individually and combined, show a close link with the spatial distributions of AOT and PWV. The spatially averaged combined effect results in a − 13.9% reduction in irradiance, while the average AOT effect is − 2.3%, and the PWV effect is − 12.1%. The temporal analysis focuses on detecting trends in the anomalies. The results show overall positive trends for AOT and PWV. Consequently, signiﬁcant negative overall trends are found for the effects. However, signiﬁcant positive trends for the individual AOT and the combined AOT-PWV effects are found in speciﬁc regions, such as the eastern United States, Europe or Asia, indicating successful emission control policies in these areas. This study contributes to a better understanding of the individual and combined effects of aerosols and water vapour on solar radiation at a global scale.


Introduction
Solar radiation reaching the Earth's surface plays a fundamental role in the surface energy balance [1]. Changes in the atmospheric composition can alter the downwelling surface solar radiation, and, therefore, the balance between incoming and outgoing radiation fluxes, resulting in changes in the temperature at the surface [2]. Aerosols and water vapour are two major atmospheric components that alter the radiation that reaches the Earth's surface on clear-sky days through scattering and absorption processes [3][4][5]. Their effects on the surface solar radiation are assorted, causing, for example, diurnal and seasonal temperature oscillations [6], changes in the surface energy balance [7] or the water cycle [8][9][10][11]. In addition, atmospheric water vapour influences the radiometric properties of the aerosol and, therefore, their effects on the downwelling surface solar radiation, impacting climate and climate change. Consequently, global and long-term assessment of the effects of aerosols and water vapour on solar radiation is essential for climatological studies [12].
In fact, since aerosols and water vapour are always present in the atmosphere and interact with each other, it is their combined effect that gives a complete description of their impact on downwelling radiation. Despite its great interest, this topic has been poorly studied so far (e.g., [29,30]). In this sense, Obregón et al. [30] took an interesting step by proposing and validating a new methodology to estimate this combined effect. Furthermore, Obregón et al. [29] applied the proposed methodology to the Mediterranean area showing it to be suitable for large areas. However, no study has been conducted at a global scale. Hence, in this study, we propose to apply this methodology to the entire globe.
This study aims to calculate and analyse the spatial and temporal variation at a global scale of the combined effect of aerosols and water vapour on the solar radiation reaching the Earth's surface. For this purpose, the new methodology proposed and validated by Obregón et al. [30] has been applied at a global scale for the first time. For this aim, aerosol and water vapour data provided by SYN1deg Ed4 CERES (Clouds and the Earth's Radiant Energy System) product are used. These data have already been used to monitor aerosols (e.g., [23]) and water vapour (e.g., [31]), but the combined effect of both components has not been addressed yet. Another important aspect in the present study is the period considered (2000-2019) that represents two decades of climate data derived from remote sensing. The long period and global scale allow performing a complete and significant analysis of the spatial and temporal variability of the combined effects of aerosol and water vapour on the downwelling SW irradiance, essential for evaluating climatic effects. The individual effects of aerosols and water vapour have also been calculated to compare the results with other studies at the global scale.
The paper is structured as follows: Section 2 describes the dataset and the methodology, the results obtained are presented and discussed in Section 3, and some conclusions are given in Section 4.

Dataset and Methodology
The combined effect of aerosols and water vapour on SW irradiance at ground level was calculated using the methodology recently developed and validated by Obregón et al. [30]. According to this methodology, aerosols are characterised by one of the most important aerosol optical quantities, the aerosol optical thickness (AOT), and water vapour by the precipitable water vapour (PWV).
In this study, global coverage AOT and PWV data provided by CERES (SYN1deg Ed4A) [32] have been used. The SYN1deg product combines Terra and Aqua CERES, MODIS and geostationary satellite observations, and Earth System model data. In particular, aerosol data came from the NASA/GSFC MODIS MOD04_L2/MYD04_L2 products and were assimilated by the Model for Atmospheric Transport and Chemistry (MATCH) [33]. The MATCH model assimilates and spatially and temporally interpolates MODIS aerosol optical thickness to the CERES resolutions. PWV data were obtained from the Global Modeling and Assimilation Office (GMAO)'s Goddard Earth Observing System Data Assimilation System (GEOS-DAS V5.4.1) product [34]. Products provided by CERES were adapted to the same temporal and spatial resolution, facilitating global studies. The acronym "SYN" (Synoptic Radiative Fluxes and Clouds) means that this version provides radiation data on clear and all-sky conditions, and "1deg" means it has a 1-degree spatial resolution. CERES samples 1-km MODIS data every four pixels and every two scan lines to result in an effective nominal resolution of 8 km 2 (2 km × 4 km) [35]. Cloud masking is performed at this spatial resolution, ensuring proper management of small clouds. At the end of the process, this high spatial resolution is then degraded to also provide products at lower spatial resolutions, such as the 1-degree products used in the present study. Despite being coarser than the initial resolution, it is suitable for global-scale studies, such as the Remote Sens. 2021, 13, 708 3 of 20 present one. In addition to the radiation flux, this version provides MODIS cloud properties and aerosols information. Regarding the cloud masking process, CERES Edition 4 (Ed4) used in this study represents a substantial improvement over the previous edition 2, mainly in terms of its MODIS cloud mask [36] The use of revised algorithms and input data to Terra and Aqua MODIS radiances greatly improve the accuracy of cloud amount. The Ed4 cloud mask employs 5-7 additional channels, new models for the clear sky, ocean and snow/ice surface radiances, and revised Terra MODIS calibrations. As a result, excellent consistency is obtained between the Aqua and Terra cloud fraction. Comparisons with CALIPSO (Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations) confirm that Ed4 cloud amounts are more or as accurate as other cloud mask systems available. Its mask product correctly identifies cloudy or clear areas 90%-96% of the time during the day and 88%-95% during the night over non-polar areas [35]. MODIS retrievals of AOT have been thoroughly validated by comparison to Aerosol Robotic Network (AERONET) data worldwide (e.g., [37][38][39][40][41]). PWV has also been validated by comparison to AERONET and GPS data [42][43][44]. Although no uncertainties for AOT and PWV as provided by the SYN1deg Ed4 CERES product have been reported, some uncertainties may be given for the original MODIS data from which they are derived. Thus, MODIS derived AOT uncertainty is estimated to be ±0.03±0.05*AOT over ocean and ±0.05±0.15*AOT over land [39]. The relative uncertainty of MODIS PWV values ranges between 5% and 10% [45]. The CERES consists of three-hourly, daily, and monthly data. In this study, daily AOT and PWV products from 1 March 2000 to 31 December 2019 were used to calculate the aerosol and water vapour combined effects on the downwelling SW irradiance. Only CERES measurements registered under clear sky conditions were selected for the calculation.
The methodology proposed by Obregón et al. [30] to estimate the combined effects of AOT and PWV is based on simulating the downwelling SW irradiance reaching the Earth's surface under cloud-free conditions and different values of AOT and PWV. These simulations were performed using the libRadtran (Version 1.7) radiative transfer model. Thus, a look-up-table with 100 simulation runs performed with AOT ranging from 0 to 1.5 and PWV from 0 to 60 mm, was obtained. The combined effect of AOT and PWV was calculated as the relative difference (Rel.Dif, %) of the SW irradiances with respect to an atmosphere with no aerosol and no water vapour (AOT = 0 and PWV = 0). The expression used to calculate the combined effect is [30]:  Similar to the case of the combined effect, two look-up tables were built for the individual AOT and PWV effects according to Equations (2) and (3), respectively. The dependence of the relative difference with AOT and PWV values is shown in Figures 2 and 3. As can be seen in Figure 2, AOT relative difference values were more negative as AOT increased (reaching values beyond −23.10%) and slightly depended on PWV values. Figure  3 shows that PWV relative difference values were also negative and were more negative as the PWV increased (reaching values beyond −18.48%). In this case, a large gradient in relative difference values for low PWV values (between 0 and 10 mm) was observed.
Obregón et al. [30] validated this methodology using irradiance measurements at nine stations scattered around the world and obtained differences lower than 3% in 84% of the cases. The sensitivity of the model to different input variables was also analysed by Obregón et al. [30], obtaining high robustness against changes in the aerosol type, vertical profiles of the thermodynamic variables, main gases and aerosols, and the surface albedo.
The combined and individual effects of aerosols and water vapour on the SW irradiance for specific values of AOT and PWV were obtained by linear interpolation of the values of the look-up-tables. This procedure was applied to a grid of 360 pixels latitude × 180 pixels longitude with daily AOT and PWV values provided by CERES. Thus, daily relative differences (Rel.Dif (AOT, PWV), Rel.Dif (AOT) and Rel.Dif (PWV)) for the globe were obtained, and their spatial and temporal variability was analysed. For temporal analysis, Sen's method [46] was applied to estimate the magnitude of the trend and the Mann-Kendall test [47,48] at 95% confidence level to estimate its statistical significance. Before trend detection, the time series were monthly-averaged and deseasonalised to eliminate the large influence of the annual cycle. Therefore, time series of AOT, PWV, AOT effects, PWV effects and AOT-PWV effects corresponded to monthly anomalies, i.e., differences between the monthly average for a particular year and the mean over all the period of study for that specific month.  Additionally, the individual effects of aerosols and water vapour were calculated to compare the results with other studies. The methodology applied to calculate the individual effects was the same as above, but, in this case, only one variable was allowed to vary. The variable fixed took the actual value of the variable in each location and time. The expression for AOT follows: where I is the simulated SW irradiance for an AOT value and I ref is the simulated SW irradiance for AOT equal to 0, while PWV value remains constant. Similarly, for PWV, the following expression is obtained: Similar to the case of the combined effect, two look-up tables were built for the individual AOT and PWV effects according to Equations (2) and (3), respectively. The dependence of the relative difference with AOT and PWV values is shown in Figures 2 and 3. As can be seen in Figure 2, AOT relative difference values were more negative as AOT increased (reaching values beyond −23.10%) and slightly depended on PWV values. Figure 3 shows that PWV relative difference values were also negative and were more negative as the PWV increased (reaching values beyond −18.48%). In this case, a large gradient in relative difference values for low PWV values (between 0 and 10 mm) was observed.

Spatial Analysis
Before analysing the combined effect of AOT-PWV on SW irradiance at ground level, the individual effects of AOT and PWV were analysed. Analysis of individual effects fa- Obregón et al. [30] validated this methodology using irradiance measurements at nine stations scattered around the world and obtained differences lower than 3% in 84% of the cases. The sensitivity of the model to different input variables was also analysed by Obregón et al. [30], obtaining high robustness against changes in the aerosol type, vertical profiles of the thermodynamic variables, main gases and aerosols, and the surface albedo.
The combined and individual effects of aerosols and water vapour on the SW irradiance for specific values of AOT and PWV were obtained by linear interpolation of the values of the look-up-tables. This procedure was applied to a grid of 360 pixels latitude × 180 pixels longitude with daily AOT and PWV values provided by CERES. Thus, daily relative differences (Rel.Dif (AOT, PWV), Rel.Dif (AOT) and Rel.Dif (PWV)) for the globe were obtained, and their spatial and temporal variability was analysed. For temporal analysis, Sen's method [46] was applied to estimate the magnitude of the trend and the Mann-Kendall test [47,48] at 95% confidence level to estimate its statistical significance. Before trend detection, the time series were monthly-averaged and deseasonalised to eliminate the large influence of the annual cycle. Therefore, time series of AOT, PWV, AOT effects, PWV effects and AOT-PWV effects corresponded to monthly anomalies, i.e., differences between the monthly average for a particular year and the mean over all the period of study for that specific month.

Spatial Analysis
Before analysing the combined effect of AOT-PWV on SW irradiance at ground level, the individual effects of AOT and PWV were analysed. Analysis of individual effects facilitates understanding of the combined effect, as well as allows comparing the results with other studies.

AOT Individual Effect
The spatial distribution of AOT and average values of the AOT effect for the study period are shown in Figure 4. The AOT values shown in Figure 4a reveal a significant spatial variability, with larger AOT over the Northern Hemisphere than over the Southern, and also over land with respect to oceanic regions. Only in oceanic areas close to land regions, with high aerosol concentrations, the AOT values became as large as over land. AOT ranged between 0 and 1.1 approximately. The values lower than 0.1 were found in 51% of the global surface, specifically in latitudes above 60 • and ocean regions. However, the highest values (higher than 0.5) were concentrated in 2% of the global surface. Specifically, they were found above east China (AOT up to 0.8), India and the Arabian and Sahara Deserts. The high concentration of aerosols in these regions is due, on the one hand, to the influence of the deserts and, on the other hand, to the rapid economic expansion of these areas, increasing fossil fuel demand and causing an increase in AOT. Although AOT values are usually lower over oceans, it should be noted that the oceanic areas surrounding these regions are affected by high concentrations of aerosols. An example is the Atlantic Ocean near Central and South America, likely due to the transport of mineral dust [49].
Moderate AOT values were found over some countries in South America, such as Brazil, due to the seasonal deforestation and sugarcane burning [50]. This spatial distribution was also described by other authors (e.g., [51][52][53]) and coincided with the spatial distribution of the AOT effect ( Figure 4b). The absolute AOT effect on SW irradiance at ground level was higher in those regions where the concentration of aerosols was higher. As can be seen in this figure, AOT effect values were negative, indicating a decrease in downwelling SW irradiance, that is, a surface radiative cooling due to the presence of aerosols in the atmosphere. This decrease ranged between 0% and −18%, with an average of −2.3%. In the regions with the highest concentration of aerosols (east China, the Arabian and Sahara Deserts and India, as shown in Figure 4a), the reduction in SW irradiance ranged between −8% and −18 %, while in the regions with AOT lower than 0.1, this reduction did not exceed −2%.
welling SW irradiance, that is, a surface radiative cooling due to the presence of aerosols in the atmosphere. This decrease ranged between 0% and −18%, with an average of −2.3%. In the regions with the highest concentration of aerosols (east China, the Arabian and Sahara Deserts and India, as shown in Figure 4a), the reduction in SW irradiance ranged between −8% and −18 %, while in the regions with AOT lower than 0.1, this reduction did not exceed −2%.   Figure 5 shows the spatial distribution of seasonally averaged AOT effects on downwelling SW irradiance during the period 2000-2019. A high seasonal variation was evident, with more negative values during spring and summer than in autumn and winter. These high AOT effect values in summer and spring were due to a combination of factors, such as an increment in dust outbreaks and forest fires or the lack of precipitation. Precipitation results in a scavenging of atmospheric aerosols reducing the aerosol load and, consequently, the AOT effect. Its impact extends for several days after the precipitation event until new aerosols are incorporated from the soil or by air mass advection. Thus, the lack of precipitation in summer and spring, along with other factors, such as the frequent dust outbreaks, the forest fires and the expansion of the atmospheric boundary layer, contributes to an increase in the AOT effect. It should be noted that the regions with the most negative AOT effect values remained in all seasons. These regions were east China, India, and the Arabian and Sahara Deserts. They coincided with the regions most affected by high AOT values when the global spatial distribution of the AOT effects averaged was analysed ( Figure  4). Unlike what happens in most regions of the planet, in the middle of South America, AOT effects were very high in autumn due to the seasonal deforestation activities [50]. In southwestern Africa, high absolute AOT effects have been found during summer and autumn due to Savanna wildfires [54].
Remote Sens. 2021, 13, 708 7 of 20 and the Arabian and Sahara Deserts. They coincided with the regions most affected by high AOT values when the global spatial distribution of the AOT effects averaged was analysed ( Figure 4). Unlike what happens in most regions of the planet, in the middle of South America, AOT effects were very high in autumn due to the seasonal deforestation activities [50]. In southwestern Africa, high absolute AOT effects have been found during summer and autumn due to Savanna wildfires [54].

PWV Individual Effect
The 20-year mean global spatial distributions of PWV and PWV effects are shown in Figure 6a,b, respectively. Both figures show predominantly zonal patterns along latitude lines. Figure 6a shows that the highest PWV values (around 35-60 mm) were located around the Intertropical Convergence Zone (ITCZ), mainly over the oceans, particularly Western Pacific and sea areas near Indonesia; over land, in central Africa and northern South America. The lowest PWV values (less than 10 mm) were found over regions with latitudes greater than 60°, such as Antarctica and Greenland, and mountainous regions, such as the Himalayas, Andes, and the Rocky Mountains. In the rest of the planet, PWV ranged between 10 and 35 mm. A similar global spatial distribution of PWV was described by Wang et al. [55]. In regard to regions previously identified to have high concentrations of aerosols (east China, the Arabian and Sahara Deserts and India), only east China and

PWV Individual Effect
The 20-year mean global spatial distributions of PWV and PWV effects are shown in Figure 6a,b, respectively. Both figures show predominantly zonal patterns along latitude lines. Figure 6a shows that the highest PWV values (around 35-60 mm) were located around the Intertropical Convergence Zone (ITCZ), mainly over the oceans, particularly Western Pacific and sea areas near Indonesia; over land, in central Africa and northern South America. The lowest PWV values (less than 10 mm) were found over regions with latitudes greater than 60 • , such as Antarctica and Greenland, and mountainous regions, such as the Himalayas, Andes, and the Rocky Mountains. In the rest of the planet, PWV ranged between 10 and 35 mm. A similar global spatial distribution of PWV was described by Wang et al. [55]. In regard to regions previously identified to have high concentrations of aerosols (east China, the Arabian and Sahara Deserts and India), only east China and India showed high PWV values. Conversely, PWV values were low in the Arabian and Sahara Deserts due to the extremely dry atmosphere in these regions.
The global spatial distribution of PWV effects is shown in Figure 6b. As can be seen, the higher the PWV, the higher its effect on solar radiation at the surface. The sign of these effects was negative, indicating a reduction in downwelling SW irradiance due to the presence of water vapour in the atmosphere. The PWV effect ranged between 0% and −18%, with an average of −12.1%. The absolute value of this average was notably higher than the average AOT value mentioned. Therefore, it can be said that, on average, the effect of vapour on downwelling SW irradiance is greater than the effect of aerosols. Figure 7 illustrates the spatial distribution of the PWV effect for different seasons, where the North-South displacement of the General Atmospheric Circulation can be observed. The PWV effect showed less variability over the ocean than over land. Equatorial regions showed very little variation linked to the seasonal displacement of the ITCZ, and the effect of PWV on solar radiation was high all year round. At middle and high latitudes, PWV over land and, consequently, the PWV effect showed a higher seasonal variability due to the relationship between PWV and surface temperature [56,57]. This fact was more evident in regions such as Siberia and Greenland. However, there were also land regions, such as Antarctica, desert and mountainous regions, where the seasonal variability of the PWV effect was low due to the extremely dry atmosphere in these regions all the year round. In general, the most negative values were found in summer, while the least negative values were obtained in winter. Spring and autumn showed similar values. Thus, the area with PWV effect was more negative than the mean value (which is −12.1%), varied from 65% of the global surface in summer to 54% in spring and autumn, and 50% in winter. India showed high PWV values. Conversely, PWV values were low in the Arabian and Sahara Deserts due to the extremely dry atmosphere in these regions. The global spatial distribution of PWV effects is shown in Figure 6b. As can be seen, the higher the PWV, the higher its effect on solar radiation at the surface. The sign of these effects was negative, indicating a reduction in downwelling SW irradiance due to the presence of water vapour in the atmosphere. The PWV effect ranged between 0% and −18%, with an average of −12.1%. The absolute value of this average was notably higher than the average AOT value mentioned. Therefore, it can be said that, on average, the effect of vapour on downwelling SW irradiance is greater than the effect of aerosols.  Figure 7 illustrates the spatial distribution of the PWV effect for different seasons, where the North-South displacement of the General Atmospheric Circulation can be observed. The PWV effect showed less variability over the ocean than over land. Equatorial regions showed very little variation linked to the seasonal displacement of the ITCZ, and the effect of PWV on solar radiation was high all year round. At middle and high latitudes, PWV over land and, consequently, the PWV effect showed a higher seasonal variability due to the relationship between PWV and surface temperature [56,57]. This fact was more evident in regions such as Siberia and Greenland. However, there were also land regions, such as Antarctica, desert and mountainous regions, where the seasonal variability of the PWV effect was low due to the extremely dry atmosphere in these regions all the year round. In general, the most negative values were found in summer, while the least negative values were obtained in winter. Spring and autumn showed similar values. Thus, the area with PWV effect was more negative than the mean value (which is −12.1%), varied from 65% of the global surface in summer to 54% in spring and autumn, and 50% in winter.

AOT-PWV Combined Effect
After analysing the individual effects of aerosols and water vapour on downwelling solar irradiance, the combined effect of both atmospheric quantities was also calculated. The study of the combined effects of AOT and PWV is of great interest because it can also help to understand how aerosols and water vapour interact, affecting each other in the atmosphere, and complement other existing studies in this field. For example, Sun et al. [58] showed the effect of dust on PWV, and Zhao et al. [59] showed the PWV effect on the growth of aerosol particles and then on the AOT. Figure 8a shows the global spatial distribution of averaged AOT-PWV combined effect during the period 2000-2019. As can be seen in this figure, the AOT-PWV effect was also negative, as in the case of the individual effects, indicating a decrease in downwelling SW irradiance. However, the reduction in downwelling SW irradiance had increased, ranging between 0% and −31%, while the reductions for individual effects were limited to the range between 0% and −18%. Likewise, the global average AOT-PWV effect was more negative (−13.9%), while the average AOT effect was −2.3% and the average PWV effect −12.1%. From these values, it is concluded that the combined effect AOT-PWV is less intense than the sum of the individual effects. The difference is due to the fact that the interaction between water vapour and aerosols was not considered when the individual effects wer4e estimated. This interaction determines the aerosol properties in the atmosphere and, therefore, their effects on the downwelling solar radiation on the surface. As mentioned before, this study aimed to calculate the combined AOT-PWV effect (Figure 8), and the individual AOT (Figures 4b and 5) and PWV (Figures 6b and 7) effects on solar radiation, and to demonstrate that the combined effect cannot be solely derived adding the individual effects. For this aim, Figure 8b illustrates the difference between the combined AOT-PWV effect, as calculated by Equation (1), and the sum of the AOT and PWV individual effects, as calculated by Equations (2) and (3), respectively. As can be seen, this difference was not equal to zero in most parts of the globe, ranging from 0 in Antarctica to −4% (highest absolute value) in areas of east China, India and the Sahara Desert. These regions combined the highest AOT values on the planet (Figure 4) with PWV values greater than 30 mm ( Figure 6). Under these conditions, aerosols with a high affinity for water adsorb moisture, increasing their size and their ability to scatter solar radiation, causing an increase in multiple scattering and, consequently, enhancing the downwelling diffuse component of the radiation. This increase in diffuse radiation would compensate for the reduction in direct solar radiation reaching the surface due to the effect of aerosols and water vapour, resulting in a less negative AOT-PWV effect. This comparison between individual and combined effects highlights the importance of analysing the combined effect of aerosols and water vapour on solar radiation, mainly in regions with high concentrations of aerosols and water vapour.
The analysis of the spatial distribution of the AOT-PWV effect showed that the highest reductions in SW irradiance (between −25% and −32%) occurred over east China, India and the surrounding water regions, north and central Africa and South America, coinciding with areas with high influence of aerosols and water vapour simultaneously. There were other regions with high concentrations of only one of them (aerosols or water vapour), and, therefore, their combined effect on solar radiation was not so noticeable. This was the case of the Middle East region (AOT-PWV effect between −15% and −20%), where, although high AOT values were found, PWV values remained very low due to its extremely dry atmosphere. The lowest reductions in SW irradiance (AOT-PWV effect between 0% and −7%) were found in latitudes greater than 60 • and mountainous regions and correspond to areas with a low influence of aerosols, water vapour or both. pour), and, therefore, their combined effect on solar radiation was not so noticeable. This was the case of the Middle East region (AOT-PWV effect between −15% and −20%), where, although high AOT values were found, PWV values remained very low due to its extremely dry atmosphere. The lowest reductions in SW irradiance (AOT-PWV effect between 0% and −7%) were found in latitudes greater than 60° and mountainous regions and correspond to areas with a low influence of aerosols, water vapour or both. Similar to individual effects, the combined effect of AOT and PWV on solar radiation showed seasonal variability (Figure 9). This seasonal variation was higher over land than over ocean. Over ocean, the main variations occurred in regions surrounding areas affected by high concentrations of aerosols, such as regions next to the Arabian and Sahara Deserts. A higher seasonal variability was found over land, with the most negative AOT-PWV effect values in summer. High AOT-PWV effect in summer is mainly due to the increase in atmospheric boundary layer height [60] and the surface temperature due to solar heating favoring the vertical mixing of aerosols (e.g., [61,62]), the increase in dust outbreaks (e.g., [63,64]) and forest fires (e.g., [65]), lack of precipitation (e.g., [66]) and the presence of water vapour in the atmosphere (e.g., [67]).
A spatial analysis of the AOT-PWV effect in Figure 9 revealed that east China was the region with the highest reductions in SW downwelling irradiance throughout the year, with values ranging between −25% and −35%. Over India and north Africa, a strong seasonal pattern was observed, with the highest effect in summer. Other regions with a marked seasonal variability were the north and northwest of China and Canada, where the AOT-PWV effect varied from −20% in summer to values close to 0 in winter.
A spatial analysis of the AOT-PWV effect in Figure 9 revealed that east China was the region with the highest reductions in SW downwelling irradiance throughout the year, with values ranging between −25% and −35%. Over India and north Africa, a strong seasonal pattern was observed, with the highest effect in summer. Other regions with a marked seasonal variability were the north and northwest of China and Canada, where the AOT-PWV effect varied from −20% in summer to values close to 0 in winter.

Temporal Analysis
A temporal analysis was conducted focusing on the detection of trends in the time series of the monthly anomalies of AOT, PWV, AOT effects, PWV effects and AOT-PWV effects, as explained in Section 2. Sen's method was used to estimate the magnitude of the trend, while the Mann-Kendall test was used to calculate its statistical significance. The term significant has been used to express that the trend is statistically significant (confidence level of 95%). Figure 10 shows the monthly evolution of these anomalies. AOT monthly anomalies showed an increase over the years, with significant slope equal to 0.0002 per year ( Figure  10a). In the case of PWV, monthly anomalies also showed an increase over the years, with

Temporal Analysis
A temporal analysis was conducted focusing on the detection of trends in the time series of the monthly anomalies of AOT, PWV, AOT effects, PWV effects and AOT-PWV effects, as explained in Section 2. Sen's method was used to estimate the magnitude of the trend, while the Mann-Kendall test was used to calculate its statistical significance. The term significant has been used to express that the trend is statistically significant (confidence level of 95%). Figure 10 shows the monthly evolution of these anomalies. AOT monthly anomalies showed an increase over the years, with significant slope equal to 0.0002 per year (Figure 10a). In the case of PWV, monthly anomalies also showed an increase over the years, with a slope equal to 0.0014 per year (Figure 10b), although the slope was not significant. Figure 10c shows the evolution of the AOT effect and PWV effect anomalies. The higher the value of AOT or PWV, the higher its effects on solar radiation at the surface, revealing a decrease throughout the study period. The slopes were significant and presented negative values (−0.0038% and −0.0020% per year) since the AOT effect and PWV effect values were negative. Figure 10c also shows the temporal evolution of AOT-PWV effect anomalies, with a significant trend equal to −0.0052% per year. This negative trend also indicated an increase in the absolute combined effect on solar radiation at the surface since the AOT-PWV effect was also negative.
In addition to the detection of the general trends in the AOT and PWV anomalies and the anomalies of their effects, it is interesting to analyse its spatial distribution throughout the globe. Figure 11 shows the spatial distribution of the AOT linear trend per decade ( Figure 11a) and its corresponding AOT effect (Figure 11b). Pixels with levels of significance p < 0.05 are marked with a black dot. The results obtained show an average of the linear trend per decade over the global of 0.0017 for AOT, during the period 2000-2019, which represents an absolute increase of almost 0.0034. A global positive trend in AOT has also been obtained by other authors (e.g., [68,69]). Their values differ depending on the data period used. Thus, for example, Zhang and Reid [68]  throughout the study period. The slopes were significant and presented negative values (−0.0038% and −0.0020% per year) since the AOT effect and PWV effect values were negative. Figure 10c also shows the temporal evolution of AOT-PWV effect anomalies, with a significant trend equal to −0.0052% per year. This negative trend also indicated an increase in the absolute combined effect on solar radiation at the surface since the AOT-PWV effect was also negative. In addition to the detection of the general trends in the AOT and PWV anomalies and the anomalies of their effects, it is interesting to analyse its spatial distribution throughout the globe. Figure 11 shows the spatial distribution of the AOT linear trend per decade ( Figure 11a) and its corresponding AOT effect (Figure 11b). Pixels with levels of significance p < 0.05 are marked with a black dot. The results obtained show an average of the linear trend per decade over the global of 0.0017 for AOT, during the period 2000-2019, which represents an absolute increase of almost 0.0034. A global positive trend in AOT has also been obtained by other authors (e.g., [68,69]). Their values differ depending on the data period used. Thus, for example, Zhang and Reid [68]   Approximately 72% of the global surface showed a positive trend of AOT, and it was significant at a 5% significance level in 50% of that surface. The highest values of significant positive trend (up to 0.13 per decade) were found in India and its surrounding oceanic areas. This finding is in line with recent studies that have reported significant positive trends of AOT over the Indian region [70,71]. This trend is associated with an increase in urban/industrial pollution due to the rapid economic expansion [72]. There was also a significant positive trend, although with lower values (up to 0.05), in some areas of the Sahara Desert, northern South America, and the region of the Philippines, Malaysia or Indonesia, as well as in the west of the Middle East region. The significant positive trend of AOT in the western Middle East is due to an increase in the frequency and intensity of dust outbreaks from the Arabian Desert [73][74][75]. In contrast, significant negative trends (up to −0.05 per decade) were found in eastern United States, Europe, the Mediterranean basin, East Asia and South America. This negative trend could be due to the effective enforcement of emission control policies in these regions [76][77][78][79][80]. This spatial distribution agrees with the spatial distribution described by other authors (e.g., [51,53,69,76]).
globe was equal to −0.029% per decade during the period 2000-2019, resulting in an absolute AOT effect increase of 0.057%. A negative trend was also estimated by Subba et al. [23] for aerosol effects on solar radiation. The spatial distribution was similar to that of AOT. The more positive the AOT trend, the more negative the AOT effect trend. The most negative AOT effect trends were found in India, likely due to the increasing trend of anthropogenic aerosols in that region [81]. In the case of the AOT effect (Figure 11b), the average of the linear trend over the globe was equal to −0.029% per decade during the period 2000-2019, resulting in an absolute AOT effect increase of 0.057%. A negative trend was also estimated by Subba et al. [23] for aerosol effects on solar radiation. The spatial distribution was similar to that of AOT. The more positive the AOT trend, the more negative the AOT effect trend. The most negative AOT effect trends were found in India, likely due to the increasing trend of anthropogenic aerosols in that region [81]. Figure 12 shows the global spatial distribution of the PWV linear trend per decade (Figure 12a) and the PWV effect (Figure 12b). An overall global trend of 0.023 mm per decade during the period 2000-2019 was found for PWV, which represented an absolute increase of almost 0.046 mm. This positive sign of the trend is in line with other studies (e.g., [31,82]) and is found in almost 60% of the world. The regions with significant positive trends were distributed all over the world. The most positive trends were found in South America. Significant negative trends were obtained in the easternmost part of South America, South Africa, the north of the Indian Ocean, northern Australia and Pacific Ocean regions. While the negative trend values for these last four regions are in line with the results reported by Chen et al. [31], the spatial distribution differs from that shown by Mieruch et al. [83], possibly due to the different data period used in the studies. Thus, while Mieruch et al. [83] analysed only ten years (1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005), in the current study, the data period was extended to two decades (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018)(2019). It is interesting to note that Polar Remote Sens. 2021, 13, 708 14 of 20 Regions behave in the opposite way, i.e., while PWV was increasing in the Arctic area, it decreased in Antarctica. This opposite pattern is consistent with other studies that indicate the Arctic is warming faster than the global warming average [84]. Thus, this increase in surface temperature in the Arctic would yield more evaporation and, consequently, more PWV [85]. Once the trends of the individual effects were analysed, the trends of the combined effect was studied. A global trend in the combined AOT-PWV effect equal to −0.053% per decade during the period 2000-2019 was found, representing a decrease of 0.11% along the period of study. This negative trend was found in 66% of the world's area, and it was significant at a 5% significance level in 38% of that surface. Figure 13 shows the global spatial distribution of the linear trend per decade in the AOT-PWV effect. The most negative significant trend values were obtained in India. These values were related to the positive trend of AOT, likely due to the rapid economic expansion. Ocean waters surrounding India, central Africa, northwest Pacific, northern South America, and the region of the Philippines, Malaysia or Indonesia showed significant negative trends, although their absolute values were lower. A negative trend was also obtained in the Arctic due to the influence of the negative trend of PWV effect. Conversely, significant positive trends were found in the eastern United States, South America, Europe, the Mediterranean Sea, southern Africa, Antarctic, Asia, Australia, and some specific areas of the Atlantic Ocean (regions where significant positive trends of AOT effect were obtained). In a recent study, Obregón et al. [29] also reported significant positive trends of the AOT-PWV effect over the Mediterranean Sea. The positive trend in these regions indicates that emission control policies are working, something that is not happening in India due to the great economic expansion the country is experiencing. In the case of the PWV effect (Figure 12b), the global trend was equal to −0.019% per decade during the period 2000-2019, resulting in an absolute PWV effect increase of 0.039%. The spatial distribution of the PWV effect was similar, in most of the study area, to that of PWV. The more positive the PWV trend, the more negative the PWV effect trend becomes. This similarity was not fulfilled in the maritime areas of Polar Regions. In these areas, the absolute PWV trends were not high; however, the PWV effect trends showed high absolute values. This lack of similarity is due to the fact that PWV in these regions was very small, and any small alteration can result in a large variation in the PWV effect, as can be seen in Figure 3.
Once the trends of the individual effects were analysed, the trends of the combined effect was studied. A global trend in the combined AOT-PWV effect equal to −0.053% per decade during the period 2000-2019 was found, representing a decrease of 0.11% along the period of study. This negative trend was found in 66% of the world's area, and it was significant at a 5% significance level in 38% of that surface. Figure 13 shows the global spatial distribution of the linear trend per decade in the AOT-PWV effect. The most negative significant trend values were obtained in India. These values were related to the positive trend of AOT, likely due to the rapid economic expansion. Ocean waters surrounding India, central Africa, northwest Pacific, northern South America, and the region of the Philippines, Malaysia or Indonesia showed significant negative trends, although their absolute values were lower. A negative trend was also obtained in the Arctic due to the influence of the negative trend of PWV effect. Conversely, significant positive trends were found in the eastern United States, South America, Europe, the Mediterranean Sea, southern Africa, Antarctic, Asia, Australia, and some specific areas of the Atlantic Ocean (regions where significant positive trends of AOT effect were obtained). In a recent study, Obregón et al. [29] also reported significant positive trends of the AOT-PWV effect over the Mediterranean Sea. The positive trend in these regions indicates that emission control policies are working, something that is not happening in India due to the great economic expansion the country is experiencing.

Conclusions
The global spatial and temporal distribution of the combined effects of AOT and PWV on solar radiation at the surface were calculated and analysed. In addition, the individual effects of aerosols and water vapour were also calculated to compare the results with other analyses at a global scale. For this purpose, the methodology recently developed and validated by Obregón et al. [30] was applied at the global scale for the first time. The analysis of the spatial distribution of individual and combined AOT and PWV effects showed a close link with spatial distributions of AOT and PWV. These effects were negative, indicating a decrease in downwelling SW irradiance due to the presence of aerosols and water vapour in the atmosphere. AOT-PWV effects ranged between 0% and −31%, while the reductions for individual effects were limited to the range between 0% and −18%. Likewise, the global average AOT-PWV effect is more negative (−13.9%), while the average AOT effect was −2.3% and of the average PWV effect −12.1%. From these values, it is concluded that the combined effect AOT-PWV is less intense than the sum of the individual effects, mainly in regions with the highest concentrations of aerosols and water vapour. The difference between the sum of the individual effects and the combined effect is due to the lack of consideration of the interaction between water vapour and aerosols when the individual effects are calculated. The most negative AOT effects were found above east China, India, and the Arabian and Sahara Desert due to the influence of the deserts and the rapid economic expansion. While the most negative PWV effect was located around the Intertropical Convergence Zone. The most negative AOT-PWV effects were found in areas with a high influence of aerosols and water vapour simultaneously. High seasonal variability was found, with the most negative values during summer due

Conclusions
The global spatial and temporal distribution of the combined effects of AOT and PWV on solar radiation at the surface were calculated and analysed. In addition, the individual effects of aerosols and water vapour were also calculated to compare the results with other analyses at a global scale. For this purpose, the methodology recently developed and validated by Obregón et al. [30] was applied at the global scale for the first time. The analysis of the spatial distribution of individual and combined AOT and PWV effects showed a close link with spatial distributions of AOT and PWV. These effects were negative, indicating a decrease in downwelling SW irradiance due to the presence of aerosols and water vapour in the atmosphere. AOT-PWV effects ranged between 0% and −31%, while the reductions for individual effects were limited to the range between 0% and −18%. Likewise, the global average AOT-PWV effect is more negative (−13.9%), while the average AOT effect was −2.3% and of the average PWV effect −12.1%. From these values, it is concluded that the combined effect AOT-PWV is less intense than the sum of the individual effects, mainly in regions with the highest concentrations of aerosols and water vapour. The difference between the sum of the individual effects and the combined effect is due to the lack of consideration of the interaction between water vapour and aerosols when the individual effects are calculated. The most negative AOT effects were found above east China, India, and the Arabian and Sahara Desert due to the influence of the deserts and the rapid economic expansion. While the most negative PWV effect was located around the Intertropical Convergence Zone. The most negative AOT-PWV effects were found in areas with a high influence of aerosols and water vapour simultaneously. High seasonal variability was found, with the most negative values during summer due to a combination of factors, such as an increment in dust outbreaks and forest fires, the lack of precipitation, and the increase in the surface temperature. This spatial analysis has contributed to a better understanding of the individual and combined effects of AOT and PWV on downwelling SW irradiance and highlights the importance of analysing the combined effects, mainly in regions with high concentrations of aerosols and water vapour.
The analysis of the temporal distribution focused on the detection of trends in the AOT and PWV anomalies and the anomalies of their effects. The results obtained show positive global trends per decade for AOT and PWV, with values equal to 0.0017 and 0.0234 mm, respectively. Consequently, significant negative global trends were found for individual AOT and PWV effects and combined effects (−0.029%, −0.019% and −0.053% per decade). These negative trends also indicate an increase in the absolute effects on solar radiation at the surface since effects were negative. The spatial distribution of AOT effect trends showed that the most significant negative values were found in India and its surrounding oceanic areas due to the rapid economic expansion that this region has experienced. In contrast, significant positive trends were found in the eastern United States, Europe, South America or East Asia due to the effective enforcement of emission control policies. In the case of the PWV effect, the most negative trends were found in South America, central Africa, western Unites States and the maritime area of Arctic, and the most positive trends were found in South Africa, northern Australia and the maritime area of the Antarctic. The opposite pattern among the Polar Regions is consistent with temperature studies that suggest the Arctic is warming faster than the global average. The global spatial distribution of the AOT-PWV effect trend shows that, in general, trends in the combined AOT-PWV effect are directly related to the trends in the individual AOT and PWV effects. However, a comprehensive description of the aerosol and water vapour effects is only provided by analysing the trend of the combined effects since both constituents are simultaneous present in the atmosphere. Thus, this combined effect is the suitable one to be considered for evaluating climatic effects. Furthermore, the study of the combined effects of AOT and PWV may be used in future research to address how aerosols and water vapour can affect each other in the atmosphere. This study focused on the individual and combined effects of aerosols and water vapour on downward SW radiation, quantifying these variables using AOT and PWV data provided by satellite products. In particular, the CERES SYN1deg Ed4 product was used. Therefore, this study is limited to cloud-free cases since AOT and PWV can only be obtained from the satellite under these conditions. However, it should be noted that the role of PWV also extends to cases with clouds, as it seems to be correlated with the likelihood of clouds, which definitely have an impact on the downward SW irradiance. In addition, aerosols contribute positively to cloud formation by acting as condensation and ice nuclei. Funding: This work was partially supported by the FCT (Fundação para a Ciência e a Tecnologia) through the grant SFRH/BPD/86498/2012 with national funding. The work is co-funded by the European Union through the European Regional Development Fund, included in the COM-PETE 2020 (Operational Program Competitiveness and Internationalization) through the ICT project (UIDB/04683/2020) with the reference POCI-01-0145-FEDER-007690 and also through TOMAQAPA (PTDC/CTAMET/29678/2017) and CILIFO (0753_CILIFO_5_E) projects. The work has been also partially funded by FEDER/Ministerio de Ciencia, Innovación y Universidades-Agencia Estatal de Investigación through project RTI 2018-097332-B-C22, and by Junta de Extremadura-FEDER through project GR18097.