Modelling of weighted-mean temperature using regional radiosonde observations in Hunan China

Precipitable water vapour (PWV) over a ground station can be estimated from the global navigation satellite systems (GNSS) signal’s zenith wet delays (ZWD) multiplying by a conversion factor that is a function of weighted-mean temperature (Tm). The commonly used Bevis Tm model (BTM) may not perform well in some regions due to its use of data from North America in the model development. In this study, radiosonde observations in 2012 from three stations Changsha, Huaihua, Chenzhou in Hunan province, China were used to establish a new regional Tm model (RTM) based on a numerical integration and the least squares estimation methods. Four seasonal RTMs were also established and assessed for 2012. The RTM-derived Tm at the three stations from 2012 2014 were validated by comparing it with radiosondederived Tm. Results showed that the accuracy of the yearly RTM was improved by 29% over the BTM, and the bias and root mean square (RMS) of all the four seasonal RTMs were slightly smaller than the yearly RTM, and the accuracy of spring, summer, autumn and winter Tm models is improved by 5, 13, 4, and 5% respectively. In addition, the bias and RMS of the differences between the GNSS-PWV resulting from the RTM-derived Tm and the radiosonde-PWV were 1.13 and 3.21 mm respectively, which are reduced by 34 and 10% respectively. Thus the seasonal RTMs are recommended for GNSS meteorology for Hunan Province. Article history: Received 18 December 2016 Revised 18 April 2017 Accepted 26 May 2017


INtrODUCtION
Water vapour is one of the main driving forces to weather changes and atmospheric circulation. The trend and dynamics of variation in water vapour have strong influence on climate and weather forecasting (Li et al. 2012;Guerova et al. 2016). However, due to its fast formation and development, it is not easy to trace the variation trend and distribution of water vapour on small and medium scale weather systems timely and accurately (Li et al. 2007;Yu and Liu 2009;Wang et al. 2011), especially for the detection and prediction of rainstorms (Bennitt and Jupp 2012;De Haan 2013;Shi et al. 2015). The challenge is that atmospheric variables observed from traditional meteorological sensors most have low temporal and spatial resolutions, e.g., radiosonde balloons are usually launched twice a day, thus tracing water vapour at a high resolution based on the traditional meteorological data is difficult (Bevis et al. 1994;Guerova et al. 2016).
Nowadays, Global Navigation Satellite Systems (GNSS) have heralded a new era for robust retrievals of atmospheric precipitable water vapour (PWV) due to its 24-hr availability, global coverage, high accuracy, high resolution, and low cost (Bevis et al. 1994;Fan et al. 2016;Yeh et al. 2016). Although the zenith tropospheric delay (ZTD) estimated from the GNSS measurements can be directly assimilated into numerical weather prediction (NWP) (Poli et al. 2007;Bennitt and Jupp 2012;De Haan 2013), GNSSderived PWV has the potential to be used for the studies of severe weather (Iwabuchi et al. 2006;Yeung et al. 2009;Li et al. 2012) and climate changes (Kruczyk 2015;Bianchi et al. 2016;Yeh et al. 2016). Previous studies (Song and Terr. Atmos. Ocean. Sci., Vol. 29, No. 2, 187-199, April 2018 Grejner-Brzezinska 2009; Li et al. 2012;Wang et al. 2015) have shown that most severe rainfall events occur in the descending trends of time series PWV over a station after a long ascending period. Moreover, it is suggested by Benevides et al. (2015) that the reliability and accuracy of severe weather event forecast could be improved by the analysis of 2D or 3D variation of PWV fields with the use of other meteorological data (Yeung et al. 2009;Shi et al. 2015;Jiang et al. 2016). GNSS-derived PWV over a station can be obtained from the zenith wet delay (ZWD) multiplying by a conversion factor Π, which is a function of the weightedmean temperature (T m ) over the station (Bevis et al. 1994;Wang et al. 2016). Thus, the accuracy of T m affects the accuracy of the resultant PWV.
The most accurate method for calculating T m is to obtain both temperature and humidity profiles from the radiosonde data (Liu et al. 2001;Wang et al. 2016). In practice, for most GNSS stations there are no co-located radiosonde data available due to the high cost of radiosonde stations. The most common method is to use the relationship between surface temperature T s and T m (Yao et al. 2014(Yao et al. , 2015, e.g., the Bevis model (T m = 0.72T s + 70.2) established in 1992. The Bevis T m model (BTM) was derived from the profiles of vapour pressure and dewpoint temperature from 8718 observations over 13 mid-latitude radiosonde stations (27 -65°N) in North America over a 2-year period. The BTM can be applied for real-time applications since it does not need to use meteorological measurements. However, due to the rapid spatial and temporal variation of atmospheric wet profiles, it is found that the BTM performs unevenly globally, such as in China, its systematic bias is generally above 4 K, with an extreme value of 8 K in some areas (Lan et al. 2016). The error in the BTM can result in a several millimetres error in PWV under an extreme weather condition with large amounts of water vapour . The new numerical relationships between T m and T s have been studied by many researchers using the linear regression method and local data for regional usage in various countries (Liu et al. 2001;Song and Boutiouta 2012;Jiang et al. 2014;Kruczyk 2015), such as the Eastern China region, Hong Kong, Beijing, Wuhan, Chengdu, Zhengzhou, and Jiangxi (Li and Mao 1998;Wang et al. 2007;Lü et al. 2008;Gong 2013).
Since 2011, the Hunan Continuously Operating Reference Stations (HNCORS) network has been established by the Meteorological Bureau and National Land Agency of Hunan Province, China. It consists of 93 stations covering the whole Hunan province (Hu and Tang 2009), and meteorological applications of the ground-based GNSS infrastructure have been put on the business agenda. There are three radiosonde stations -Changsha, Huaihua, and Chenzhou are available in the whole Hunan region, a new T m model based on these local radiosonde data may perform better.
The outline of this paper is as follows. The methodology of computing T m and PWV from radiosonde and ground-based GNSS-ZWD will be given in the second section, followed by one year radiosonde data in 2012 collected at the above three stations were used to establish a regional T m model (RTM) after a brief introduction to the dataset and methodologies. Then the accuracy of the RTM is evaluated using the following two years' true T m derived from radiosonde data by comparing three sets of T m derived from radiosondes, the BTM and the new RTM. Furthermore, the PWVs resulting from the Chenzhou radiosonde data are also compared for investigation of the effect of T m on PWV. Conclusions are given in the last section.

Calculating t m
The Constant, Bevis formula, numerical integration and approximate integration methods are mainly used to calculate T m . The data required, realisation difficulty and the accuracy level of these four methods are shown in Table 1 (Askne and Nordius 1987;Duan et al. 1996;Liu et al. 2001).
Among these methods, the Constant method is the simplest but has the lowest accuracy; the Bevis formula is most widely used currently and can be used for real-time applications, but its accuracy may vary with location and season; the approximate integration needs temperature lapse rate α and vapour pressure decline rate λ, which cannot be easily obtained; the numerical integration is the most accurate method and is also easily implemented (Liu et al. 2001). In E temperature lapse rate α; vapour pressure decline rate λ; surface temperature T s ; gas constant R; gravitational acceleration g hardest low Table 1. Comparison of four typical T m calculation methods (Askne and Nordius 1987;Duan et al. 1996;Liu et al. 2001).
this study, the numerical integration is adopted and its mathematical expression can be formulated as follows (the same as the one shown in Table 1) (Askne and Nordius 1987;Liu et al. 2001).
In fact, vapour pressure is a measure of the amount of moisture in air. When air reaches saturation vapour pressure, water vapour in air will condense. In this case, dew point temperature is the same as air temperature. Vapour pressure under the condition of saturation is the saturated vapour pressure (e), which is a function of dewpoint temperature (T d ) (Gong 2013;Bianchi et al. 2016 Where T d is in Celsius. Equation (4) was given by the World Meteorological Organization in 2008 (Jarraud 2008).

Obtaining PWV
Total precipitable water is the amount of water contained in a vertical column of unit cross section extending from the earth's surface to the "top" of the atmosphere. PWV can be defined as the amount of liquid water if the entire atmospheric vapour within the vertical column were compressed to the point of condensation (Bevis et al. 1992). PWV is measured in millimetres (mm) and four methods can be used to obtain PWV (Chen and Liu 2014). In this study two methods are adopted and they will be elaborated below.

From radiosondes
where g is the gravitational constant (cm s -2 ), P Z and P 0 are the atmospheric pressures at the surface and top of the atmosphere (hPa) respectively, and q, the specific humidity (g/kg), is the mass of dry air (the atmosphere minus water vapour) divided by the total mass of the atmosphere, and is usually expressed as where e (vapour pressure) is calculated by Where T d is dewpoint temperature ( The reason for i being restricted to 16 here is that no water vapour exists at altitudes higher than the 16 km from surface of earth (Jarraud 2008).

From Conversion of GNss-ZWD
Generally, the ZTD of GNSS signals can be estimated using un-differenced precise point positioning (PPP) or differential strategies. In this study, the un-differenced PPP strategy is adopted using the precise satellite orbits (SP3), satellite clock biases and ERPs from the IGS Final/Rapid products with 5-min intervals Kouba and Heroux 2001;Li et al. 2012). The daily GNSS measurements sampled with 30 seconds are used with an elevation cut-off angle of 7°. Epoch-by-epoch Kalman filter was used to estimate the receivers' ZTDs along with their coordinates, clock biases, tropospheric gradients and float ionosphere-free ambiguities. IERS standards are used to the models' corrections of relativistic effects, the phase centre biases and variations of the satellite and the receiver, solid earth and ocean tide loading, phase wide-up, polar motion and nutation McCarthy and Petit 2004). A-priori modelling of the tropospheric delay is done using the global mapping function (GMF) (Boehm et al. 2006).
The ZTD can be divided into two parts: the ZWD and zenith hydrostatical delay (ZHD), and 90% of the ZTD is induced by dry air in the atmosphere. The ZHD has a positive correlation with ground pressure. The ZWD is mainly caused by the atmospheric water vapour which varies rapidly in both spatial and temporal domains, so it can only be estimated by an empirical model at an accuracy of 10 -20%.
The accuracy of the ZHD calculated using the most commonly used Saastamoinen model can be at an accuracy level of millimetres. The Saastamoinen model can be expressed as (Saastamoinen 1972) . .
where p c , c { , and h c are the surface pressure (hPa), geographic latitude and altitude of the station, respectively (km).
The accuracy of surface pressure measurements is generally 0.2 -0.5 hPa, and the accuracy of the ZHD calculated by the Saastamoinen model can be at millimetres. The accuracy of GNSS-PPP-derived ZTD can be also at a level of millimetres. As a result, the accuracy of the ZWD, calculated by ZWD = ZTD -ZHD, is also at the level of millimetres (Bevis et al. 1992;Rocken et al. 1995;Li et al. 2015). The ZWD can be converted to PWV by the following formula, where П is the conversion factor; w t is the water density; R = 461 (J kg -1 K -1 ); k = (3.776 ± 0.014) × 10 5 (K 2 hPa -1 ); and kl = 16.48 (K hPa -1 ). From Eq. (11), we can see that the accuracy of PWV is determined by the accuracy of ZWD and T m .

Data source
In this study, radiosonde data collected from balloonborne instrument platforms with radio transmitting twice a day from the aforementioned three radiosonde stations in Hunan province, Changsha, Huaihua, and Chenzhou (Fig. 1), were used to calculate the time series of T m over the three stations using Eq. (1) during the period of 2012 -2014. Since there are no co-located GNSS data both at Changsha and Huaihua radiosonde stations, PWV from the Chenzhou radiosonde station in 2015 was also calculated using Eq. (9). The radiosonde-derived T m from the three stations during the period of 2012 -2014 and radiosonde-derived PWV from the Chenzhou station in 2015 were used as the ground truth in the performance assessment of BTM and RTM.
The time series of T m (two T m observing values each day) and its corresponding surface temperature T s from the three stations in 2012 were used to establish a new regional T m -T s model or relationship, and the T m resulting from the radiosonde in the period of 2012 -2014 were used to assess the performance of the new RTM. Moreover, PWV over the co-located Chenzhou GNSS station in 2015 were calculated using T m derived from both the Bevis and the RTM. These two sets of GNSS-derived PWVs were then compared with the radiosonde-derived PWV (as the truth) at the Chenzhou station to assess the performance of these two T m models and the impact of the errors in the two models on the resultant PWV.
The performance was measured by the bias and root mean square (RMS) of the time series as defined below where X i is the model-derived (or predicted) T m value from the RTM, and X true is the radiosonde-derived T m . Fig. 1. Location of the three radiosonde stations in Hunan Province, China.

bevis t m Model
In order to analyse the systematic deviation of T m and the accuracy of the BTM in Hunan province, the time series T m (again, two results each day) and the difference between BTM-derived (using the T s of the three years) and radiosonde-derived T m during the period of 2012 -2014 at the aforementioned three GNSS and radiosonde co-located stations are shown in Fig. 2.
As shown in the left-side sub-figure in Fig. 2, a seasonal trend in both the T m time series and the difference can be found. Both the mean value and the amplitude of the T m at the three stations are different. The right-side sub-figures show a certain positive systemic bias in the Bevis-derived T m from the truth. From the statistical result of the above time series, i.e., bias and RMS listed in Table 2, we can see that the biases of the three stations' results are 2.58, 1.81, and 2.10 K; and their RMSs are 3.01, 2.95, and 3.59 K respectively. The RMS indicates the accuracy of the model. The mean bias and mean RMS of all the three stations is 2.16 and 3.18 K respectively. The positive bias value means that the BTM overestimates the T m in Hunan region.

Modelling of regional t m and Its Accuracy
The two-dimensional linear fitting method adopted to model the regional T m has the same expression as the BTM, i.e., T m = a × T s + b, and its observation equation matrix is expressed as where V is the residual vector. In this study, radiosonde-derived T m and T s from all the aforementioned three stations obtained from different amounts of data (4, 8, 12, and 24 months) during the period of 2012 -2013 are substituted into Eq. (13). Then the coefficients a and b are estimated using the least squares estimation method. The RTM results are listed in Table 3. The comparison of the bias and RMS of the time series T m derived from these four new RTM models during the period of 2012 -2014 against the radiosonde-derived T m at the three stations are also listed in Table 3. We can see that the RMSs resulting from eight months or more are almost at the same level (~2.17 mm), while the bias is reduced from 0.36 to 0.17 mm when 12 months or more data are used. These imply that 12 months data are sufficient for establishing a RTM in the Hunan region at a desired accuracy, so the radiosonde-derived T m and T s from all the three aforementioned stations during year 2012 will be used to establish the RTM as expressed by . .
The comparison of the bias and RMS of the time series T m derived from this new model at the three stations during the period of 2012 -2014 against the radiosonde-derived T m at the same stations and same period are listed in Table 4. The accuracy of the RTM is better than that of the BTM in terms of both mean bias and mean RMS of the three stations, 0.18 and 2.17 K against 2.16 and 3.18 K, respectively.
In fact, many RTMs have been developed using the same expression as BTM, the coefficients a and b are rounded or truncated to a two-decimal level. The four-decimal coefficients of Eq. (14) can be truncated to two decimals as following .
. Table 5 indicates the accuracy difference resulting from the two-decimal and four-decimal RTMs in the bias and RMS of the two models at the three stations in the period of 2012 -2014. The four-decimal model result is better than the twodecimal model's, especially its bias is only 0.18 K, which is significantly smaller than the two-decimal's -1.36 K. The bias and RMS of four-decimal T m model are improved by 87 and 14%, respectively over the two decimals. This suggests that the truncation error cannot be neglected, so our following calculations and discussions are all based on the fourdecimal model.

Performance Comparison of the two Models
In Fig. 3a, the red points are the new RTM-derived T m during the period of 2012 -2014 and the black line is the T m fitting result expressed by Eq. (14), the correlation coefficient between T m and T s reaches 0.9425, which means a strong positive correlation between T m and T s . The distribution of the differences between the RTM-derived T m and their truth (from radiosondes) are most in the range of about -2 -2 K (Fig. 3b). Figure 4 is for comparisons between the time series of T m obtained from the new RTM and the BTM and Tables 2 -4 shows the difference between the two models' statistical results. It can be seen that the RMS of the RTM (2.17 K) is smaller than that of the BTM (3.18 K). As for the bias, the RTM's (0.18 K) is much smaller than the BTM's (2.16 K). These two sets results indicate improvements of 91 and 32% in bias and RMS, respectively.

Comparison of GNss-PWV resulting from two t m Models and PWV from radiosonde
The Chenzhou CORS station (named CZSQ) has a colocated radiosonde station with a horizontal distance of about   10 m. The GTS1 digital radiosonde produced by Shanghai Changwang Meteotech Co. Ltd. is used at the Chenzhou station. Generally, the GTS1 digital radiosonde can provide highly reproducible temperature measurements in the troposphere, but it has a warm bias during the daytime (< 0.5 K) and a relatively large cool bias at nighttime (-0.2 to -1.6 K) in the stratosphere. It can also measure atmospheric pressure with an 1 -2 hPa accuracy from ground up to 30 km altitude (Bian et al. 2010).
The GNSS-ZTD of CZSQ is calculated using the PPP strategy as mentioned in section 2 (Li et al. 2012(Li et al. , 2015. The ZHD is determined using Eq. (10) with pressure values from a pressure sensor mounted at this station. The two sets of PWVs converted from the GNSS-ZTDs together with T m from both the BTM and RTM are compared against the true PWV (derived from the GTS1 digital radiosonde data) in 2015 -2016, as shown in Fig. 5a. Figure 5b indicates the accuracy of the two sets of PWVs. However, some GNSS data are missing during the period of 2015 -2016 due to the instrument or network failures, thus the PWVs during the period of data missing cannot be obtained. Figure 6 shows the correlation between the two sets of PWVs and their true values, and Table 6 is their statistical result. The comparison between the two models' results shown in the table, the bias and RMS of the RTM-derived PWV are improved by 34 and 10%, respectively. Figure 2 in section 3.2 shows a seasonal trend in the T m time series. In this section, the correlation of T m and T s across different seasons is investigated. Due to the RMS results are almost at the same level (~2.17 k) when the data are from 8 months or longer, four seasonal RTMs are established using the radiosonde-derived T m and T s from all the aforementioned three stations during the period of 2012 -2014 (9 months data for each season) in this study. Figure 7 shows the seasonal model results and the coeffi-cients of correlation between T m and T s in each season. The correlation coefficients of the spring, summer and autumn seasons all are above 0.86, but the winter season has a much lower value 0.675, which may be due to low temperature and high humidity.

Accuracy of seasonal t m Models
The time series of T m from each of the four seasonal models and the yearly RTM obtained from radiosonde data in 2012 -2014 are compared, and the difference between these models-derived T m and the reference value (truth) of T m (derived from radiosondes) is a measure of their accuracy, see the right column in Fig. 8. It is shown that all the seasonal RTMs are slightly better than the yearly RTM.
The statistics of the results in Fig. 8 are listed in Table 7. We can see that both bias and RMS of all the seasonal RTMs are smaller than that of the yearly RTM and their improvements are 5, 13, 4, and 5%.

CONCLUsION
In this study, a new RTM is developed and its performance is tested based on both numerical integration and least squares estimation methods, and radiosonde observations collected at three stations in Hunan province, Changsha, Huaihua, and Chenzhou in 2012. Results show that the performance of the RTM (in terms of bias and RMS) is better than the BTM (bias: 0.18 vs 2.16 K; RMS: 2.17 vs 3.18 K; improvements are of 91 and 32%). This indicates the new RTM, which is based on local radiosonde data, outperforms the BTM in the same region. With respective to the reference of PWV (derived from radiosonde data), the bias and RMS of GNSS-PWV resulting from the new RTM are 1.03 and 3.21 mm respectively, which leads to 34 and 10% improvements over the BTM results. Four seasonal RTMs are also investigated using each of the seasonal data in 2012 -2014 and their performance is slightly better than the yearly RTM, with the improvements of 5, 13, 4, and 5%, in terms of RMS.   Table 6. Bias, RMS, and correlation coefficient of the PWV resulting from the Bevis and RTM (relative to the truth -radiosonde-derived PWV) at the Chenzhou station (mm).  Therefore, the seasonal RTMs are recommended for studies in the area of GNSS meteorology in the Hunan region. Generally, the bias of the BTM relative to radiosonde-derived T m is above 4 K, and its maximum value may reach 8 K or even more. Theoretically, the relative error of PWV resulting from the error in T m can be approximated by the relative error of T m . Therefore, our newly established RTM can improve the accuracy of PWV by more than 3 mm under an extreme weather condition caused by large amounts of water vapour, which is very common in summer of the Hunan region, where PWV may reach as large as 70 mm.