Improved Modeling of Global Ionospheric Total Electron Content Using Prior Information

The Ionosphere Working Group of the International GNSS Service (IGS) has been a reliable source of global ionospheric maps (GIMs) since 1998. Modeling of the global ionospheric total electron content (TEC) is performed daily by several Ionosphere Associate Analysis Centers (IAACs). Four IAACs (CODE, ESA, CAS and WHU) use the spherical harmonic (SH) expansion as their primary method for modeling GIMs. The IAACs generally solve a normal equation to obtain the SH coefficients and Differential Code Biases (DCBs) of satellites and receivers by traditional least-squares estimation (LSE) without any prior knowledge. In this contribution, an improved method is proposed and developed for global ionospheric modeling based on utilizing prior knowledge. Prior values of SH coefficients and DCBs of satellites and receivers, as well as the variance factor and covariance matrix, could be obtained from the ionospheric modeling on the previous day. The parameters can subsequently be updated through GNSS measurements to achieve higher accuracy. Comparisons are carried out between WHU products based either on priori information or original LSE and IGS final products, other IAAC products, and JASON data for the year 2014. The results indicate that there is improved consistency between WHU GIMs and IGS final GIMs, other IAAC products, and JASON data, particularly in comparison with ESA and UPC products, with the probabilities of achieving better consistency with these products exceeding 95%. Moreover, WHU-produced DCBs of satellites also have slightly improved consistency with IGS final GIMs and IAAC products.

The Ionosphere Working Group of the International GNSS Service (IGS) [27] has been a reliable source of global ionospheric maps (GIMs) since 1998 [28,29].The group consists of several Ionosphere Associate Analysis Centers (IAACs), such as the Center for Orbit Determination in Europe (CODE) [30], the European Space Operations Center (ESOC) of the European Space Agency (ESA) [2], the Jet Propulsion Laboratory (JPL) [31], and the Technical University of Catalonia (UPC) [32].At its workshop in February 2016, the IGS announced that three additional members had joined the Ionosphere Working Group, namely Natural Resources Canada (NRCan), the Chinese Academy of Sciences (CAS), and Wuhan University (WHU).Based on a general assessment of the new GIMs, they appear to provide similar or better precision than the existing IGS GIMs [33].It is anticipated that the final IGS GIMs culminating from more products generated by the additional IAACs would be more robust and have better accuracy.
Methods for modeling the global ionospheric TEC have previously been developed and are maturing, though improved modeling is still an open area for research.The further improvement of global vertical TEC (VTEC) maps remains as an important goal.Several researchers have proposed various approaches to overcome the gaps in VTEC maps on the global scale that occur due to the scarcity of ionospheric pierce points (IPPs) over the oceans and in southern latitudes.The Kriging interpolation technique has previously been used to update VTEC maps and provides better UPC GIMs that exhibit approximately 12% improvement in a self-consistency test [34].Zhang proposed the inequality-constrained least squares (ICLS) method to eliminate non-physical negative values in GIMs [35].In our previous study, ionospheric empirical model IRI 2012 was used to provide VTEC values as virtual measurements for global ionospheric modeling and the release of improved GIMs [10].Additional a priori VTEC values obtained from IRI or other empirical models can also be used for modeling; this has resulted in slightly improved precision [36,37].
Presently, IAACs process GNSS measurements for the release of daily GIMs.The daily data processing is performed separately.However, additional a priori information could actually be available before the next daily modeling of global ionospheric VTEC values, such as the coefficients of the model, DCBs of satellites and receivers, and the standard deviation of and information about the normal equation from the ionospheric modeling on the previous day.In this contribution, the priori information is used to provide the initial values for modeling, particularly the coefficients of the model and the DCBs of satellites.The remaining manuscript is organized as follows.In Section 2, a methodology for global ionospheric modeling based on utilizing priori information is outlined.In Section 3, the improved results of the modeling using this methodology are presented and analyzed.Finally, the study's conclusions are presented in the last section.

Basic Methodology
Since 1998, IAACs have been releasing daily GIMs independently using different approaches.CODE uses a spherical harmonic (SH) expansion referring to a solar geomagnetic reference frame for representing GIMs [30].UPC uses a two-layer tomographic model for the TEC estimation, then presents improved GIMs using a Kriging interpolation technique [32,34].ESA has developed a three-dimensional mathematical ionosphere model based on a simple Chapman profile approach; an approach similar to CODE's SH expansion is then followed.JPL's approach is based on interpolating TEC within triangular tiles; the approach has been extended to include climatological models as simulated data, so that VTEC maps can be generated without gaps [31].WHU also uses SH functions for global ionospheric modeling, along with CODE and ESA.The main equation used for modeling is presented as follows [30,35,38]: where P 1 and P 2 are the smoothed dual-frequency code measurements; f 1 and f 2 are the carrier frequencies of the L1 and L2 signals, respectively; m f is the ionospheric mapping function, which depends on the zenith distance z at the receiver's location; VTEC is the vertical TEC at the IPP; c is the speed of light; DCB s and DCB r are the differential code biases of satellites and receivers, respectively.Code measurements (P 1 and P 2 ) are smoothed by the carrier-phase measurements to obtain high-precision code observables.The code observables are actually replaced by the carrier-phases, shifted by the average value of code minus the phase in a continuous arc [39].We followed the widely used thin-shell approximation of the ionosphere, and used the same MLSM mapping function [30] as CODE to transform Slant TEC to VTEC.SH functions were used for modeling VTEC with reference to a solar geomagnetic reference frame according to the following Equation (2): where ϕ is the geomagnetic latitude of the IPP; λ is the sun-fixed longitude of the IPP; n and m are the degree and order of the model, respectively; P nm is the normalized associated Legendre function of degree n and order m; a nm and b nm are the unknown SH coefficients and GIM parameters, respectively.

Estimation Using Priori Information
Suppose L is a vector of observations.Moreover, suppose the probability density function of vector L is dependent on a vector of unknown parameters X.The random parameter vector X has the probability density function p(X).The density function of L is introduced as the conditional density function p(L|X) , which may be called a likelihood function.The probability density function p(X|L) of parameters X can be described as follows: where ∞ denotes proportionality.The detailed theoretical derivation is available in Koch [40] and Berger, et al. [41].The density function p(X) is known for parameters X before observations L have been taken.The function therefore represents the prior distribution of parameters X.Once data L have been obtained, p(X|L) represents the posterior distribution of parameters X.Thus, the posterior distribution combines the prior information already available with the updated information obtained from the observations.Bayes theorem states that through the likelihood function, the observations modify the prior values of the parameters, thus leading to the posterior density function of the parameters [40].
If the prior information for the parameters vector X is known, the prior mathematic expectation is X, the prior covariance matrix is Σ X , and because X is supposed to obey the normal distribution function N, the prior distribution of parameters X is X ∼ N(X, Σ X ).Moreover, if the observations vector L obeys a normal distribution and the variance factor is σ 2 , the prior distribution of parameters X is a conjugate prior.Hence, the posterior density of X is also normally distributed [40], , where A is the design matrix, P is the weight matrix, and X 0 is the posterior mathematic expectation.Therefore, the least squares Bayesian estimation of X is given as [42]: Furthermore, the variance factor is estimated as: where To consider whether the prior density of X has a normal distribution, a quantile-quantile (QQ) plot is used to investigate the gaussianity of the distribution of prior parameters X. Figure 1 shows the QQ plot, the linear fitting, and a 95% confidence interval.Firstly, the value of the coefficient of determination (R 2 = 0.984) indicated that the standardized SH coefficient data was highly related to the normal theoretical quantiles.Indeed, the QQ plot showed only minimal deviation from a straight-line pattern.Most of the QQ points were within the 95% confidence interval, except at the extreme tails.Therefore, this data provides evidence supporting the assumption that the prior parameters X have a normal distribution (Figure 1).
To consider whether the prior density of X has a normal distribution, a quantile-quantile (QQ) plot is used to investigate the gaussianity of the distribution of prior parameters X. Figure 1 shows the QQ plot, the linear fitting, and a 95% confidence interval.Firstly, the value of the coefficient of determination (R 2 = 0.984) indicated that the standardized SH coefficient data was highly related to the normal theoretical quantiles.Indeed, the QQ plot showed only minimal deviation from a straightline pattern.Most of the QQ points were within the 95% confidence interval, except at the extreme tails.Therefore, this data provides evidence supporting the assumption that the prior parameters X have a normal distribution (Figure 1).The blue and green dashed lines present the lower limit and upper limit of the 95% confidence interval, respectively.
In practical daily modeling of the global ionospheric TEC, the prior values of unknown parameters X (including SH coefficients and DCBs of satellites and receivers) can be ascertained from the ionospheric modeling on the previous day, which involves calculation of the prior mathematic expectation, the prior covariance matrix and the prior variance factor.Then, once the GNSS measurements are available for modeling on the current day, the parameters vector X can be updated by Bayesian estimation with priori information.

Results and Analysis
In order to evaluate the performance of the proposed method, GPS data from approximately 350 IGS stations in 2014 are collected for daily global ionospheric modeling by both original LSE without prior knowledge and Bayesian estimation with prior knowledge of SH coefficients and DCBs of satellites and receivers.Two types of comparison outcomes were investigated: the daily average (bias) and the root mean square (RMS) of the differences between the modeled data and the IGS final GIMs or products from other IAACs, as shown in Equations ( 6) and (7).In these equations, n is the total number of grid points of a daily GIMs product, which is equal to 67,379 (71 × 73 × 13) [30,43], and and are WHU VTEC values and VTEC values from final IGS GIMs or the other IAAC products, respectively.The blue and green dashed lines present the lower limit and upper limit of the 95% confidence interval, respectively.
In practical daily modeling of the global ionospheric TEC, the prior values of unknown parameters X (including SH coefficients and DCBs of satellites and receivers) can be ascertained from the ionospheric modeling on the previous day, which involves calculation of the prior mathematic expectation, the prior covariance matrix and the prior variance factor.Then, once the GNSS measurements are available for modeling on the current day, the parameters vector X can be updated by Bayesian estimation with priori information.

Results and Analysis
In order to evaluate the performance of the proposed method, GPS data from approximately 350 IGS stations in 2014 are collected for daily global ionospheric modeling by both original LSE without prior knowledge and Bayesian estimation with prior knowledge of SH coefficients and DCBs of satellites and receivers.Two types of comparison outcomes were investigated: the daily average (bias) and the root mean square (RMS) of the differences between the modeled data and the IGS final GIMs or products from other IAACs, as shown in Equations ( 6) and (7).In these equations, n is the total number of grid points of a daily GIMs product, which is equal to 67,379 (71 × 73 × 13) [30,43], and VTEC w and VTEC r are WHU VTEC values and VTEC values from final IGS GIMs or the other IAAC products, respectively.
Firstly, the differences among the original four IAAC products are shown in Figures 2-4.The remainder of this section subsequently presents details on the differences in the VTEC maps and DCBs of satellites and receivers analyzed using bias and RMS calculations for the differences between the modeled data and the combined IGS final products or the original four IAAC products.
Figure 2 shows the differences in the VTEC maps among the IAAC products in 2014.Generally, the bias among IAAC products was approximately within ±2 TECU.Most RMS values were within 2~8 TECU, except for days 95 and 96, likely because inadequate data were used by UPC for modeling on those days, as ascertained through the respective UPC GIMs.The typical number of IGS stations used for UPC modeling was approximately 250, in contrast to the 133 and 145 stations used on days 95 and 96, respectively.However, the VTEC maps of ESA and UPC had the strongest consistency, with an annual mean RMS of 3.99 TECU.The VTEC maps of CODE and the other IAACs also have good consistency, with their relatively low RMS values.Next, the differences in satellite DCBs among IAAC products is presented in Figure 3.The biases of the satellite DCB comparisons were approximately within ±0.2 ns.Moreover, annual means of the biases were nearly zero.Additionally, most RMS values were approximately 0.2 ns or less.There was a small number of RMS values beyond 0. DCBs of satellites and receivers analyzed using bias and RMS calculations for the differences between the modeled data and the combined IGS final products or the original four IAAC products.
Figure 2 shows the differences in the VTEC maps among the IAAC products in 2014.Generally, the bias among IAAC products was approximately within ±2 TECU.Most RMS values were within 2~8 TECU, except for days 95 and 96, likely because inadequate data were used by UPC for modeling on those days, as ascertained through the respective UPC GIMs.The typical number of IGS stations used for UPC modeling was approximately 250, in contrast to the 133 and 145 stations used on days 95 and 96, respectively.However, the VTEC maps of ESA and UPC had the strongest consistency, with an annual mean RMS of 3.99 TECU.The VTEC maps of CODE and the other IAACs also have good consistency, with their relatively low RMS values.Next, the differences in satellite DCBs among IAAC products is presented in Figure 3.The biases of the satellite DCB comparisons were approximately within ±0.2 ns.Moreover, annual means of the biases were nearly zero.Additionally, most RMS values were approximately 0.2 ns or less.There was a small number of RMS values beyond 0.2 ns and several RMS values beyond 0.4 ns.The set of RMS values beyond 0.2 ns included the RMS of differences of DCB on pseudorandom noise (PRN) 6 and PRN 10 between the CODE and JPL products.The set of RMS values beyond 0.4 ns included the RMS of differences of DCB on PRN 6, PRN 10, and PRN 30 between the CODE and UPC products.It was apparent that the CODE and ESA satellite DCBs had the strongest consistency in terms of RMS values.Next, Figure 4 showed the differences in receiver DCBs among IAAC products.The bias was approximately zero for the comparisons in 2014, except for those in spring and early summer.The majority of RMS values were close to or less than 1 ns, particularly for the differences between CODE and JPL.There were a few RMS values of approximately 2 ns or more, as shown in Figure 4.Among these larger RMS values, the majority were for differences between CODE and UPC or for those between JPL and UPC.Nevertheless, the low annual mean RMS values as shown in Figure 4 indicated that the receiver DCBs provided by IAACs were not significantly different.

Comparison of VTEC Maps
The IGS Ionosphere Working Group has been providing reliable ionospheric products by combining global ionospheric VTEC maps released by CODE, JPL, ESA and UPC since 1998.Two kinds of VTEC maps from WHU, i.e., those generated either by LSE or Bayesian estimation, are investigated to test the performance of the modeling.Figures 5 and 6 present the differences in the new VTEC maps compared to the IGS final GIMs and other IAAC products, respectively.Also, the annual mean of biases and RMS values obtained through the comparisons are presented in Table 1.

Comparison of VTEC Maps
The IGS Ionosphere Working Group has been providing reliable ionospheric products by combining global ionospheric VTEC maps released by CODE, JPL, ESA and UPC since 1998.Two kinds of VTEC maps from WHU, i.e., those generated either by LSE or Bayesian estimation, are investigated to test the performance of the modeling.Figures 5 and 6         For further analysis, the performance of the new approach with Bayesian estimation as compared to original LSE was investigated and presented in terms of decrement of RMS values in Figures 7 and 8.The RMS decrease value was calculated by the difference of RMS values between ORG and BYS solutions compared with IAAC GIMs.Most RMS decrease values were positive, whereas a few are negative.These results indicated that there was a strong probability of improved consistency between WHU GIMs based on Bayesian estimation and IGS final GIMs, as well as other IAAC products.Although the RMS decrease values for the JPL comparison appeared more negative than those for the IGS, CODE, ESA and UPC comparisons, the values were probably positive in most cases.Table 2 presents the proportion of positive values, which reflects the probability of achieving better consistency, as well as the minimum, maximum and mean of RMS decrease values in units of TECU.From the probabilities depicted in Table 2, it was apparent that for the most part, WHU GIMs based on Bayesian estimation have improved consistency with IGS final GIMs and other IAAC products, particularly ESA and UPC, with probabilities of achieving better consistency exceeding 95%.Although it was difficult to state for certain that a real improvement has occurred based on the comparison of two results that are at the same level of performance, the improved consistency was clear in comparison with both the IGS final GIMs and all other IAAC products.Thus, a real and significant improvement of VTEC maps by modeling based on Bayesian estimation is indicated.For further analysis, the performance of the new approach with Bayesian estimation as compared to original LSE was investigated and presented in terms of decrement of RMS values in Figures 7 and 8.The RMS decrease value was calculated by the difference of RMS values between ORG and BYS solutions compared with IAAC GIMs.Most RMS decrease values were positive, whereas a few are negative.These results indicated that there was a strong probability of improved consistency between WHU GIMs based on Bayesian estimation and IGS final GIMs, as well as other IAAC products.Although the RMS decrease values for the JPL comparison appeared more negative than those for the IGS, CODE, ESA and UPC comparisons, the values were probably positive in most cases.Table 2 presents the proportion of positive values, which reflects the probability of achieving better consistency, as well as the minimum, maximum and mean of RMS decrease values in units of TECU.From the probabilities depicted in Table 2, it was apparent that for the most part, WHU GIMs based on Bayesian estimation have improved consistency with IGS final GIMs and other IAAC products, particularly ESA and UPC, with probabilities of achieving better consistency exceeding 95%.Although it was difficult to state for certain that a real improvement has occurred based on the comparison of two results that are at the same level of performance, the improved consistency was clear in comparison with both the IGS final GIMs and all other IAAC products.Thus, a real and significant improvement of VTEC maps by modeling based on Bayesian estimation is indicated.

Comparison of Satellites DCBs
In addition to the VTEC maps, satellite DCBs were also compared.Figures 9 and 10 present the bias and RMS values, respectively, of the differences of 32 GPS satellite DCBs between WHU products and IGS final products or other IAAC products.Also, the average of the biases and RMS values of the differences in 32 GPS satellite DCBs are presented in Table 3.Most bias values for satellite DCBs were approximately within ±0.2 ns, which shows that WHU products had nearly the same performance as IAAC products.However, the bias values for satellite DCBs in comparisons involving the WHU BYS solution changed little compared to those involving the WHU ORG solution.For the comparison with IGS final products, most RMS values were approximately 0.1 ns, except those for the PRN 1 and PRN 30 satellites, which had RMS values of approximately 0.4 ns and 0.3 ns, respectively.The situation of differences in satellite DCBs for the comparison with CODE products was similar to that of the IGS comparison.Most RMS values were approximately 0.1 ns, with the exception of those for certain satellites, particularly PRN 1, 6 and 30, which had RMS values beyond 0.3 ns.The results were similar for the JPL, ESA and UPC comparisons, as shown in Figure 10.Most RMS values were approximately 0.2 ns, except for a few satellites; these values were similar to RMS values for the differences in satellites DCBs among IAAC products.Furthermore, the RMS values for the differences with most of the satellite DCBs when comparing against both IGS final products and other IAAC products were lower for the new solutions based on Bayesian estimation than those based on original LSE.Therefore, the proposed method for global ionospheric modeling based on Bayesian

Comparison of Satellites DCBs
In addition to the VTEC maps, satellite DCBs were also compared.Figures 9 and 10 present the bias and RMS values, respectively, of the differences of 32 GPS satellite DCBs between WHU products and IGS final products or other IAAC products.Also, the average of the biases and RMS values of the differences in 32 GPS satellite DCBs are presented in Table 3.Most bias values for satellite DCBs were approximately within ±0.2 ns, which shows that WHU products had nearly the same performance as IAAC products.However, the bias values for satellite DCBs in comparisons involving the WHU BYS solution changed little compared to those involving the WHU ORG solution.For the comparison with IGS final products, most RMS values were approximately 0.1 ns, except those for the PRN 1 and PRN 30 satellites, which had RMS values of approximately 0.4 ns and 0.3 ns, respectively.The situation of differences in satellite DCBs for the comparison with CODE products was similar to that of the IGS comparison.Most RMS values were approximately 0.1 ns, with the exception of those for certain satellites, particularly PRN 1, 6 and 30, which had RMS values beyond 0.3 ns.The results were similar for the JPL, ESA and UPC comparisons, as shown in Figure 10.Most RMS values were approximately 0.2 ns, except for a few satellites; these values were similar to RMS values for the differences in satellites DCBs among IAAC products.Furthermore, the RMS values for the differences with most of the satellite DCBs when comparing against both IGS final products and other IAAC products were lower for the new solutions based on Bayesian estimation than those based on original LSE.Therefore, the proposed method for global ionospheric modeling based on Bayesian estimation with prior knowledge could improve not only VTEC maps, but also satellites DCBs.Also, the differences of satellites DCBs between the two approaches were very small.This is because there are so many GNSS measurements (over 3 million observed by around 300 IGS stations per day) for global ionosphere modeling.These rich data are quite enough to estimate 32 GPS satellite DCBs with high accuracy.Thus, the BYS approach brings no obvious improvement.estimation with prior knowledge could improve not only VTEC maps, but also satellites DCBs.Also, the differences of satellites DCBs between the two approaches were very small.This is because there are so many GNSS measurements (over 3 million observed by around 300 IGS stations per day) for global ionosphere modeling.These rich data are quite enough to estimate 32 GPS satellite DCBs with high accuracy.Thus, the BYS approach brings no obvious improvement.

Comparison of Receiver DCBs
Because the prior values of SH coefficients and DCBs of satellites and receivers from the previous day's modeling are used for Bayesian estimation, a comparison of receiver DCBs was also performed.Figures 11 and 12 present the bias and RMS values, respectively, of the differences in receiver DCBs between WHU products and IGS final products or other IAAC products.Also, the annual mean of estimation with prior knowledge could improve not only VTEC maps, but also satellites DCBs.Also, the differences of satellites DCBs between the two approaches were very small.This is because there are so many GNSS measurements (over 3 million observed by around 300 IGS stations per day) for global ionosphere modeling.These rich data are quite enough to estimate 32 GPS satellite DCBs with high accuracy.Thus, the BYS approach brings no obvious improvement.

Comparison of Receiver DCBs
Because the prior values of SH coefficients and DCBs of satellites and receivers from the previous day's modeling are used for Bayesian estimation, a comparison of receiver DCBs was also performed.Figures 11 and 12 present the bias and RMS values, respectively, of the differences in receiver DCBs between WHU products and IGS final products or other IAAC products.Also, the annual mean of

Comparison of Receiver DCBs
Because the prior values of SH coefficients and DCBs of satellites and receivers from the previous day's modeling are used for Bayesian estimation, a comparison of receiver DCBs was also performed.Figures 11 and 12 present the bias and RMS values, respectively, of the differences in receiver DCBs between WHU products and IGS final products or other IAAC products.Also, the annual mean of the biases and RMS values of the differences in receiver DCBs are presented in Table 4.The bias of receiver DCBs was mostly distributed in an interval from −1 ns to 0 ns, with the exception of some bias values smaller than −1 ns seen with the UPC comparison for spring days.Most RMS values were approximately 1 ns, which indicates that WHU products had the same level of performance as other IAAC products.Notably, the values were beyond 1 ns on most of the days in 2014 for the UPC comparison, particularly in spring.Moreover, the annual means of RMS values were slightly larger based on Bayesian estimation than those based on the original LSE.However, the relative value of this increase was fairly small compared to the absolute values of receivers DCB (general mean value beyond 10 ns).This was due to the receiver DCB being considered as a constant parameter during a 1-day period in the ionosphere modeling.The modeling errors of receivers DCB may have been larger than the differences between two approaches.Thus, the BYS approach gives approximately the same performance as the ORG approach.
the biases and RMS values of the differences in receiver DCBs are presented in Table 4.The bias of receiver DCBs was mostly distributed in an interval from −1 ns to 0 ns, with the exception of some bias values smaller than −1 ns seen with the UPC comparison for spring days.Most RMS values were approximately 1 ns, which indicates that WHU products had the same level of performance as other IAAC products.Notably, the values were beyond 1 ns on most of the days in 2014 for the UPC comparison, particularly in spring.Moreover, the annual means of RMS values were slightly larger based on Bayesian estimation than those based on the original LSE.However, the relative value of this increase was fairly small compared to the absolute values of receivers DCB (general mean value beyond 10 ns).This was due to the receiver DCB being considered as a constant parameter during a 1-day period in the ionosphere modeling.The modeling errors of receivers DCB may have been larger than the differences between two approaches.Thus, the BYS approach gives approximately the same performance as the ORG approach.the biases and RMS values of the differences in receiver DCBs are presented in Table 4.The bias of receiver DCBs was mostly distributed in an interval from −1 ns to 0 ns, with the exception of some bias values smaller than −1 ns seen with the UPC comparison for spring days.Most RMS values were approximately 1 ns, which indicates that WHU products had the same level of performance as other IAAC products.Notably, the values were beyond 1 ns on most of the days in 2014 for the UPC comparison, particularly in spring.Moreover, the annual means of RMS values were slightly larger based on Bayesian estimation than those based on the original LSE.However, the relative value of this increase was fairly small compared to the absolute values of receivers DCB (general mean value beyond 10 ns).This was due to the receiver DCB being considered as a constant parameter during a 1-day period in the ionosphere modeling.The modeling errors of receivers DCB may have been larger than the differences between two approaches.Thus, the BYS approach gives approximately the same performance as the ORG approach.Generally, the comparisons indicated that WHU VTEC maps and satellite DCBs calculated based on Bayesian estimation with prior knowledge had improved consistency with IGS final combined products and other IAAC products.Moreover, the receiver DCBs appeared to exhibit approximately the same performance when calculated based on either of the two approaches.According to the results presented in this section, it is apparent that WHU products have the same level of performance as other IAAC products.

Validation with JASON Data
An independent source of VTEC measurements is the data observed by a dual-frequency altimeter instrument on a JASON satellite; this data has previously been used to validate final IGS ionospheric products and IAAC products [28].In this study, JASON VTEC (hereafter J2TEC) measurements were used to validate the VTEC values of GNSS-derived VTEC maps (including IGS GIMs, WHU products based on ORG and BYS solutions, and other IAAC products).Figures 13-16 show bias and RMS values for the differences between GNSS-derived VTEC maps and J2TEC in 2014.
As depicted in Figure 13, the values of bias at low latitudes were larger than those at mid-high latitudes.Other scholars have obtained results similarly showing that J2TEC values are larger than the VTEC values derived from GNSS measurements at mid-high latitudes [28,37,44].In general, the bias among GNSS-derived VTEC maps and J2TEC presented a similar trend, with an average difference of approximately 2 TECU. Figure 14 shows that the RMS values of differences ranged from 4 TECU to more than 12 TECU.The differences among RMS values were smaller for mid-low latitudes and middle latitudes than those for low latitudes, and were especially small for middle latitudes of the Northern Hemisphere.There were obvious differences among the RMS values for low latitudes of the Northern Hemisphere, with a maximum difference of more than 4 TECU.Figures 15 and 16 showed that WHU products perform similarly to IGS GIMs and other IAAC products with respect to J2TEC.Bias for most latitudes decreased with the BYS solution, particularly so for low latitudes.Additionally, as shown in Figure 16, RMS values also decreased for low latitudes and mid-low latitudes of the Southern Hemisphere.These results indicated that there was a greater consistency between WHU products based on Bayesian estimation and JASON VTEC measurements.

Conclusions
Global ionospheric VTEC maps provided by the IGS Ionospheric Working Group are important resources for studies of the ionosphere and for satellite positioning.Different methods and data processing strategies for ionospheric modeling are performed by IAACs independently, which is the main reason of the differences among the IAACs' products.Although the IGS final combined GIMs have been a reliable source of ionospheric products, further improvement of VTEC maps remains an open area for research.IAACs have a responsibility to provide better and more robust ionospheric products.To this end, this contribution proposes an improved method for global ionospheric modeling based on Bayesian estimation with prior knowledge.The estimated SH coefficients and DCBs of satellites and receivers obtained from the previous day's modeling could be considered as prior knowledge, thereby providing approximate initial values for unknown parameters.The initial values can then be updated through real GNSS measurements such that the final derived VTEC maps

Conclusions
Global ionospheric VTEC maps provided by the IGS Ionospheric Working Group are important resources for studies of the ionosphere and for satellite positioning.Different methods and data processing strategies for ionospheric modeling are performed by IAACs independently, which is the main reason of the differences among the IAACs' products.Although the IGS final combined GIMs have been a reliable source of ionospheric products, further improvement of VTEC maps remains an open area for research.IAACs have a responsibility to provide better and more robust ionospheric products.To this end, this contribution proposes an improved method for global ionospheric modeling based on Bayesian estimation with prior knowledge.The estimated SH coefficients and DCBs of satellites and receivers obtained from the previous day's modeling could be considered as prior knowledge, thereby providing approximate initial values for unknown parameters.The initial values can then be updated through real GNSS measurements such that final derived VTEC maps have better consistency with IAACs GIMs than those by using original LSE without prior knowledge.VTEC maps, satellite DCBs and receiver DCBs from WHU products based on either Bayesian estimation or original LSE are compared with those from IGS final products and other IAAC products.The results show that there is strong agreement between WHU GIMs derived by Bayesian estimation and IGS final GIMs or other IAAC products, particularly those of ESA and UPC, with probabilities of improved consistency exceeding 95%.Satellite DCBs derived from the Bayesian estimation also have improved consistency with IGS final products and IAAC products.At the same time, the Bayesian approach gives approximately the same performance of receiver DCBs as that by original LSE approach.Additionally, a validation of VTEC maps against external independent JASON data is performed.The investigation shows that there is improved consistency between WHU products derived by Bayesian estimation and J2TEC measurements relative to those derived by the original LSE solution.These results also indicate that WHU's products have the same level of performance as the other IAAC products.However, this study investigates only a limited dataset covering a one-year timeframe.Moreover, high levels of solar activity and geomagnetic storms would also impact the ionosphere.Therefore, it would be worthwhile for further studies to investigate the performance of the proposed method under different space weather conditions.

Figure 1 .
Figure 1.Quantile-quantile (QQ) plot with black points representing standardized spherical harmonic (SH) coefficients versus a normal distribution.The red line displays the linearity of the QQ points.The blue and green dashed lines present the lower limit and upper limit of the 95% confidence interval, respectively.

Figure 1 .
Figure 1.Quantile-quantile (QQ) plot with black points representing standardized spherical harmonic (SH) coefficients versus a normal distribution.The red line displays the linearity of the QQ points.The blue and green dashed lines present the lower limit and upper limit of the 95% confidence interval, respectively.
2 ns and several RMS values beyond 0.4 ns.The set of RMS values beyond 0.2 ns included the RMS of differences of DCB on pseudorandom noise (PRN) 6 and PRN 10 between the CODE and JPL products.The set of RMS values beyond 0.4 ns included the RMS of differences of DCB on PRN 6, PRN 10, and PRN 30 between the CODE and UPC products.It was apparent that the CODE and ESA satellite DCBs had the strongest consistency in terms of RMS values.Next, Figure 4 showed the differences in receiver DCBs among IAAC products.The bias was approximately zero for the comparisons in 2014, except for those in spring and early summer.The majority of RMS values were close to or less than 1 ns, particularly for the differences between CODE and JPL.There were a few RMS values of approximately 2 ns or more, as shown in Figure 4.Among these larger RMS values, the majority were for differences between CODE and UPC or for those between JPL and UPC.Nevertheless, the low annual mean RMS values as shown in Figure 4 indicated that the receiver DCBs provided by IAACs were not significantly different.Remote Sens. 2018, 10, 63 5 of 19

Figure 2 .Figure 2 .Figure 3 .
Figure 2. Differences in VTEC maps among IAAC products in 2014.(a) Biases of global vertical total electron content (VTEC) maps comparisons among Ionosphere Associate Analysis Center (IAAC) products.(b) Root Mean Square (RMS) of differences in VTEC maps among IAAC products.(b) Figure 2. Differences in VTEC maps among IAAC products in 2014.(a) Biases of global vertical total electron content (VTEC) maps comparisons among Ionosphere Associate Analysis Center (IAAC) products.Root Mean Square (RMS) of differences in VTEC maps among IAAC products.

Figure 3 .
Figure 3. Differences in satellite Differential Code Biases (DCBs) among IAAC products in 2014.(a) Biases of satellites DCB comparisons among IAAC products.(b) RMS of differences in satellites DCB among IAAC products.

Figure 4 .
Figure 4. Differences in receiver DCBs among IAAC products in 2014.(a) Biases of receivers DCB comparisons among IAAC products.(b) RMS of differences in receivers DCB among IAAC products.
Each figure includes the biases and RMS values of the differences between the VTEC maps from WHU based on original LSE (ORG) or Bayesian estimation (BYS) and an IAAC product or the IGS final GIMs.The figures show that there is a minor increase in the absolute values of biases for the Bayesian estimation compared to those for the original LSE solution.At the same time, most RMS values decrease slightly in the comparisons made with the Bayesian estimation relative to those made with the original LSE solution.According to the annual means of the RMS values, the WHU GIMs based on Bayesian estimation have better consistency with the IGS final GIMs as well as the other four IAAC products.

Figure 4 .
Figure 4. Differences in receiver DCBs among IAAC products in 2014.(a) Biases of receivers DCB comparisons among IAAC products.(b) RMS of differences in receivers DCB among IAAC products.
present the differences in the new VTEC maps compared to the IGS final GIMs and other IAAC products, respectively.Also, the annual mean of biases and RMS values obtained through the comparisons are presented in Table1.Each figure includes the biases and RMS values of the differences between the VTEC maps from WHU based on original LSE (ORG) or Bayesian estimation (BYS) and an IAAC product or the IGS final GIMs.The figures show that there is a minor increase in the absolute values of biases for the Bayesian estimation compared to those for the original LSE solution.At the same time, most RMS values decrease slightly in the comparisons made with the Bayesian estimation relative to those made with the original LSE solution.According to the annual means of the RMS values, the WHU GIMs based on Bayesian estimation have better consistency with the IGS final GIMs as well as the other four IAAC products.

Figure 5 .Figure 6 .
Figure 5. Differences between VTEC maps and The Ionosphere Working Group of the International GNSS Service (IGS) final global ionospheric maps (GIMs).The red and blue points show the RMSs of differences in VTEC maps based on original least-squares estimation (LSE) or Bayesian estimation, respectively, compared with IGS final GIMs.

Figure 5 .
Figure 5. Differences between VTEC maps and The Ionosphere Working Group of the International GNSS Service (IGS) final global ionospheric maps (GIMs).The red and blue points show the RMSs of differences in VTEC maps based on original least-squares estimation (LSE) or Bayesian estimation, respectively, compared with IGS final GIMs.

Figure 5 .Figure 6 .
Figure 5. Differences between VTEC maps and The Ionosphere Working Group of the International GNSS Service (IGS) final global ionospheric maps (GIMs).The red and blue points show the RMSs of differences in VTEC maps based on original least-squares estimation (LSE) or Bayesian estimation, respectively, compared with IGS final GIMs.

Figure 6 .
Figure 6.Differences between VTEC maps and IAAC GIMs.The red points and blue points show the RMSs of differences in VTEC maps based on original LSE and Bayesian, respectively, compared with IAAC GIMs.(a) Compared with Center for Orbit Determination in Europe (CODE).(b) Compared with Jet Propulsion Laboratory (JPL).(c) Compared with European Space Agency (ESA).(d) Compared with Technical University of Catalonia (UPC).

Figure 7 .
Figure 7. RMS decrease of differences in VTEC maps between Wuhan University (WHU) GIMs based on Bayesian estimation and the IGS final GIMs compared to those between GIMs based on original LSE and the IGS final GIMs.

Figure 7 .
Figure 7. RMS decrease of differences in VTEC maps between Wuhan University (WHU) GIMs based on Bayesian estimation and the IGS final GIMs compared to those between GIMs based on original LSE and the IGS final GIMs.

Figure 8 .
Figure 8. RMS decrease of differences in VTEC maps between WHU GIMs based on Bayesian estimation and IAAC GIMs compared to those between GIMs based on original LSE and the IGS final GIMs.(a) Compared with CODE.(b) Compared with JPL.(c) Compared with ESA.(d) Compared with UPC

Figure 8 .
Figure 8. RMS decrease of differences in VTEC maps between WHU GIMs based on Bayesian estimation and IAAC GIMs compared to those between GIMs based on original LSE and the IGS final GIMs.(a) Compared with CODE.(b) Compared with JPL.(c) Compared with ESA.(d) Compared with UPC.

Figure 9 .
Figure 9. Biases of satellites DCBs from the comparisons between WHU products based on original LSE or Bayesian estimation and the IGS final combined products or IAAC products.

Figure 10 .
Figure 10.RMS of the differences in satellite DCBs between WHU products based on original LSE or Bayesian estimation, and the IGS final combined products or IAAC products.

Figure 9 .
Figure 9. Biases of satellites DCBs from the comparisons between WHU products based on original LSE or Bayesian estimation and the IGS final combined products or IAAC products.

Figure 9 .
Figure 9. Biases of satellites DCBs from the comparisons between WHU products based on original LSE or Bayesian estimation and the IGS final combined products or IAAC products.

Figure 10 .
Figure 10.RMS of the differences in satellite DCBs between WHU products based on original LSE or Bayesian estimation, and the IGS final combined products or IAAC products.

Figure 10 .
Figure 10.RMS of the differences in satellite DCBs between WHU products based on original LSE or Bayesian estimation, and the IGS final combined products or IAAC products.

Figure 11 .Figure 11 .
Figure 11.Differences in receiver DCBs between WHU products based on original LSE or Bayesian estimation and the IGS final combined products.

Figure 11 .Figure 12 .
Figure 11.Differences in receiver DCBs between WHU products based on original LSE or Bayesian estimation and the IGS final combined products.

Figure 13 .
Figure 13.Bias values for comparison between IGS GIMs or IAAC products and J2TEC in 2014.

Figure 13 .
Figure 13.Bias values for comparison between IGS GIMs or IAAC products and J2TEC in 2014.

Figure 14 .
Figure 14.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.

Figure 15 .
Figure 15.Bias values for comparison between WHU products and J2TEC in 2014.

Figure 14 .
Figure 14.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.

Figure 13 .
Figure 13.Bias values for comparison between IGS GIMs or IAAC products and J2TEC in 2014.

Figure 14 .
Figure 14.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.

Figure 15 .
Figure 15.Bias values for comparison between WHU products and J2TEC in 2014.Figure 15.Bias values for comparison between WHU products and J2TEC in 2014.

Figure 16 .
Figure 16.RMS values of differences between WHU products and J2TEC in 2014.

Figure 16 .
Figure 16.RMS values of differences between WHU products and J2TEC in 2014.

Figure A1 .
Figure A1.Bias values for comparison between IGS GIMs or IAAC products and J2TEC in 2014.

Figure A2 .
Figure A2.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.Figure A2.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.

Figure A2 .
Figure A2.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.Figure A2.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.

Figure A1 .
Figure A1.Bias values for comparison between IGS GIMs or IAAC products and J2TEC in 2014.

Figure A2 .
Figure A2.RMS values of differences between IGS GIMs or IAAC products and J2TEC in 2014.

Figure A4 .
Figure A4.Differences in satellite DCBs between ORG and BYS approaches.Figure A4.Differences in satellite DCBs between ORG and BYS approaches.

Figure A4 .
Figure A4.Differences in satellite DCBs between ORG and BYS approaches.Figure A4.Differences in satellite DCBs between ORG and BYS approaches.

Figure A4 .
Figure A4.Differences in satellite DCBs between ORG and BYS approaches.

Figure A5 .
Figure A5.Differences in receiver DCBs between ORG and BYS approaches.
It should be noted that the denominator in Equation (5) is n − t (where t is the number of unknown parameters) if the equation is solved by LSE.Conversely, if it is solved by Bayesian estimation with prior knowledge, which involves t virtual measurements of parameters in addition to the number of observations L, the degrees of freedom would be (n + t) − t = n.
and n is the number of degrees of freedom, which is equal to total number of the observations L in the Bayesian estimation.

Table 1 .
Differences in VTEC maps based on original LSE and Bayesian approaches compared with IAAC GIMs (unit: TECU).

Table 1 .
Differences in VTEC maps based on original LSE and Bayesian approaches compared with IAAC GIMs (unit: TECU).

Table 1 .
Differences in VTEC maps based on original LSE and Bayesian approaches compared with IAAC GIMs (unit: TECU).

Table 3 .
Differences in satellites DCBs based on original LSE and Bayesian approaches compared with IAAC products (unit: TECU).

Table 3 .
Differences in satellites DCBs based on original LSE and Bayesian approaches compared with IAAC products (unit: TECU).

Table 3 .
Differences in satellites DCBs based on original LSE and Bayesian approaches compared with IAAC products (unit: TECU).

Table 4 .
Differences in receiver DCBs based on original LSE and Bayesian approaches compared with IAAC products (unit: TECU).