Evaluation of Precipitable Water Vapor from Five Reanalysis Products with Ground-Based GNSS Observations

At present, the global reliability and accuracy of Precipitable Water Vapor (PWV) from different reanalysis products have not been comprehensively evaluated. In this study, PWV values derived by 268 Global Navigation Satellite Systems (GNSS) stations around the world covering the period from 2016 to 2018 are used to evaluate the accuracies of PWV values from five reanalysis products. The temporal and spatial evolution is not taken into account in this analysis, although the temporal and spatial evolution of atmospheric flows is one of the most important information elements available in numerical weather prediction products. The evaluation results present that five reanalysis products with PWV accuracy from high to low are in the order of the fifth generation of European Centre for Medium-Range Weather Forecasts (ECMWF) Reanalysis (ERA5), ERA-Interim, Japanese 55-year Reanalysis (JRA-55), National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR), and NCEP/DOE (Department of Energy) according to root mean square error (RMSE), bias and correlation coefficient. The ERA5 has the smallest RMSE value of 1.84 mm, while NCEP/NCAR and NCEP/DOE have bigger RMSE values of 3.34 mm and 3.51 mm, respectively. The findings demonstrate that ERA5 and two NCEP reanalysis products have the best and worst performance, respectively, among five reanalysis products. The differences in the accuracy of the five reanalysis products are mainly attributed to the differences in the spatial resolution of reanalysis products. There are some large absolute biases greater than 4 mm between GNSS PWV values and the PWV values of five reanalysis products in the southwest of South America and western China due to the limit of terrains and fewer observations. The accuracies of five reanalysis products are compared in different climatic zones. The results indicate that the absolute accuracies of five reanalysis products are highest in the polar regions and lowest in the tropics. Furthermore, the effects of different seasons on the accuracies of five reanalysis products are also analyzed, which indicates that RMSE values of five reanalysis products in summer and in winter are the largest and the smallest in the temperate regions. Evaluation results from five reanalysis products can help us to learn more about the advantages and disadvantages of the five released water vapor products and promote their applications.


Introduction
The change of atmospheric water vapor content has an important impact on weather prediction and global climate change [1,2]. In general, each 1 K rise in temperature can make the atmospheric

Global Navigation Satellite Systems (GNSS) Precipitable Water Vapor (PWV)
To obtain the PWV value, Zenith Wet Delay (ZWD) can be calculated by subtracting Zenith Hydrostatic Delay (ZHD) from Zenith Tropospheric Delay (ZTD). There is no way to provide accurate meteorological data with respect to most GNSS stations because we find that about 65% of GNSS stations do not install meteorological observation instruments or the equipment and instruments of a few GNSS stations are damaged according to meteorological data available provided by the CDDIS (Crustal Dynamics Data Information System). Therefore, ZHD and VMF1 (Vienna Mapping Function) coefficients can be obtained by the gridded VMF1 from the analysis data of the ECMWF. ZWD can be calculated by using Bernese5.2 software [44].
The conversion formula from ZWD to PWV is as follows: where Π is the water vapor conversion coefficient, which is a function of weighted mean temperature Tm directly obtained from the Global Pressure and Temperature 3 (GPT3) model by inputting the coordinate value and time [45].

Reanalysis of PWV Products
The accuracies of PWV from five global reanalysis products are verified in this paper, including ERA5, ERA-Interim, Japanese 55-year Reanalysis (JRA-55), NCEP/National Center for Atmospheric Research (NCAR), and NCEP/DOE reanalysis. The information of five global reanalysis products is summarized in Table 1 and the details are described below.

ERA-Interim
The ERA-Interim is based on the reanalysis of ECMWF ERA-40 from 1979 to 31st August 2019. Compared to ERA-40 reanalysis, there are many improvements in the assimilation system and atmospheric model. The ERA-Interim reanalysis products are available every 6 hours and provide atmospheric and surface parameters [46].

ERA5
ERA5 is the latest reanalysis product of ECMWF, providing hourly atmospheric reanalysis products. ERA5 reanalysis products are produced using a high-resolution of 0.25 • × 0.25 • (atmosphere) and are available in the Climate Data Store (CDS). The ERA5 reanalysis will be completed by mid-2020, when temporal coverage is from 1950 to the present. ERA5 reanalysis products use HadISST.2, reprocess the climate data records of ECMWF and implement RTTOV11 radiative transfer, which improves the accuracy of ERA5 reanalysis products compared with ERA-Interim reanalysis products [47].

Japanese 55-year Reanalysis (JRA-55)
JRA-55 is the second reanalysis project and uses a more sophisticated Data Assimilation (DA) system to perform data assimilation [48]. The period of JRA-55 reanalysis is covered from 1958 to the present. The JRA-55 reanalysis project alleviated many of the deficiencies exhibited by JRA-25 using an improved DA system to produce high-quality reanalysis products [49].

National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR)
The National Centers for Environmental Prediction (NCEP)/NCAR reanalysis project is carried out by using a typical analysis/forecast system to assimilate observation data from 1948 to the present. The NCEP/NCAR reanalysis products, with a temporal resolution of four times per day, are available from Physical Sciences Laboratory (PSL) [50].

NCEP/Department of Energy (DOE)
The successor to NCEP/NCAR reanalysis, NCEP/DOE reanalysis (Reanalysis-2), is available, going back to 1979 from the present. The main goal of NCEP/NCAR reanalysis is to improve the accuracy of NCEP/NCAR reanalysis products by updating physical processes parameters and correcting the errors [51].

Comparison between GNSS PWV Derived from Different Pressure Values from Different Reanalysis Pressure; Products and PWV Derived from Different Reanalysis Water Vapor Products
We calculate the ZHD with pressure values with respect to five reanalysis products, and convert the ZTD from Bernese into PWV. GNSS PWV values derived from ZTD and different pressure values with respect to five reanalysis products are compared with PWV values derived from respective reanalysis products, which is shown in Figure 1.  Figure 1. depicts that there is little difference in the accuracy evaluation of PWV values from five reanalysis water vapor products with GNSS PWV values from different pressure values with respect to five reanalysis products. GNSS PWV values derived from Bernese ZTD and ZHD of pressure values from ECMWF (ERA5 and Interim) reanalysis products are in better agreement with PWV values derived from five reanalysis water vapor products. Therefore, in this study, pressure values from ECMWF reanalysis pressure products are used to calculate GNSS PWV values.

Elevation Correction
ERA5 is a 1-hour temporal resolution reanalysis product, and the other four types of reanalysis products are 6-hour temporal resolution data. However, in this study, a 2-hour temporal resolution GNSS data is used. When the temporal resolution is interpolated, the accuracies of the reanalysis products may be reduced, affecting the accuracy evaluation of reanalysis products. Therefore, in order to better evaluate the accuracies of the reanalysis products, we do not interpolate in terms of temporal resolution. PWV data from reanalysis and from GNSS at corresponding time are selected. Since altitude has a great influence on PWV value, it is necessary to correct the PWV value from reanalysis products to the PWV value of GNSS station altitude, which can better evaluate the accuracy of PWV from reanalysis products using GNSS PWV [52].
The 0 PWV is water vapor value of the reanalysis products height, PWV is water vapor value from 0 PWV of reanalysis products height correction to GNSS station height , 2 C is a constant of 0.439, and h is the difference of the height of the reanalysis products minus the height of the GNSS station.

Statistical Indicators
GNSS PWV values are used to evaluate the PWV accuracy of five reanalysis products. Statistical indicators to measure the quality of reanalysis products include correlation coefficient R, bias, root

Elevation Correction
ERA5 is a 1-h temporal resolution reanalysis product, and the other four types of reanalysis products are 6-h temporal resolution data. However, in this study, a 2-h temporal resolution GNSS data is used. When the temporal resolution is interpolated, the accuracies of the reanalysis products may be reduced, affecting the accuracy evaluation of reanalysis products. Therefore, in order to better evaluate the accuracies of the reanalysis products, we do not interpolate in terms of temporal resolution. PWV data from reanalysis and from GNSS at corresponding time are selected. Since altitude has a great influence on PWV value, it is necessary to correct the PWV value from reanalysis products to the PWV value of GNSS station altitude, which can better evaluate the accuracy of PWV from reanalysis products using GNSS PWV [52].
The PWV 0 is water vapor value of the reanalysis products height, PWV is water vapor value from PWV 0 of reanalysis products height correction to GNSS station height, C 2 is a constant of 0.439, and ∆h is the difference of the height of the reanalysis products minus the height of the GNSS station.

Statistical Indicators
GNSS PWV values are used to evaluate the PWV accuracy of five reanalysis products. Statistical indicators to measure the quality of reanalysis products include correlation coefficient R, bias, root Remote Sens. 2020, 12, 1817 6 of 18 mean square error (RMSE), relative bias (R-bias), and relative root mean square error (R-RMSE) [34,53]. The corresponding formulas are as follows: where N is the number of PWV pairs, PWV gi is PWV value of GNSS, PWV ri is PWV value of reanalysis products, PWV g is the mean PWV value of GNSS, and PWV r is the mean PWV value of reanalysis products.

The Influence of Height Difference on the Accuracy Evaluation of Five Reanalysis Products
Although Equation 2 is used to correct the influence of height difference between GNSS and reanalysis products on the accuracy evaluation of five reanalysis products, the height difference may still have an important influence on the accuracy evaluation, which is shown in Figure 2. mean square error (RMSE), relative bias (R-bias), and relative root mean square error (R-RMSE) [34,53]. The corresponding formulas are as follows: where N is the number of PWV pairs, gi PWV is PWV value of GNSS, ri PWV is PWV value of reanalysis products, g PWV is the mean PWV value of GNSS, and r PWV is the mean PWV value of reanalysis products.

The Influence of Height Difference on the Accuracy Evaluation of Five Reanalysis Products
Although Equation 2 is used to correct the influence of height difference between GNSS and reanalysis products on the accuracy evaluation of five reanalysis products, the height difference may still have an important influence on the accuracy evaluation, which is shown in Figure 2.  a big influence on accuracy evaluation of five reanalysis water vapor products. The finding indicates that large height difference between GNSS height and reanalysis products height has an important influence on accuracy evaluation of five reanalysis water vapor products.

Evaluation of PWV Values From Five Reanalysis Products Around the World
The accuracies of five reanalysis products around the world are evaluated with respect to 268 GNSS stations, covering the period from 2016 to 2018. The information of accuracies of five reanalysis products around the world is shown in Table 2, and R-bias refers to mean absolute relative bias. What is more, the distributions of biasand RMSE of five reanalysis products are depicted in Figures 3 and 4, and the column graph of RMSE of five reanalysis products around the world is presented in Figure 5.  Figure 2 depicts that the height difference less than 1000 m has a small influence on accuracy evaluation of five reanalysis water vapor products, while the height difference greater than 1000 m has a big influence on accuracy evaluation of five reanalysis water vapor products. The finding indicates that large height difference between GNSS height and reanalysis products height has an important influence on accuracy evaluation of five reanalysis water vapor products.

Evaluation of PWV Values From Five Reanalysis Products Around the World
The accuracies of five reanalysis products around the world are evaluated with respect to 268 GNSS stations, covering the period from 2016 to 2018. The information of accuracies of five reanalysis products around the world is shown in Table 2, and R-bias refers to mean absolute relative bias. What is more, the distributions of biasand RMSE of five reanalysis products are depicted in Figure 3 and Figure 4, and the column graph of RMSE of five reanalysis products around the world is presented in Figure 5.       Figure 3 shows the number of the large absolute bias values greater than 3 mm of two NCEP reanalysis products are more than those of the other three kinds of reanalysis products, whereas the number of the large absolute bias values greater than 3 mm of ERA5 is the least among five reanalysis products. Figure 4 depicts that the distributions of RMSE of five reanalysis products, the larger RMSE values of which are mainly located at low latitudes, in southwestern South America and western China, and are similar to those of bias around the world.   As shown in Table 2, five reanalysis products with PWV accuracy from high to low are in the order of ERA5, ERA-Interim, JRA-55, NCEP/NCAR, and NCEP/DOE around the world in accordance with RMSE, bias, R-RMSE, R-bias, and R. The accuracy difference among five reanalysis products may be directly related to the difference of spatial resolution. The higher the spatial resolutions of reanalysis products are, the higher are the accuracies of reanalysis products. The ERA5 has the smallest RMSE, mean absolute bias, R-RMSE, R-bias, and the largest correlation coefficient, which are 1.84 mm, 0.67 mm, 12.01%, 5.21%, and 0.97, which indicate that ERA5 are generally in better agreement with GNSS compared with the other reanalysis products. ERA5 assimilates more recent satellite observations, with the latest data assimilation method adopted, which improves the performance of ERA5.
The RMSE values, mean absolute bias values, and R-RMSE values of NCEP/NCAR and NCEP/DOE are greater than 3 mm, 1 mm, and 20%, respectively. The correlation coefficients of NCEP/NCAR and NCEP/DOE are 0.91 and 0.90, but those of the other reanalysis products are beyond 0.95. The results indicate that the accuracies of two NCEP reanalysis products are poorer than those of the other reanalysis products, which may be due to the use of old assimilation systems and atmospheric physical models, and fewer recent observations being assimilated.  According to Figure 5, the percentages of ERA5 and ERA-Interim, the RMSE values of which are less than 3 mm, exceed 80%. The percentage of ERA5 whose RMSE values are less than 3 mm is more than 95%. The percentages of two NCEP reanalysis products with the RMSE of greater than 5   Figure 3 shows the number of the large absolute bias values greater than 3 mm of two NCEP reanalysis products are more than those of the other three kinds of reanalysis products, whereas the number of the large absolute bias values greater than 3 mm of ERA5 is the least among five reanalysis products. Figure 4 depicts that the distributions of RMSE of five reanalysis products, the larger RMSE values of which are mainly located at low latitudes, in southwestern South America and western China, and are similar to those of bias around the world.
As shown in Table 2, five reanalysis products with PWV accuracy from high to low are in the order of ERA5, ERA-Interim, JRA-55, NCEP/NCAR, and NCEP/DOE around the world in accordance with RMSE, bias, R-RMSE, R-bias, and R. The accuracy difference among five reanalysis products may be directly related to the difference of spatial resolution. The higher the spatial resolutions of reanalysis products are, the higher are the accuracies of reanalysis products. The ERA5 has the smallest RMSE, mean absolute bias, R-RMSE, R-bias, and the largest correlation coefficient, which are 1.84 mm, 0.67 mm, 12.01%, 5.21%, and 0.97, which indicate that ERA5 are generally in better agreement with GNSS compared with the other reanalysis products. ERA5 assimilates more recent satellite observations, with the latest data assimilation method adopted, which improves the performance of ERA5.
The RMSE values, mean absolute bias values, and R-RMSE values of NCEP/NCAR and NCEP/DOE are greater than 3 mm, 1 mm, and 20%, respectively. The correlation coefficients of NCEP/NCAR and NCEP/DOE are 0.91 and 0.90, but those of the other reanalysis products are beyond 0.95. The results indicate that the accuracies of two NCEP reanalysis products are poorer than those of the other reanalysis products, which may be due to the use of old assimilation systems and atmospheric physical models, and fewer recent observations being assimilated.
According to Figure 5, the percentages of ERA5 and ERA-Interim, the RMSE values of which are less than 3 mm, exceed 80%. The percentage of ERA5 whose RMSE values are less than 3 mm is more than 95%. The percentages of two NCEP reanalysis products with the RMSE of greater than 5 mm exceeds 10%. The results indicate that the accuracies of ERA5, ERA-Interim, and JRA55 are higher than those of two NCEP reanalysis products, and ERA5 has the highest accuracy among five reanalysis products.

Evaluation of PWV Values from Five Reanalysis Products in Different Climatic Zones
As can be seen from Figures 3 and 4, biasand RMSE of five reanalysis products are different at different latitudes. According to the distribution of solar heat on the Earth's surface, we have divided the globe into three different climatic zones around the world, which include tropical regions, temperate regions, and polar regions. The information of accuracies of five reanalysis products in different climatic zones is shown in Table 3, and R-bias refers to mean absolute relative bias. As shown in Table 3, in the polar regions, higher absolute accuracy illustrates that the PWV values of five reanalysis products agree well with those of GNSS, but the relative accuracies of five reanalysis products are poor. In the tropical regions, the PWV values of five reanalysis products have a large difference compared to those of GNSS. The possible reasons for low accuracy of the five reanalysis products in the tropical regions can be explained as follows. On the one hand, in the tropical regions, the number of radiosondes is limited. On the other hand, the quality of satellite remote-sensing observation is low due to the thick cloud cover, and the super refraction, which is an abnormally rapid decrease in the refractive index with height in the atmosphere resulting in the anomalous propagation of radio waves, affecting the quality of radio occultation data observation in the tropical regions. Table 3 presents five reanalysis products with RMSE values from small to large are in the order of ERA5, ERA-Interim, JRA-55, NCEP/NCAR, and NCEP/DOE even in different climatic zones. The RMSE values of ERA5 are smaller than those of the other four reanalysis products in different climatic zones. In addition, the RMSE values of NCEP/NCAR and NCEP/DOE reanalysis products are both less than 1.8 mm in the polar regions, which demonstrates the two NCEP reanalysis products in the polar regions have a better performance compared with the other climatic zones. However, the RMSE values and mean absolute bias values of NCEP/NCAR and NCEP/DOE reanalysis products are both larger than 4.5 mm and 1.6 mm, respectively, and the correlation coefficients of NCEP/NCAR and NCEP/DOE reanalysis products are both less than 0.9 in the tropical regions. The results show that the accuracies of two NCEP reanalysis products are far lower than those of the other three reanalysis products in the tropical regions, which is mainly attributed to the lower spatial resolution.

Seasonal Variations of the Accuracy of Five Reanalysis Water Vapor Products in the Different Climatic Zones
In order to better analyze seasonal variations of the accuracy of five reanalysis water vapor products in the different climatic zones, we take ERA5 as an example and analyze the variations of RMSE values with the change of seasons in different climatic zones. As we know, there are four typical seasons in the temperate regions including spring, summer, autumn and winter. However, seasons in the tropical regions which are controlled by precipitation variability are divided into the dry season and wet season [54], and seasons in the polar regions mainly include the warm season and cold season [55].
Due to the small number of GNSS stations in the polar regions, all GNSS stations are selected in the polar regions to analyze the seasonality of RMSE for reanalysis products. There are more GNSS stations in the tropical and temperate regions, GNSS stations are selected in these regions according to the latitudinal belts. The information of GNSS station location selected is shown in Figure 6. The seasonality of RMSE from ERA5 is shown in Figure 7, and the seasonality of RMSE from the other four kinds of reanalysis products in the temperate regions is depicted in Figure 8.
In order to better analyze seasonal variations of the accuracy of five reanalysis water vapor products in the different climatic zones, we take ERA5 as an example and analyze the variations of RMSE values with the change of seasons in different climatic zones. As we know, there are four typical seasons in the temperate regions including spring, summer, autumn and winter. However, seasons in the tropical regions which are controlled by precipitation variability are divided into the dry season and wet season [54], and seasons in the polar regions mainly include the warm season and cold season [55].
Due to the small number of GNSS stations in the polar regions, all GNSS stations are selected in the polar regions to analyze the seasonality of RMSE for reanalysis products. There are more GNSS stations in the tropical and temperate regions, GNSS stations are selected in these regions according to the latitudinal belts. The information of GNSS station location selected is shown in Figure 6. The seasonality of RMSE from ERA5 is shown in Figure 7, and the seasonality of RMSE from the other four kinds of reanalysis products in the temperate regions is depicted in Figure 8.    KELY  KIRU  KIR0  INVK  NRIL  TRO1  SCOR  TIXG  TIXI  THU2  THU3  NYA1  NYAL  ALRT  SCTB  MCM4  SYOG  DAV1  MAW1 KELY  KIRU  KIR0  INVK  NRIL  TRO1  SCOR  TIXG  TIXI  THU2  THU3  NYA1  NYAL  ALRT  SCTB  MCM4  SYOG  DAV1  MAW1     As can be seen from Figure 7, RMSE values of ERA5 in the polar regions have no obvious seasonal variations. The RMSE values of ERA5 in the wet season are slightly greater than those in the dry season in some GNSS stations, which demonstrates that some seasonal variations are not significant in the tropical regions.
However, RMSE values of ERA5 in the temperate regions have obvious seasonality. Figures 8 and 7b show that RMSE values of five reanalysis products in summer and in winter are the largest and the smallest in the temperate regions, respectively. The main reason could be that the PWV value and its uncertainty are relatively large in summer and small in winter in the temperate regions.

Conclusions and Discussion
Firstly, PWV values derived from 268 GNSS stations around the world from 2016 to 2018 are used to verify the accuracies of the five reanalysis water vapor products on a global scale. We can draw the conclusion that the five reanalysis products with PWV accuracy from high to low are ERA5, ERA-Interim, JRA-55, NCEP/NCAR, and NCEP/DOE. RMSE, bias, and correlation coefficient of ERA5 are 1.84 mm, 0.67 mm, and 0.97 respectively, and the percentage of RMSE greater than 3 mm is over 95%, which indicates that ERA5 has the best performance among the five reanalysis products. The reasons for high accuracy of ERA5 are attributed to two respects. On the one hand, ERA5 has higher spatial and temporal resolution. On the other hand, ERA5 assimilates more recent satellite observations with an updated assimilation method, which improves its performance. Two NCEP reanalysis products have larger RMSE and bias, which are beyond 3 mm and 1 mm, respectively, and the percentage of two NCEP reanalysis products, whose RMSE values are larger than 3 mm, exceeds 10%. These results show that the accuracies of two NCEP reanalysis products are lower than those of the other three reanalysis products. The performance differences between two NCEP reanalysis products and the other reanalysis may be due to the differences in spatial resolution of reanalysis products. In addition, there are some large deviations with respect to five reanalysis products in the southwest of South America and western China. The main reason why there are large RMSE values in southwestern South America and western China-where the terrain is complex and mountainous and where few observations are available-is that reanalysis products cannot represent the strong gradients that sometimes are caused by the variable terrain.
Secondly, the global accuracy distribution of five reanalysis products varies with different climatic zones, including tropical, temperate, and polar regions. The RMSE value, bias value, and R-RMSE value of ERA5 are smaller than the other reanalysis products in each climatic zone, which shows the accuracy of ERA5 is the highest even in different climatic zones. RMSE values of five reanalysis products are less than 2 mm in the polar regions, and those are beyond 3 mm in the tropical regions except for the ERA5 whose RMSE value is 2.45 mm. The results show the absolute accuracies of five reanalysis products are higher in the polar regions and lower in the tropical regions. The number and spatiotemporal resolution of radiosonde stations are limited, the accuracy of water vapor retrieval with respect to remote-sensing satellites is low due to the lower signal-to-noise ratios of the measured spectra and the limitation of weather conditions, and the super refraction affects the quality of GNSS radio occultation data observation in the tropical regions. All these factors result in the lower absolute accuracies of five reanalysis water vapor products in the tropical regions compared with the other regions.
Lastly, the effects of seasonality on the accuracies of the five reanalysis products are analyzed. The results indicate that the RMSE values of ERA5 in the polar regions do not have obvious seasonality. However, RMSE values of five reanalysis products in the temperate regions vary significantly with the seasons, and those in summer and in winter are the largest and the smallest in the temperate regions, respectively. Seasons in the tropical regions which are controlled by precipitation variability are divided into the dry season and wet season. RMSE values with respect to reanalysis products in wet season are slightly greater than those in dry season in some GNSS stations. The results demonstrate that there are certain seasonal variations which are not very significant in terms of RMSE values with respect to reanalysis products in the tropical regions. Our future research will focus on assimilating GNSS PWV into global reanalysis products and improving the accuracy of the reanalysis products. In addition, it is also necessary to analyze the PWV accuracy of the five reanalysis products over the ocean regions.