A mean-squared-error condition for weighting ionospheric delays in GNSS baselines

Abstract

Although ionosphere-weighted GNSS parameter estimation is a popular technique for strengthening estimator performance in the presence of ionospheric delays, no provable rules yet exist that specify the needed weighting in dependence on ionospheric circumstances. The goal of the present contribution is therefore to develop and present the ionospheric conditions that need to be satisfied in order for the ionosphere-weighted solution to be mean squared error (MSE) superior to the ionosphere-float solution. When satisfied, the presented conditions guarantee from an MSE performance view, when (a) the ionosphere-fixed solution can be used, (b) the ionosphere-float solution must be used, or (c) an ionosphere-weighted solution can be used.

1 Introduction

Ionosphere-weighted GNSS parameter estimation is a popular and very flexible technique for strengthening estimator performance in the presence of ionospheric model delays. It is used in a wide variety of GNSS applications, ranging from the use of external ionospheric models (Schaer 1999; Memarzadeh 2009; Feltens et al. 2011) to the incorporation of ionospheric corrections from reference networks, for instance for PPP or PPP-RTK (Odijk 2000; Paziewski 2016; Tomaszewski et al. 2020; Teunissen 2021), and the strengthening of medium to long baseline models (Teunissen 1997; Bock 1998; Odolinski and Teunissen 2017; Brack et al. 2021). However, no clear and provable rules have been established that specify the needed weighting in dependence on ionospheric circumstances. Although for some applications, heuristic rules are available that describe the weighting in dependence on baseline length, the circumstances under which they guarantee an improved solution have not been established.

Figure 1 shows an ionosphere-float and an ionosphere-fixed baseline-scatterplot of the same data from two receivers 176 km apart. Normally one would not think of hard-constraining the between-receiver differential ionospheric delays to zero for baselines of this length. However, as the results of Fig. 1 surprisingly show, in this case the ionosphere-fixed scatterplot is clearly to be preferred over its ionosphere-float counterpart. Although the ionosphere-fixed solution is biased, its mean squared error (MSE) is significantly smaller than that of the ionosphere-float solution. This example thus shows that we need means to be able to recognize and predict the occurrence of such circumstances. The goal of the present contribution is therefore to develop and present the ionospheric condition that needs to be satisfied in order for the ionosphere-weighted solution to be MSE-superior to the ionosphere-float solution. At present no such objective condition exists.

Fig. 1

Ionosphere-float versus ionosphere-fixed single-epoch code-only position scatterplots of baseline CUT0–KELN ($\sim$ 176 km) on 5 January 2021 (12:00–15:00 Perth time, GPS L1/L2, 1 Hz). Left: Ionosphere float (the DD ionospheric delays are assumed completely unknown); Right: Ionosphere fixed (the DD ionospheric delays are hard-constrained to zero)

In order to show how the ionospheric circumstances affect the mean-squared-error performance of the ionosphere-weighted estimators, we will use the following simple metric to capture the variability of the ionosphere,

$$\begin{aligned} {\mathbf {\mathsf{{d}}}}^{2} = \sum _{s=1}^{m} (\texttt {i}_{t}^{s}-\bar{\texttt {i}}_{t})^{2}/(m-1) \end{aligned}$$

(1)

where $\texttt {i}_{t}^{s}$ is the between-receiver single-difference ionospheric delay of satellite s at epoch t and $\bar{\texttt {i}}_{t}$ the average of these delays over all m satellites. In this contribution we develop the conditions that (1) has to satisfy in order to guarantee that from a mean-squared-error performance view, (a) the ionosphere-fixed solution can be used, (b) the ionosphere-float solution must be used, or (c) an ionosphere-weighted solution can be used.

This contribution is organized as follows. In Sect. 2 we revisit the ionosphere-weighted model and provide an alternative view of constructing its estimators. In particular, we take a deterministic view of the ionospheric constraints, while weighting them in the stochastic settings of the model. With the so obtained interpretation of the ionosphere-weighted estimator, its MSE can be evaluated, thus making it possible in Sect. 3 to formulate the general conditions the ionosphere-weighting need to obey in order for the ionosphere-weighted estimator to be MSE-superior to the ionosphere-float estimator. To develop these conditions further and to be able to express them in the metric of (1), we introduce in Sect. 4 the multi-frequency GNSS baseline model that forms the basis of our further MSE-analyses. As a result, we are able in Sect. 5 to find and formulate the various specific conditions the ionospheric variability needs to satisfy in order to guarantee improved MSE-performance of the ionosphere-weighted GNSS parameters over its ionosphere-float counterparts. Finally, a summary with conclusions is provided in Sect. 6.

The following notation will be used. We indicate random variables by means of an underscore. Thus, ${\underline{x}}$ is random, while x is not. ${\mathsf {E}}(\cdot )$ and ${\mathsf {D}}(\cdot )$ denote the expectation and dispersion operators, respectively. The square-weighted-norm is denoted as $||\cdot ||_{M}^{2}=(\cdot )^{\mathrm{T}}M^{-1}(\cdot )$. The matrix inequality $A \le B$ means that matrix $B-A$ is positive semi-definite, and $A \otimes B$ denotes the Kronecker product of A and B.

2 The iono-weighted model revisited

The ionosphere-weighted system of GNSS observation equations is given in linear(ized) model form as (Bock 1998; Teunissen 1998; Liu 2001; Odijk 2002),

$$\begin{aligned} {\mathsf {E}}\left[ \begin{array}{c} {\underline{y}} \\ {\underline{i}} \end{array} \right] = \left[ \begin{array}{c} A \\ A_{i} \end{array} \right] x \;\;,\;\; {\mathsf {D}}\left[ \begin{array}{c} {\underline{y}} \\ {\underline{i}} \end{array} \right] = \left[ \begin{array}{cc} Q_{yy} &{} 0\\ 0 &{} Q_{ii} \end{array} \right] \end{aligned}$$

(2)

with ${\underline{y}}\in {\mathbb {R}}^{m}$ the vector of GNSS pseudorange and carrier-phase observables, $A \in {\mathbb {R}}^{m \times n}$ the GNSS design matrix, $x=(b^{\mathrm{T}}, i^{\mathrm{T}})^{\mathrm{T}} \in {\mathbb {R}}^{n}$ the vector of unknown parameters, partitioned into the non-ionospheric parameter vector $b \in {\mathbb {R}}^{n-p}$ and the ionospheric delay vector $i \in {\mathbb {R}}^{p}$, and $Q_{yy}$ the variance-covariance (vc) matrix of ${\underline{y}}$. The p equations of ${\mathsf {E}}({\underline{i}})=A_{i}x$ form the ionospheric observation equations, with $Q_{ii}$ the vc-matrix of the ionospheric observable vector ${\underline{i}}$. The matrix $A_{i}$ can take various forms. Usually it will be given as $A_{i}=(0,I_{p})$. However, when the ionospheric equations pertain only to the time-variation of the ionospheric delays for instance, the identity-matrix $I_{p}$ will be replaced by a time-differencing matrix. Here and in the following we assume $\mathrm{rank}(A)=n$ and that both $Q_{yy}$ and $Q_{ii}$ are positive definite.

The ionosphere-weighted formulation is flexible and in use in a wide range of different GNSS applications. It is a form of regularized least-squares (Gunst and Mason 1980; Toutenburg 1982; Tarantola 2005), of which Tikhonov regularization (Tihonov 1963), ridge-regression (Hoerl and Kennard 1970), and principal component estimation (Jolliffe 2002) are prime examples. The ionosphere-weighted formulation is used in combination with external ionospheric models, see, e.g., (Schaer 1999; Memarzadeh 2009; Jee et al. 2010; Feltens et al. 2011; Olivares-Pulido et al. 2019), and to incorporate ionospheric corrections from reference networks, for instance for PPP or PPP-RTK, see, e.g., (Odijk 2000, 2002; Collins et al. 2012; Paziewski 2016; Odijk et al. 2016; Wang et al. 2018; Psychas and Verhagen 2020; Tomaszewski et al. 2020; Teunissen 2021). As the ionosphere is known to decorrelate as a function of baseline length, the ionosphere-weighted formulation is also used for strengthening medium to long baseline models, see, e.g., (Schaffrin and Bock 1988; Goad and Yang 1994; Bock 1998; Teunissen 1998; Odijk 2002; Odolinski and Teunissen 2017, 2019; Brack et al. 2021).

The flexibility of the ionosphere-weighted model can be made explicit when its best linear unbiased estimator (BLUE) of x, denoted as ${\hat{\underline{{\hat{x}}}}}$, is compared with its ionosphere-float counterpart $\underline{{\hat{x}}}$, which excludes the ionosphere-weighting and is thus solely based on ${\mathsf {E}}({\underline{y}})=Ax$, ${\mathsf {D}}({\underline{y}})=Q_{yy}$. As we have

$$\begin{aligned} \begin{array}{lcl} \hat{{\underline{x}}} &{}=&{} Q_{{\hat{x}}{\hat{x}}}[A^{\mathrm{T}}Q_{yy}^{-1}{\underline{y}}]\\ \hat{\hat{{\underline{x}}}} &{}=&{} Q_{\hat{{\hat{x}}}\hat{{\hat{x}}}}[A^{\mathrm{T}}Q_{yy}^{-1}{\underline{y}}+A_{i}^{\mathrm{T}}Q_{ii}^{-1}{\underline{i}}] \end{array} \end{aligned}$$

(3)

with vc-matrices $Q_{{\hat{x}}{\hat{x}}}\,=\,[A^{\mathrm{T}}Q_{yy}^{-1}A]^{-1}$ and $Q_{\hat{{\hat{x}}}\hat{{\hat{x}}}}\,=\,[A^{\mathrm{T}}\,Q_{yy}^{-1}A +A_{i}^{\mathrm{T}}Q_{ii}^{-1}A_{i}]^{-1}$, we may also write

$$\begin{aligned} \begin{array}{lcl} \hat{\hat{{\underline{x}}}} &{}\overset{(\mathrm{i})}{=}&{} \hat{{\underline{x}}} +[Q_{{\hat{x}}{\hat{x}}}^{-1}+A_{i}^{\mathrm{T}}Q_{ii}^{-1}A_{i}]^{-1}A_{i}^{\mathrm{T}}Q_{ii}^{-1}[{\underline{i}}-A_{i}\hat{{\underline{x}}}] \\ &{}\overset{(\mathrm{v})}{=}&{} \hat{{\underline{x}}} +Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}[Q_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}[{\underline{i}}-A_{i}\hat{{\underline{x}}}] \end{array} \end{aligned}$$

thus showing how $\hat{{\underline{x}}}$ gets updated as a consequence of the extra ionospheric equation ${\mathsf {E}}({\underline{i}})=A_{i}x$. Representation $(\mathrm{i})$ is the information form and $(\mathrm{v})$ the variance form of the update equation. Although both estimators, $\hat{\hat{{\underline{x}}}}$ and $\hat{{\underline{x}}}$, are BLUEs of x, and therefore unbiased, ${\mathsf {E}}(\hat{\hat{{\underline{x}}}})={\mathsf {E}}(\hat{{\underline{x}}})=x$, their precision differs. An application of the variance propagation law to the above variance-form gives

$$\begin{aligned} Q_{\hat{{\hat{x}}}\hat{{\hat{x}}}}=Q_{{\hat{x}}{\hat{x}}}-Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}Q_{vv}^{-1}A_{i}Q_{{\hat{x}}{\hat{x}}} \end{aligned}$$

(4)

with $Q_{vv}=Q_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}$ the vc-matrix of the predicted ionospheric residual ${\underline{v}}={\underline{i}}-A_{i}\underline{{\hat{x}}}$. Expression (4) shows how the precision of the ionosphere-weighted BLUE ${\hat{\underline{{\hat{x}}}}}$ lies in between that of the ionosphere-fixed solution ($Q_{ii}=0$) and the ionosphere-float solution ($Q_{ii}=\infty$). The precision improvement gets larger for $Q_{ii} \downarrow 0$, with in the limit providing that of the ionosphere-fixed solution, in which case the extra ionospheric equations have become hard constraints, as is used in case of short baselines for instance. The precision improvement of the estimator gets less, the less precise the extra ionospheric information is. In the limit of $Q_{ii} \uparrow$, one obtains the precision of the ionosphere-float solution, $Q_{\hat{{\hat{x}}}\hat{{\hat{x}}}} \rightarrow Q_{{\hat{x}}{\hat{x}}}$, as used for very long baselines for instance.

A basic assumption of the above ionosphere-weighted formulation is that

$$\begin{aligned} {\mathsf {E}}({\underline{i}})=A_{i}x\;\;\mathrm{and}\;\; {\mathsf {D}}({\underline{i}})=Q_{ii} \end{aligned}$$

(5)

must hold. If these assumptions fail to hold, then also the BLUE’s properties that are assigned to ${\hat{\underline{{\hat{x}}}}}$ fail to hold. As there are situations in practice for which it may be difficult to justify the strict assumptions of (5), it is important to understand to what extent these assumptions can be relaxed. For instance, it may not always be justified to treat the externally provided ionospheric information as a random variable ${\underline{i}}$, or if it is justified, one may have difficulty to justify that it is unbiased. The practice, for instance, of treating the differential ionospheric delay, in the absence of external information, as a pseudo-observable for which always the value zero is taken as sample-value, indicates that the assumption of randomness and unbiasedness could sometimes be considered a bit of a stretch.

To relax the assumptions of (5), we will refrain in the following from working with the random variable ${\underline{i}}$. Instead we assume to have a deterministic approximation, say $i_{0}$, of $A_{i}x$ available, the contribution of which will be weighted in the solution. That is, as one can expect $i_{0} \ne A_{i}x$, we will not enforce $i_{0}=A_{i}x$ as a hard constraint, but instead use a weight-matrix ${\bar{Q}}_{ii}^{-1}$ to weigh the difference $d=i_{0}-A_{i}x$ in a least-squares sense when solving for the unknown parameter vector x. As a result, the ionosphere-weighting will be achieved by working with an estimator having the same structure as ${\hat{\underline{{\hat{x}}}}}$ of (3), but one in which the random variable ${\underline{i}}$ is replaced by the deterministic approximation $i_{0}$ and the vc-matrix $Q_{ii}$ is replaced by the inverted weight-matrix ${\bar{Q}}_{ii}$. The approximation $i_{0}$ can often be obtained from an a priori ionospheric model, $i_{0}=A_{i}x_{0}$, such that it is the difference $d=-A_{i}(x-x_{0})$ that is weighted in the solution.

3 Improving precision at the cost of bias

3.1 Precision and bias

When replacing ${\underline{i}}$ and $Q_{ii}$ in the information- and variance-form of (3) by $i_{0}$ and ${\bar{Q}}_{ii}$, respectively, we obtain the ionosphere-weighted estimator

$$\begin{aligned} \begin{array}{lcl} \hat{{\underline{x}}}_{w} &{}\overset{(\mathrm{i})}{=}&{} \hat{{\underline{x}}} +[Q_{{\hat{x}}{\hat{x}}}^{-1}+A_{i}^{\mathrm{T}}{\bar{Q}}_{ii}^{-1}A_{i}]^{-1}A_{i}^{\mathrm{T}}{\bar{Q}}_{ii}^{-1}[i_{0}-A_{i}\hat{{\underline{x}}}] \\ &{}\overset{(\mathrm{v})}{=}&{} \hat{{\underline{x}}} +Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}[{\bar{Q}}_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}[i_{0}-A_{i}\hat{{\underline{x}}}] \end{array} \end{aligned}$$

(6)

When applying the propagation laws for means and variances to (6), we now have to recognize that $i_{0}$ is nonrandom (i.e., the only random vector in (6) is $\underline{{\hat{x}}}$) and that $i_{0}$ is not necessarily equal to the expectation of $A_{i}\underline{{\hat{x}}}$. Application of these propagation laws to (6) gives (Appendix)

$$\begin{aligned} \begin{array}{lcl} {\mathsf {E}}(\hat{{\underline{x}}}_{w}) &{}=&{} x+M d \\ {\mathsf {D}}(\hat{{\underline{x}}}_{w}) &{}=&{}Q_{{\hat{x}}{\hat{x}}}-M(2{\bar{Q}}_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}})M^{\mathrm{T}} \end{array} \end{aligned}$$

(7)

with $d=i_{0}-A_{i}x$ and $M=Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}[{\bar{Q}}_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}$. This result shows that $\hat{{\underline{x}}}_{w}$ has a better precision than $\underline{{\hat{x}}}$. In fact, since $i_{0}$ is nonrandom, it even has a better precision than the ionosphere-weighted estimator ${\hat{\underline{{\hat{x}}}}}$. However, unlike ${\hat{\underline{{\hat{x}}}}}$, the estimator $\underline{{\hat{x}}}_{w}$ is biased, having as bias, $b_{{\hat{x}}_{w}}={\mathsf {E}}(\hat{{\underline{x}}}_{w})-x=M d$. Hence, in order to judge the performance of $\underline{{\hat{x}}}_{w}$, we cannot rely solely on its precision. As we need to include, next to its variance, also its bias, we will use the mean-squared-error criterion for comparing the performance of estimators.

Recall that for an estimator $\hat{{\underline{z}}}$ of $z \in {\mathbb {R}}^{n}$, the mean-squared-error matrix (MSEM) and mean squared error (MSE) are defined as

$$\begin{aligned} \begin{array}{lcl} \mathrm{MSEM}(\hat{{\underline{z}}})&{}=&{}{\mathsf {E}}((\hat{{\underline{z}}}-z)(\hat{{\underline{z}}}-z)^{\mathrm{T}})\\ \mathrm{MSE}(\hat{{\underline{z}}})&{}=&{}{\mathsf {E}}((\hat{{\underline{z}}}-z)^{\mathrm{T}}(\hat{{\underline{z}}}-z)). \end{array} \end{aligned}$$

(8)

Thus, $\mathrm{MSE}(\hat{{\underline{z}}})={\mathsf {E}}(||\hat{{\underline{z}}}-z||^{2})=\mathrm{trace}(\mathrm{MSEM}(\hat{{\underline{z}}}))$. The MSE can be decomposed into a variance-plus-$\hbox {bias}^{2}$ contribution as

$$\begin{aligned} \begin{array}{lcl} \mathrm{MSEM}(\hat{{\underline{z}}})&{}=&{}Q_{{\hat{z}}{\hat{z}}}+b_{{\hat{z}}}b_{{\hat{z}}}^{\mathrm{T}}\\ \mathrm{MSE}(\hat{{\underline{z}}}) &{}=&{} \mathrm{trace}(Q_{{\hat{z}}{\hat{z}}})+||b_{{\hat{z}}}||^{2}. \end{array} \end{aligned}$$

(9)

Now let $\hat{{\underline{z}}}_{1}$ and $\hat{{\underline{z}}}_{2}$ be two estimators of $z \in {\mathbb {R}}^{n}$, and let us be interested in estimating the scalar $\theta =f^{\mathrm{T}}z$. Then, using the MSE-criterion, one would prefer $\hat{{\underline{\theta }}}_{1}$ over $\hat{{\underline{\theta }}}_{2}$ if $\mathrm{MSE}(\hat{{\underline{\theta }}}_{1}) \le \mathrm{MSE}(\hat{{\underline{\theta }}}_{2})$, and thus if $f^{\mathrm{T}}\mathrm{MSEM}(\hat{{\underline{z}}}_{1})f \le f^{\mathrm{T}}\mathrm{MSEM}(\hat{{\underline{z}}}_{2})f$. If this inequality is fulfilled for every $f \in {\mathbb {R}}^{n}$, then

$$\begin{aligned} \mathrm{MSEM}(\hat{{\underline{z}}}_{1}) \le \mathrm{MSEM}(\hat{{\underline{z}}}_{2}) \end{aligned}$$

(10)

meaning that the matrix $\mathrm{MSEM}(\hat{{\underline{z}}}_{2})-\mathrm{MSEM}(\hat{{\underline{z}}}_{1})$ is non-negative definite. The estimator $\hat{{\underline{z}}}_{1}$ is said to be MSE-superior over $\hat{{\underline{z}}}_{2}$ if their MSEMs satisfy the matrix condition (10). Thus, every linear function of an MSE-superior vector estimator has an MSE that is not larger than the one with which it is compared.

3.2 Iono-weighted MSE

We can now use the MSE-criterion to compare the MSE-performance of the ionosphere-weighted estimator $\hat{{\underline{x}}}_{w}$ with that of the ionosphere-float BLUE $\hat{{\underline{x}}}$.

Theorem 1a

(Iono-weighted MSE) The biased ionosphere-weighted estimator ${\hat{x}}_{w}$ is MSE-superior to the ionosphere-float BLUE ${\hat{x}}$, i.e., $\mathrm{MSEM}(\hat{{\underline{x}}}_{w}) \le \mathrm{MSEM}(\hat{{\underline{x}}})$, if and only if

$$\begin{aligned} d^{\mathrm{T}}[2 {\bar{Q}}_{ii}+Q_{{\hat{i}}{\hat{i}}}]^{-1}d \le 1 \end{aligned}$$

(11)

with $d=i_{0}-A_{i}x$ and $Q_{{\hat{i}}{\hat{i}}}=A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}$.

Proof

See Appendix.$\square$

This result provides the condition that the ionospheric approximation $i_{0}$ of $A_{i}x$ and the used weight-matrix ${\bar{Q}}_{ii}^{-1}$ need to satisfy in order to guarantee that the ionosphere-weighted estimator $\underline{{\hat{x}}}_{w}$ is MSE-superior to the ionosphere-float solution $\underline{{\hat{x}}}$.

Condition (11) describes an origin-centered ellipsoidal region in which the difference $d=i_{0}-A_{i}x$ has to lie. Next to $i_{0}$, one has the opportunity to choose an appropriate inverted weight-matrix ${\bar{Q}}_{ii}$, to tune the size and orientation of this ellipsoidal region. It is most strict if ${\bar{Q}}_{ii}=0$, which corresponds with a hard constraint enforcing $i_{0}=A_{i}\hat{{\underline{x}}}_{w}$, thus producing the ionosphere-fixed solution. The larger ${\bar{Q}}_{ii}$ is chosen, the larger the size of the ellipsoid and the larger the bias that is permitted for $\hat{{\underline{x}}}_{w}$ to be still superior to $\underline{{\hat{x}}}$. And in the limit, if ${\bar{Q}}_{ii}=\infty$ or ${\bar{Q}}_{ii}^{-1}=0$, we have $\hat{{\underline{x}}}_{w}=\hat{{\underline{x}}}$, being the ionosphere-float solution.

One may wonder, since the difference $d=i_{0}-A_{i}x$ is causing $\underline{{\hat{x}}}_{w}$ to be biased, how the ionosphere-weighted estimator ${\hat{\underline{{\hat{x}}}}}$ compares with $\underline{{\hat{x}}}_{w}$ and $\underline{{\hat{x}}}$, if the same bias $d={\mathsf {E}}({\underline{i}})-A_{i}x$ would be allowed to impact ${\hat{\underline{{\hat{x}}}}}$. Application of the mean and variance propagation laws to the variance-form of ${\hat{\underline{{\hat{x}}}}}$ gives then

$$\begin{aligned} \begin{array}{lcl} {\mathsf {E}}({\hat{\underline{{\hat{x}}}}}) &{}=&{} x+M d \\ {\mathsf {D}}({\hat{\underline{{\hat{x}}}}}) &{}=&{}Q_{{\hat{x}}{\hat{x}}}-M(Q_{ii}+Q_{{\hat{i}}{\hat{i}}})M^{\mathrm{T}} \end{array} \end{aligned}$$

(12)

with $d={\mathsf {E}}({\underline{i}})-A_{i}x$. If we compare this result with that of (7), we see that ${\hat{\underline{{\hat{x}}}}}$ and $\underline{{\hat{x}}}_{w}$ now have the same bias. Their precision differs however. The vc-matrix of ${\hat{\underline{{\hat{x}}}}}$ is obtained from the vc-matrix expression of $\underline{{\hat{x}}}_{w}$ if we replace $2{\bar{Q}}_{ii}$ by $Q_{ii}$. This implies that the condition for ${\hat{\underline{{\hat{x}}}}}$ to be MSE-superior to $\underline{{\hat{x}}}$ becomes $d^{\mathrm{T}}[Q_{ii}+Q_{{\hat{i}}{\hat{i}}}]^{-1}d \le 1$, with $d={\mathsf {E}}({\underline{i}})-A_{i}x$. This is clearly a more stringent condition than condition (11) for $\underline{{\hat{x}}}_{w}$. Moreover, with (11) one is not restricted to have $Q_{ii}$ be the vc-matrix of ${\underline{i}}$, i.e., one is still free to choose an appropriate inverted weight-matrix ${\bar{Q}}_{ii}$.

As condition (11) achieves MSE-superiority for all functions of $\underline{{\hat{x}}}_{w}$, thus also for functions that one may not necessarily be interested in, it may be considered too restrictive for all purposes. To relax the condition, we now consider only functions that one may be interested in, say $\theta =F^{\mathrm{T}}x$, with $F \in {\mathbb {R}}^{n \times q}$, and instead of requiring $\mathrm{MSEM}(\underline{{\hat{x}}}_{w}) \mathrm{MSEM}(\underline{{\hat{x}}})^{-1} \le I_{n}$, as in Theorem 1a, we now only require ${\mathsf {E}}(||\underline{{\hat{\theta }}}_{w}-\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})=\mathrm{trace}[\mathrm{MSEM}(\underline{{\hat{\theta }}}_{w}) \mathrm{MSEM}(\underline{{\hat{\theta }}})^{-1}] \le q$. We have the following result.

Theorem 1b

(Iono-weighted MSE) The MSE of the ionosphere-weighted estimator $\underline{{\hat{\theta }}}_{w}=F^{\mathrm{T}}\underline{{\hat{x}}}_{w}$ is not larger than that of the ionosphere-float $\underline{{\hat{\theta }}}=F^{\mathrm{T}}\underline{{\hat{x}}}$, i.e., $\mathrm{MSE}(\hat{{\underline{\theta }}}_{w}) \le \mathrm{MSE}(\hat{{\underline{\theta }}})$, with $\mathrm{MSE}(\cdot )={\mathsf {E}}(||\cdot -\,\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})$, if and only if

$$\begin{aligned} \frac{d^{\mathrm{T}}Xd}{\mathrm{trace}([2{\bar{Q}}_{ii}+Q_{{\hat{i}}{\hat{i}}}]X)} \le 1\;\;\mathrm{or}\;\;X=0 \end{aligned}$$

(13)

with $X=[{\bar{Q}}_{ii}+Q_{{\hat{i}}{\hat{i}}}]^{-1}D_{{\hat{i}}{\hat{i}}}[{\bar{Q}}_{ii}+Q_{{\hat{i}}{\hat{i}}}]^{-1}$, $Q_{{\hat{i}}{\hat{i}}}=A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}$, $Q_{{\hat{i}}{\hat{i}}|F^{\mathrm{T}}x}=A_{i}Q_{{\hat{x}}{\hat{x}}|F^{\mathrm{T}}x}A_{i}^{\mathrm{T}}$, $D_{{\hat{i}}{\hat{i}}}=Q_{{\hat{i}}{\hat{i}}}-Q_{{\hat{i}}{\hat{i}}|F^{\mathrm{T}}x}$, and $d=i_{0}-A_{i}x$.

Proof

See Appendix.$\square$

Compare this result with (11). Condition (13) is a more relaxed condition than (11). The latter follows if the inequality of (13) would hold for all $X\ne 0$. To operationalize condition (13), we need to make choices for F, $A_{i}$ and ${\bar{Q}}_{ii}$. By choosing F, one selects the parameters one is interested in. Its influence on the above condition is only felt through $D_{{\hat{i}}{\hat{i}}}=Q_{{\hat{i}}{\hat{i}}}-A_{i}Q_{{\hat{x}}{\hat{x}}|F^{\mathrm{T}}x}A_{i}^{\mathrm{T}}$, with $Q_{{\hat{x}}{\hat{x}}|F^{\mathrm{T}}x}$ being the conditional vc-matrix, i.e., the vc-matrix that follows when the solution is constrained by $F^{\mathrm{T}}x=0$. Thus, when all parameters are selected, we have $F=I_{n}$, and therefore $Q_{{\hat{x}}{\hat{x}}|F^{\mathrm{T}}x}=0$, which gives $D_{{\hat{i}}{\hat{i}}}=Q_{{\hat{i}}{\hat{i}}}$. This same result follows when F only selects the ionospheric parameters, since then $A_{i}Q_{{\hat{x}}{\hat{x}}|F^{\mathrm{T}}x}A_{i}^{\mathrm{T}}=0$. It also follows from (13) that the condition is trivially fulfilled when $F^{\mathrm{T}}\underline{{\hat{x}}}$ and $A_{i}\underline{{\hat{x}}}$ would be uncorrelated, in which case $A_{i}Q_{{\hat{x}}{\hat{x}}}F=0$ and thus $D_{{\hat{i}}{\hat{i}}}=0$. In order to show now how the above condition works out for different choices of F, $A_{i}$ and ${\bar{Q}}_{ii}$, we first formulate our GNSS baseline model and determine the vc-matrix of its ionospheric float solution.

4 The GNSS baseline model

In this section, we describe the multi-frequency GNSS baseline model that will form the basis of our MSE analysis of the ionosphere-weighted estimator.

4.1 Multi-epoch, multi-frequency GNSS

The m-satellite, f-frequency GNSS baseline model of epoch t reads (Khodabandeh and Teunissen 2015; Odijk et al. 2016)

$$\begin{aligned} \begin{array}{l} {\mathsf {E}}\left[ \begin{array}{c} {\underline{\phi }}_t \\ {\underline{p}}_t \end{array}\right] = \left[ \begin{array}{lr} A_{\phi _{t}} &{} \;\; - \mu \otimes I_{m-1}\\ A_{p_{t}} &{} \;\; +\mu \otimes I_{m-1}\end{array}\right] \left[ \begin{array}{c} x \\ i_t \end{array}\right] \quad \\ {\mathsf {D}}\left[ \begin{array}{c} {\underline{\phi }}_t \\ {\underline{p}}_t \end{array}\right] = \left[ \begin{array}{cc} Q_\phi &{} 0\\ 0 &{} Q_p\end{array}\right] \otimes R_t \end{array} \end{aligned}$$

(14)

with

$$\begin{aligned} \left[ \begin{array}{c} A_{\phi _{t}}\\ A_{p_{t}} \end{array} \right] x = \left[ \begin{array}{cc} e_{f}\otimes D^{\mathrm{T}}G_t &{}\; \Lambda \otimes I_{m-1}\\ e_{f}\otimes D^{\mathrm{T}}G_t &{}\; 0 \end{array} \right] \left[ \begin{array}{c} b\\ a \end{array} \right] \end{aligned}$$

and where ${\underline{\phi }}_{t}, {\underline{p}}_{t} \in {\mathbb {R}}^{f(m-1)}$ denote the vectors of double differenced (DD) phase and code observables, $\Lambda =\mathrm{diag}(\lambda _1,\ldots ,\lambda _f)$ is the diagonal wavelength matrix, $\mu =(\mu _{1}, \ldots , \mu _{f})^{\mathrm{T}}$ with the square wavelength-ratio $\mu _{j}=(\lambda _{j}/\lambda _{1})^{2}$, $e_{f}=(1, \ldots , 1)^{\mathrm{T}}$ the f-vector of ones, $D^{\mathrm{T}}=(-e_{m-1}, I_{m-1}) \in {\mathbb {R}}^{(m-1) \times m}$, with $e_{m-1}=(1, \ldots , 1)^{\mathrm{T}}$, the differencing matrix that forms satellite-differences from between-receiver differences, $G_{t} \in {\mathbb {R}}^{m \times 3}$ the receiver-satellite geometry matrix, b the between-receiver baseline vector, $i_{t} \in {\mathbb {R}}^{m-1}$ the vector of DD slant ionospheric delays at epoch t, $a\in \mathbb {Z}^{f(m-1)}$ the vector of DD ambiguities, $Q_{\phi }, Q_{p} \in {\mathbb {R}}^{f \times f}$ the f-frequency co-factor matrices of phase and code, $R_t = D^{\mathrm{T}}W_t^{-1}D \in {\mathbb {R}}^{(m-1)\times (m-1)}$ the co-factor matrix, with $W_t=\mathrm{diag}(w_t^1, \ldots ,w_t^m)$, the diagonal matrix that captures the satellite elevation-dependency of the GNSS data through the positive weights $w_t^s$ ($s=1,\ldots ,m$). In all time-dependent variables, we use the index t to emphasize the time dependency of the quantities. It is absent in b and a, since we assume a stationary baseline and ambiguities that are time-constant.

In this and the next section, model (14) will be used to study the ionosphere-weighted MSE properties. As we will be focused on short time-spans and since the receiver-satellite geometry of GNSS changes rather slowly with time, we will assume $W_t\approx W$ and $G_t \approx G$ ($t=1,\ldots ,k$), with W and G being the ‘time-averaged’ versions of $W_t$ and $G_t$, respectively. To achieve a compact description and aid the interpretation of our results, we use in the remainder of this contribution the following notation,

$$\begin{aligned} \begin{array}{l} \left[ \begin{array}{ccc} \sigma _{{\hat{\rho }}}^2&{}&{} \sigma _{{\hat{\rho }}{\hat{i}}}\quad \\ \sigma _{{\hat{\rho }}{\hat{i}}} &{}&{} \sigma _{{\hat{i}}}^2 \end{array}\right] =[[e_{f},\mu ]^{\mathrm{T}}Q_p^{-1}[e_{f},\mu ]]^{-1}\\ \left[ \begin{array}{ccc} \sigma _{{\check{\rho }}}^2 &{}&{} \sigma _{{\check{\rho }}{\check{i}}}\quad \\ \sigma _{{\check{\rho }}{\check{i}}} &{}&{} \sigma _{{\check{i}}}^2 \end{array}\right] =\left[ \left[ \begin{array}{cc}e_{f} &{} -\mu \quad \\ e_{f} &{} +\mu \end{array}\right] ^{\mathrm{T}} \!\!\left[ \begin{array}{cc}Q_{\phi }^{-1} &{} 0 \quad \\ 0 &{} Q_p^{-1} \end{array}\right] \left[ \begin{array}{cc}e_{f} &{} -\mu \quad \\ e_{f} &{} +\mu \end{array}\right] \right] ^{-1}. \end{array} \end{aligned}$$

(15)

One will recognize these variance and co-variance factors as the ones that belong to the code-only geometry-free model and the ambiguity-fixed geometry-free model, respectively.

4.2 Precision of estimated ionospheric delays

As the vc-matrix of the estimated ionospheric delays is the driving force in (11) and (13), we now first consider the precision with which the ionospheric delays can be estimated.

Lemma

(VC-matrix ionosphere) The $k(m-1) \times k(m-1)$ vc-matrices of the ambiguity-float and ambiguity-fixed BLUEs of $i_{[k]}=[i_1^{\mathrm{T}},\ldots ,i_k^{\mathrm{T}}]^{\mathrm{T}}$ in model (14), are given as

$$\begin{aligned} \begin{array}{lcl} Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}} &{}=&{}P_{k}\otimes (\alpha -\beta )\,{\mathcal {P}}_{\!{\tilde{G}}} R +Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}|b}\\ Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}|b}&{}=&{} P_{k}\otimes \beta \,R+P_{k}^{\perp }\otimes \,\sigma _{{\check{i}}|\rho }^2\,R \end{array} \end{aligned}$$

(16)

with

$$\begin{aligned} \begin{array}{ll} \alpha = \sigma _{{\hat{i}}}^2\;\;\mathrm{and}\;\;\beta =\sigma _{{\hat{i}}|\rho }^2 &{} (\mathrm{ambiguity\;float})\\ \alpha =\sigma _{{\check{i}}}^2\;\;\mathrm{and}\;\;\beta =\sigma _{{\check{i}}|\rho }^2&{} (\mathrm{ambiguity fixed}) \end{array} \end{aligned}$$

(17)

the time-averaging projector $P_{k}=\frac{1}{k}e_ke_k^{\mathrm{T}}$ ($P_{k}^{\perp }=I_k-P_{k}$), the satellite-geometry projector ${\mathcal {P}}_{\,{\tilde{G}}} = {\tilde{G}}({\tilde{G}}^{\mathrm{T}}R^{-1}{\tilde{G}})^{-1} {\tilde{G}}^{\mathrm{T}}R^{-1}$ ($P_{{\tilde{G}}}^{\perp }=I-P_{{\tilde{G}}}$), the receiver-satellite geometry matrix ${\tilde{G}} = D^{\mathrm{T}}G$, the satellite-elevation weighting DD matrix $R=D^{\mathrm{T}}W^{-1}D$ and the conditional variance-factors $\sigma _{{\hat{i}}|\rho }^{2}=\sigma _{{\hat{i}}}^{2}-\sigma _{{\hat{i}}{\hat{\rho }}}^{2}/\sigma _{{\hat{\rho }}}^{2}$ and $\sigma _{{\check{i}}|\rho }^{2}=\sigma _{{\check{i}}}^{2}-\sigma _{{\check{i}}{\check{\rho }}}^{2}/\sigma _{{\check{\rho }}}^{2}$.

Proof

See Appendix. $\square$

The above result shows that the vc-matrix of the estimated ionospheric delays has a very patterned structure. Matrix $Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}|b}$ is the ionospheric vc-matrix when the baseline b is assumed known. When $Q_{\phi }=\sigma _{\phi }^{2}I_{f}$ and $Q_{p}=\sigma _{p}^{2}I_{f}$, the variance-factors in (16) and (17) simplify to

$$\begin{aligned} \sigma _{{\hat{i}}}^2=\frac{\sigma _{p}^{2}}{||\mu -{\bar{\mu }}||^{2}}, \sigma _{{\hat{i}}|\rho }^2=\frac{\sigma _{p}^{2}}{||\mu ||^{2}}, \sigma _{{\check{i}}|\rho }^2=\frac{\sigma _{\phi }^{2}}{(1+\epsilon )||\mu ||^{2}} \end{aligned}$$

where ${\bar{\mu }}= e_{f} \tfrac{1}{f}\sum _{j=1}^{f}\mu _{j}$ is the f-vector having the average value of the $\mu _{j}$’s as its entries and $\epsilon = \sigma _{\phi }^{2}/\sigma _{p}^{2}$ is the very small phase-code variance-ratio. The phase and code variance factors $\sigma _{\phi }^{2}$ and $\sigma _{p}^{2}$ are between-receiver single-differenced values.

The above analytical expression of the ionospheric vc-matrix shows how the various factors contribute to the precision of the estimated ionospheric delays. This becomes even clearer when we consider the vc-matrix of the ionospheric delays at a specific epoch t. From (16), the vc-matrix of $\hat{{\underline{i}}}_t$ follows as

$$\begin{aligned} Q_{{\hat{i}}_t{\hat{i}}_t} = \underbrace{\tfrac{1}{k}(\sigma _{{\hat{i}}}^2\,{\mathcal {P}}_{\!{\tilde{G}}} R+\sigma _{{\hat{i}}|\rho }^2\,{\mathcal {P}}_{\!{\tilde{G}}}^\bot R)}_{\mathrm{time-averaged}} +\underbrace{\tfrac{k-1}{k}\sigma _{{\check{i}}|\rho }^2\,R}_{\mathrm{time-differenced}}. \end{aligned}$$

This shows that the vc-matrix is composed of two components, a time-averaged component and a time-differenced component. The time-averaged component is driven, through $\sigma _{{\hat{i}}}^2$ and $\sigma _{{\hat{i}}|\rho }^2$, by the code-precision, and through ${\mathcal {P}}_{\!{\tilde{G}}}$ and R, by the receiver-satellite geometry. The time-differenced component however, is driven through $\sigma _{{\check{i}}|\rho }^2$ by the phase-precision and only through the satellite-elevation dependency of R by the receiver-satellite geometry. Also note that in the geometry-free case or when $m=4$, the ionospheric vc-matrix simplifies to $Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}} =P_{k}\otimes \alpha \,R + P_{k}^{\perp }\otimes \,\sigma _{{\check{i}}|\rho }^2\,R$, since then $P_{{\tilde{G}}}^{\perp }=0$ and $P_{{\tilde{G}}}=I_{m-1}$.

5 The MSE ionosphere conditions

In this section, we will use (11) and (13) to further develop the variability conditions the ionosphere has to satisfy in order for the ionosphere-weighted estimator (6) to have a smaller mean squared error than its ionosphere-float counterpart.

5.1 Iono-weighting as a convex combination

In order to evaluate the performance of the ionosphere-weighted estimator $\underline{{\hat{x}}}_{w}$ we need to make a choice for the inverted weight-matrix ${\bar{Q}}_{ii}$. It is through the choice of this matrix that one is able to influence the mean squared error of $\underline{{\hat{x}}}_{w}$ (cf. 6). Although we have in principle complete freedom in choosing the entries of matrix ${\bar{Q}}_{ii}$, the simplest choice would be to take it as a scaled version of the ionospheric vc-matrix,

$$\begin{aligned} {\bar{Q}}_{ii}=\lambda Q_{{\hat{i}}{\hat{i}}}=\lambda A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}. \end{aligned}$$

(18)

This choice also provides for an easy interpretation of how $\underline{{\hat{x}}}_{w}$ is constructed from the data. From the information-form of (6) follows then, since $Q_{{\hat{i}}{\hat{i}}}^{-1}$ is the reduced normal matrix for the ionospheric delays, that a scaled version of this reduced matrix is added to the normal matrix of the ionospheric delays so as to compute the ionosphere-weighted solution. From the variance-form of (6) follows that the choice ${\bar{Q}}_{ii}=\lambda Q_{{\hat{i}}{\hat{i}}}$ gives $\underline{{\hat{x}}}_{w}(\lambda )=\underline{{\hat{x}}}+\tfrac{1}{\lambda +1}Q_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}[A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}[i_{0}-A_{i}\underline{{\hat{x}}}]$, thus showing how the contribution to $\underline{{\hat{x}}}_{w}(\lambda )$ of the difference $i_{0}-A_{i}\underline{{\hat{x}}}$ is scaled through the nonnegative scale-factor $\lambda$. It therefore follows, since $Q_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}[A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}[i_{0}-A_{i}\underline{{\hat{x}}}]$ equals the difference between the ionosphere-float and the ionosphere-fixed solutions, that with the choice ${\bar{Q}}_{ii}=\lambda Q_{{\hat{i}}{\hat{i}}}$, the ionosphere-weighted solution equals their convex combination,

$$\begin{aligned} \underline{{\hat{x}}}_{w}(\lambda )= \tfrac{\lambda }{\lambda +1} \underline{{\hat{x}}}+ \tfrac{1}{\lambda +1} \underline{{\hat{x}}}_{fixed}. \end{aligned}$$

(19)

Hence, it is through $\lambda$ that one can have $\underline{{\hat{x}}}_{w}(\lambda )$ inherit properties from both $\underline{{\hat{x}}}$ and $\underline{{\hat{x}}}_{fixed}$. For smaller values of $\lambda$, $\underline{{\hat{x}}}_{w}(\lambda )$ will tend toward the ionosphere-fixed solution and therefore be more precise at the expense of some bias. For larger values of $\lambda$, $\underline{{\hat{x}}}_{w}(\lambda )$ will tend toward the ionosphere-float solution and therefore become less biased at the expense of a poorer precision. The following result shows what values for $\lambda$ need to be chosen in order to ensure that $\underline{{\hat{x}}}_{w}(\lambda )$ has the best mean squared error or one that is at least better than that of the ionosphere-float solution.

Theorem 2

(Iono-weighted minimum MSE) Let ${\bar{Q}}_{ii}=\lambda Q_{{\hat{i}}{\hat{i}}}$,  $\hat{{\underline{\theta }}}_{w}(\lambda )=F^{\mathrm{T}}\underline{{\hat{x}}}_{w}(\lambda )$ and $\mathrm{MSE}(\lambda )={\mathsf {E}}(||\hat{{\underline{\theta }}}_{w}(\lambda ) -\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})$. Then, $\lambda _{\mathrm{min}}=\arg \min _{\lambda } \mathrm{MSE}(\lambda )$ is given as

$$\begin{aligned} \lambda _{\mathrm{min}} = \frac{d^{\mathrm{T}}[Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}}Q_{{\hat{i}}{\hat{i}}}^{-1}]d}{\mathrm{trace}[Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}}]} \end{aligned}$$

(20)

and $\mathrm{MSE}(\lambda ) \le \mathrm{MSE}(\infty )$ is satisfied for

$$\begin{aligned} \lambda \ge \tfrac{1}{2}\left( \lambda _{\mathrm{min}} -1 \right) \end{aligned}$$

(21)

Proof

See Appendix.$\square$

This result shows that if we are able to evaluate the ratio (20), we would know how to scale $Q_{{\hat{i}}{\hat{i}}}$ such that the ionosphere-weighted solution has the best possible mean squared error. The result also shows the range of values we are permitted to choose for $\lambda$ so as to still get an estimator with a mean squared error that is not larger than that of the ionosphere-float solution. Furthermore, it follows from (21) that if $\lambda _{\mathrm{min}}<1$, then $\lambda$ can even be chosen equal to zero. Hence, $\lambda _{\mathrm{min}} <1$ is the important condition under which the hard-constrained ionosphere-fixed solution would be MSE-better than the ionosphere-float solution.

5.2 Ionospheric dispersion and MSE

In order to show how the mean squared error of the ionosphere-weighted estimator $\underline{{\hat{x}}}_{w}(\lambda )$ is driven by the actual variability in the between-receiver ionospheric delays, we need to provide a proper metric by which this variability is measured. For this purpose, we define the spatial and temporal variability or dispersion of the ionosphere as follows.

Fig. 2

The between-receiver spatial ionospheric dispersion ${\mathbf {\mathsf{{d}}}}_{s}^{2}=\sum _{s=1}^{m} w^{s}(\texttt {i}_{{\bar{t}}}^{s}-\bar{\texttt {i}}_{{\bar{t}}})^{2}/(m-1)$ as a mean sum-of-squares of residual between-receiver satellite ionospheric delays

Definition

(Ionospheric dispersion) Let $i_{t} \in {\mathbb {R}}^{m-1}$, for epochs $t=1, \ldots , k$, be the vectors of DD ionospheric delays and let $i_{{\bar{t}}}= \tfrac{1}{k} \sum _{t=1}^{k} i_{t}$ be their time average. Then, the spatial and temporal ionospheric dispersion, ${\mathbf {\mathsf{{d}}}}_{s}^{2}$ and ${\mathbf {\mathsf{{d}}}}_{k}^{2}$, are defined as

$$\begin{aligned} {\mathbf {\mathsf{{d}}}}_{s}^{2}= \frac{||i_{{\bar{t}}}||_{R}^{2}}{m-1}\;\;\mathrm{and}\;\;{\mathbf {\mathsf{{d}}}}_{k}^{2}= \frac{\sum _{t=1}^{k}||i_{t}-i_{{\bar{t}}}||_{R}^{2}}{(k-1)(m-1)}. \end{aligned}$$

(22)

Although perhaps not directly apparent, also the expression for ${\mathbf {\mathsf{{d}}}}_{s}^{2}$ can be written in the familiar dispersion form. To see this, let $\texttt {i}_{{\bar{t}}}^{s}$ be the between-receiver single-difference ionospheric delay of satellite s at epoch ${\bar{t}}$, and let $\bar{\texttt {i}}_{{\bar{t}}}=\sum _{s=1}^{m} w^{s}{} \texttt {i}_{{\bar{t}}}^{s}/[\sum _{s=1}^{m}w^{s}]$ be the weighted average of these delays over all m satellites. Then, with $\texttt {i}_{{\bar{t}}}=[\texttt {i}_{{\bar{t}}}^{1}, \ldots , \texttt {i}_{{\bar{t}}}^{m}]^{\mathrm{T}}$, we have $\texttt {i}_{{\bar{t}}}-e_{m}\bar{\texttt {i}}_{{\bar{t}}}=P_{e_{m}}^{\perp } \texttt {i}_{{\bar{t}}}$, where $P_{e_{m}}^{\perp }=I_{m}-e_{m}(e_{m}^{\mathrm{T}}We_{m})^{-1}e_{m}^{\mathrm{T}}W=W^{-1}D(D^{\mathrm{T}}W^{-1}D)^{-1}D^{\mathrm{T}}$, and therefore, since the DD ionospheric delay vector is expressed in its single-differenced counterpart as $i_{{\bar{t}}}=D^{\mathrm{T}}{} \texttt {i}_{{\bar{t}}}$, $\sum _{s=1}^{m} w^{s}(\texttt {i}_{{\bar{t}}}^{s}-\bar{\texttt {i}}_{{\bar{t}}})^{2}= (\texttt {i}_{{\bar{t}}}-e_{m}\bar{\texttt {i}}_{{\bar{t}}})^{\mathrm{T}}W(\texttt {i}_{{\bar{t}}}-e_{m}\bar{\texttt {i}}_{{\bar{t}}})=\texttt {i}_{{\bar{t}}}^{\mathrm{T}}D(D^{\mathrm{T}}W^{-1}D)^{-1}D^{\mathrm{T}}{} \texttt {i}_{{\bar{t}}}=i_{{\bar{t}}}^{\mathrm{T}}R^{-1}i_{{\bar{t}}}=||i_{{\bar{t}}}||_{R}^{2}$, thus showing that ${\mathbf {\mathsf{{d}}}}_{s}^{2}$ of (22) can indeed be written in the familiar dispersion form (see Fig. 2),

$$\begin{aligned} {\mathbf {\mathsf{{d}}}}_{s}^{2} = \sum _{s=1}^{m} w^{s}(\texttt {i}_{{\bar{t}}}^{s}-\bar{\texttt {i}}_{{\bar{t}}})^{2}/(m-1). \end{aligned}$$

(23)

Note that, with the weight $g=\mathrm{trace}(P_{{\tilde{G}}})/(m-1)$, the ionospheric dispersion can be decomposed further into the convex combination

$$\begin{aligned} {\mathbf {\mathsf{{d}}}}_{s}^{2} = g {\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}+(1-g){\mathbf {\mathsf{{d}}}}_{s}^{\perp 2} \end{aligned}$$

(24)

where ${\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}=\tfrac{||P_{{\tilde{G}}}i_{{\bar{t}}}||_{R}^{2}}{\mathrm{trace}(P_{{\tilde{G}}})}$ and ${\mathbf {\mathsf{{d}}}}_{s}^{\perp 2}=\tfrac{||P_{{\tilde{G}}}^{\perp }i_{{\bar{t}}}||_{R}^{2}}{\mathrm{trace}(P_{{\tilde{G}}}^{\perp })}$ describe the dispersion of the ionosphere parallel ($\Vert$) and orthogonal ($\perp$) to the receiver-satellite geometry, respectively. In the geometry-free case, we have ${\mathbf {\mathsf{{d}}}}_{s}^{2} = {\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}$, since then $P_{{\tilde{G}}}^{\perp }=0$.

With the above given definitions of ionospheric dispersion, we are now in a position to apply the results of Theorem 2 to our GNSS model. We consider the case $F=I$ (i.e., all parameters are included) for three different types of ionospheric constraints. The three considered ionospheric constraints are: (1) weighted constraining the temporal variability of the ionospheric delays $(D_{k}^{\mathrm{T}} \otimes I_{m-1})i_{[k]}$, in which $D_{k}^{\mathrm{T}}=(-e_{k-1}, I_{k-1}) \in {\mathbb {R}}^{(k-1) \times k}$ is the time-differencing matrix, (2) weighted constraining the time-average of the ionospheric delays $(\frac{1}{k}e_{k}^{\mathrm{T}} \otimes I_{m-1})i_{[k]}$, and (3) weighted constraining all ionospheric delays $(I_{k} \otimes I_{m-1})i_{[k]}$. The following result shows by how much, in each of these three cases, one would need to scale the vc-matrix $Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}}$, so as to obtain the smallest mean squared error for the ionosphere-weighted estimator $\underline{{\hat{x}}}_{w}$.

Corollary 1a

(All parameters) Let $F=I$ and let the three ionospheric constraints be given as $A_{i}x \overset{(1)}{\rightarrow } (D_{k}^{\mathrm{T}} \otimes I_{m-1})i_{[k]}$, $A_{i}x \overset{(2)}{\rightarrow } (\frac{1}{k}e_{k}^{\mathrm{T}} \otimes I_{m-1})i_{[k]}$, and $A_{i}x \overset{(3)}{\rightarrow } (I_{k} \otimes I_{m-1})i_{[k]}$, respectively. Then, the corresponding minimizer of (20), expressed in the ionospheric dispersions, is given as

$$\begin{aligned} \begin{array}{lll} (1) &{} \lambda _{\mathrm{min}}^{(1)}= {\mathbf {\mathsf{{d}}}}_{k}^{2}/\sigma _{{\check{i}}|\rho }^{2}&{}(\mathrm{temporal})\quad \\ (2) &{} \lambda _{\mathrm{min}}^{(2)}= \dfrac{g{\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}}{\sigma _{{\hat{i}}}^{2}/k}+\dfrac{(1-g){\mathbf {\mathsf{{d}}}}_{s}^{\perp 2}}{\sigma _{{\hat{i}}|\rho }^{2}/k}&{}(\mathrm{spatial})\quad \\ (3) &{} \lambda _{\mathrm{min}}^{(3)}= \frac{1}{k} \lambda _{\mathrm{min}}^{(2)}+\frac{k-1}{k} \lambda _{\mathrm{min}}^{(1)}.&{} \end{array} \end{aligned}$$

Proof

See Appendix.$\square$

These results clearly show how the spatial and temporal ionospheric variability affect the optimal choice for the scaling factor $\lambda$. The larger the variability, the larger the scale factor needs to be chosen and the less weight is then given to the ionospheric constraint. Note, in case an a priori model approximation $i_{0}=A_{i}x_{0} \ne 0$ is chosen, that the ionospheric delays, of which the above dispersions are taken, are with respect to the chosen model values. As the above results concern the overall mean squared error of all parameters, including those of the ionospheric delays, we now consider only the MSE of two separate set of parameters, the ambiguities, and the baseline.

Corollary 1b

(Ambiguities or baseline) Let F select either the ambiguities or the baseline and let $A_{i}x \rightarrow (I_{k} \otimes I_{m-1})i_{[k]}$. Then, the corresponding minimizer of (20), for the ambiguities denoted as $\lambda _{\mathrm{min}}^{(a)}$ and for the baseline as $\lambda _{\mathrm{min}}^{(b)}$, is given as

$$\begin{aligned} \begin{array}{lcl} \lambda _{\mathrm{min}}^{(a)} &{}=&{} c [\frac{\gamma _{1}{\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}}{\sigma _{{\hat{i}}}^{2}/k}+\frac{\gamma _{2}{\mathbf {\mathsf{{d}}}}_{s}^{\perp 2}}{\sigma _{{\hat{i}}|\rho }^{2}/k}] \\ \lambda _{\mathrm{min}}^{(b)}&{}=&{} \frac{{\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}}{\sigma _{{\hat{i}}}^{2}/k} \end{array} \end{aligned}$$

(25)

with $c=(\gamma _{1}+\gamma _{2})^{-1}$, $\gamma _{1}=(1-\epsilon _{1})g$, $\gamma _{2}=(1-\epsilon _{2})(1-g)$, $\epsilon _{1}=\sigma _{{\check{i}}}^{2}/\sigma _{{\hat{i}}}^{2}$, and $\epsilon _{2}=\sigma _{{\check{i}}|\rho }^{2}/\sigma _{{\hat{i}}|\rho }^{2}$.

Proof

See Appendix. $\square$

As both $\epsilon _{1}$ and $\epsilon _{2}$ are in the order of the very small phase-code variance ratio, we get with $\epsilon _{1} \approx \epsilon _{2} \approx 0$,

$$\begin{aligned} \lambda _{\mathrm{min}}^{(a)} \approx \frac{g{\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}}{\sigma _{{\hat{i}}}^{2}/k}+\frac{(1-g){\mathbf {\mathsf{{d}}}}_{s}^{\perp 2}}{\sigma _{{\hat{i}}|\rho }^{2}/k} \le \frac{{\mathbf {\mathsf{{d}}}}_{s}^{2}}{\sigma _{{\hat{i}}|\rho }^{2}/k} \end{aligned}$$

(26)

Note that the upper-bound is independent of the actual receiver-satellite geometry and therefore easier to work with.

Fig. 3

Baseline MSE-curve $\mathrm{MSE}^{(b)}(\lambda )$ as function of the scale-factor $\lambda$ for GPS [L1/L2], single- (top) and quadruple-epoch (bottom) cases. Each curve corresponds to a different level of the ionospheric spatial-dispersion ${\mathbf {\mathsf{{d}}}}_{s}^{2}$ in [$\hbox {m}^2$]. The code standard deviation is set to $\sigma _p=\sqrt{2}\times 0.2$ [m], with four visible satellites. The parts, highlighted in red, indicate $\mathrm{MSE}^{(b)}(\lambda ) > \mathrm{MSE}^{(b)}(\infty )$

Fig. 4

Estimated values of the baseline MSE-minimizer $\lambda _{\mathrm{min}}^{(b)}=({\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}/\sigma _{{\hat{i}}}^{2})$, with $k=1$, corresponding to the scatterplots of baseline CUT0–KELN (cf. Fig. 1)

The results of Corollary 1 show the best choice one can make for the scale-factor $\lambda$ under the given circumstances, and together with (21), they provide the range of values that one can choose from for the scale factor to still have a better MSE-performance than the ionosphere-float solution. As an example, consider the baseline MSE-curve as function of $\lambda$,

$$\begin{aligned} \mathrm{MSE}^{(b)}(\lambda )=\mathrm{MSE}^{(b)}(\infty )-\frac{2a[\lambda -\tfrac{1}{2}(\lambda _{\mathrm{min}}^{(b)}-1)]}{(\lambda +1)^{2}} \end{aligned}$$

with $a=(1-\tfrac{\sigma _{{\hat{i}}|\rho }^{2}}{\sigma _{{\hat{i}}}^{2}})\mathrm{trace}(P_{{\tilde{G}}})$. Its curves are shown in Fig. 3 for a range of different ionospheric dispersion values ${\mathbf {\mathsf{{d}}}}_{s}^{2}$. The range of $\lambda$-values for which ionosphere-weighting makes sense are clearly visible, as well as how the interval $[\tfrac{1}{2}(\lambda _{\mathrm{min}}^{(b)}-1), \infty )$ shifts to the right with increasing ionospheric dispersion. To know whether one would be allowed to use the ionosphere-fixed solution (i.e., $\lambda =0$), one needs to verify whether or not the condition $\lambda _{\mathrm{min}}^{(b)}<1$ holds true. As the above corollary shows this is the case when ${\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2} \le \sigma _{{\hat{i}}}^{2}/k$. Figure 4 presents values of $\lambda _{\mathrm{min}}^{(b)}$ corresponding to the positioning results of the baseline in Fig. 1. These values have been computed on the basis of the a priori known coordinates of the stations CUT0 and KELN. The values are evidently shown to be way smaller than one, i.e., $\lambda _{\mathrm{min}}^{(b)}<1$, thereby addressing why the corresponding ionosphere-fixed baseline solutions are MSE-superior to their ionosphere-float counterparts. We note that in practice one does not need a priori precise coordinates of the baseline receivers to determine the ionospheric dispersion. One can make use of neighboring stations in existing CORS networks, of which nowadays there are a multitude available. One can also make use of the geometry-free model, to determine the ionospheric dispersion, or one can use, even simultaneously, algorithms of variance component estimation to determine the ionospheric dispersion.

Figure 3 also highlights the role taken by the number of epochs k in weighing the difference $d=i_{0}-A_{i}x$. The more the number of epochs, the lower weight the difference d should be given so as to obtain an MSE-superior baseline estimator. This is further observed in the results of Table 1. The table presents root mean squared errors (RMSE) of the single-epoch ($k=1$) and 10-epoch ($k=10$) baseline solutions of two different baselines where, next to the ionosphere-float and -fixed cases, an ionosphere-weighted case ($\lambda =0.1$) is also considered. When only one epoch of data is considered, the ionosphere-fixed solution outperforms its ionosphere-weighted version in an MSE-sense. For the 10-epoch, long-baseline case (CUT0–EXMT) however, it is the ionosphere-weighted solution that has the smallest MSE.

Table 1 Role of the number of epochs (k) and baseline length: root mean squared errors (RMSE) of the code-only ionosphere-float, -weighted [$\lambda =0.1$], and -fixed positioning solutions [m] of baselines CUT0–MRO1 ($\sim$ 592 km) and CUT0–EXMT ($\sim$ 1125 km) on 5 January 2021 (12:00–15:00 Perth time, GPS L1/L2, 1 Hz)

The above results have been obtained by choosing the inverted weight-matrix ${\bar{Q}}_{ii}$ simply as a scaled version of the ionospheric vc-matrix. But other choices are of course still possible as well. Some of these choices may even give the same conditions on the ionospheric dispersion. For instance, when using a similar derivation as used above, it can be shown if the inverted weight-matrix is chosen as

$$\begin{aligned} {\bar{Q}}_{ii}={\bar{\lambda }} (\sigma _{{\hat{i}}}^{2} I_{k} \otimes R) \end{aligned}$$

(27)

that the scale-factor ${\bar{\lambda }}$ needs to lie within the same range $[\tfrac{1}{2}(\tfrac{{\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}}{\sigma _{{\hat{i}}}^{2}/k}-1), \infty )$ for the MSE of the ionosphere-weighted baseline to be smaller than its ionosphere-float counterpart.

6 Summary and conclusions

Although ionosphere-weighted GNSS parameter estimation is a popular technique for strengthening estimator performance in the presence of ionospheric delays, no provable rules currently exist that specify the needed weighting in dependence on ionospheric circumstances. The goal of the present contribution has therefore been to develop and present the ionospheric conditions that need to be satisfied in order to guarantee that the ionosphere-weighted solution is mean squared error superior to the ionosphere-float solution.

We developed general conditions that apply to all parameters of the model, as well as specific conditions that apply when one is interested in specific functions of the model parameters, such as the baseline or ambiguities, for example. We have shown how the required conditions can be expressed in the metric $\sum _{s=1}^{m} w^{s}(\texttt {i}_{t}^{s}-\bar{\texttt {i}}_{t})^{2}/(m-1)$, capturing the ionospheric variability of the between-receiver single-difference ionospheric delays. Next to the choice of parameters, we also presented the conditions for different types of ionospheric constraints.

When satisfied, the conditions provided guarantee the improved mean-squared-error performance of the estimator. More specifically, they show when (a) the ionosphere-fixed solution can be used, (b) the ionosphere-float solution must be used, or (c) an ionosphere-weighted solution can be used. For the baseline, for instance, it was shown that the scale factor (cf. 18 or 27) has to be larger than $\tfrac{1}{2}(\tfrac{{\mathbf {\mathsf{{d}}}}_{s}^{\Vert 2}}{\sigma _{{\hat{i}}}^{2}/k}-1)$ in order to have a guaranteed improved mean squared error.

Finally, we remark, although attention was focused in this contribution on ionospheric delays, that the presented methodology is general. This implies that the same approach can be used to develop equivalent conditions for other GNSS parameters of interest, such as for instance for the tropospheric delays or instrumental biases.

Appendix

Proof of (7)

Application of the propagation laws for means and variances to (6) gives ${\mathsf {E}}(\hat{{\underline{x}}}_{w}) = x+M {\mathsf {E}}[i_{0}-A_{i}\hat{{\underline{x}}}]$ and ${\mathsf {D}}(\hat{{\underline{x}}}_{w})=Q_{{\hat{x}}{\hat{x}}}+MA_iQ_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}M^{\mathrm{T}}-MA_iQ_{{\hat{x}}{\hat{x}}}-Q_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}M^{\mathrm{T}}$, where $M=Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}[{\bar{Q}}_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}$. The proof follows then from the equalities ${\mathsf {E}}[i_{0}-A_{i}\hat{{\underline{x}}}]=d$, $MA_iQ_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}M^{\mathrm{T}}\!=\!MA_iQ_{{\hat{x}}{\hat{x}}}-M{\bar{Q}}_{ii}M^{\mathrm{T}}$, and $Q_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}M^{\mathrm{T}}\!=\!M[{\bar{Q}}_{ii}+A_iQ_{{\hat{x}}{\hat{x}}}A_i^{\mathrm{T}}]M^{\mathrm{T}}$. $\square$

Proof Theorem 1a

$\mathrm{MSEM}(\hat{{\underline{x}}})\!=\!Q_{{\hat{x}}{\hat{x}}}$ and $\mathrm{MSEM}(\hat{{\underline{x}}}_{w})\!=\!Q_{{\hat{x}}{\hat{x}}}+M(dd^{\mathrm{T}}\!-2{\bar{Q}}_{ii}\!-A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}} )M^{\mathrm{T}}$. Thus, $\mathrm{MSEM}(\hat{{\underline{x}}}_{w}) \le \mathrm{MSEM}(\hat{{\underline{x}}})$ if and only if $dd^{\mathrm{T}}-2{\bar{Q}}_{ii}-A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}} \le 0$ or $d^{\mathrm{T}}[2 {\bar{Q}}_{ii}+A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}]^{-1}d \le 1$. The last step follows from the equivalence: $dd^{\mathrm{T}}\le N \Leftrightarrow d^{\mathrm{T}}N^{-1}d \le 1$ for positive-definite N. ($\Rightarrow$) $dd^{\mathrm{T}} \le N$ means $f^{\mathrm{T}}(dd^{\mathrm{T}}-N)f\le 0,~\forall f\in {\mathbb {R}}^p$. Substitution of $f\!=\!N^{-1}d$ gives $(d^{\mathrm{T}}N^{-1}d)^2 - (d^{\mathrm{T}}N^{-1}d)\le 0$ or $d^{\mathrm{T}}N^{-1}d \le 1$. ($\Leftarrow$) Every $f\in {\mathbb {R}}^p$ can be expressed as $f = [d^{+T},d^{\bot }]q$ for some vector $q\in {\mathbb {R}}^{p}$, as $[d^{+T}, d^\bot ]$ is a square and invertible matrix ($d^+=(d^{\mathrm{T}}N^{-1}d)^{-1}d^{\mathrm{T}}N^{-1}$, $d^{\mathrm{T}}d^\bot =0$). This gives $f^{\mathrm{T}}\!(dd^{\mathrm{T}}-N)\!f = -q^{\mathrm{T}}\mathrm{blkdiag}\left( (d^{\mathrm{T}}N^{-1}d)^{-1}\!-\!1,d^{\bot T}Nd^\bot \right) \!q \le 0$, as $(d^{\mathrm{T}}N^{-1}d)^{-1} \ge 1$ and $d^{\bot T}Nd^\bot \ge 0$. $\square$

Proof of Theorem 1b

The $\mathrm{MSEM}$ of $\underline{{\hat{\theta }}}_{w}=F^{\mathrm{T}}{\hat{x}}_{w}$ follows from $\mathrm{MSEM}(\underline{{\hat{x}}}_{w})=Q_{{\hat{x}}{\hat{x}}}+M(dd^{\mathrm{T}}-2{\bar{Q}}_{ii}-A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}})M^{\mathrm{T}}$ as $\mathrm{MSEM}(\underline{{\hat{\theta }}}_{w})=Q_{{\hat{\theta }}{\hat{\theta }}}+F^{\mathrm{T}}M(dd^{\mathrm{T}}-2{\bar{Q}}_{ii}-A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}})M^{\mathrm{T}}F$. This gives for ${\mathsf {E}}(||\underline{{\hat{\theta }}}_{w}-\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})=\mathrm{trace}[\mathrm{MSEM}(\underline{{\hat{\theta }}}_{w}) \mathrm{MSEM}(\underline{{\hat{\theta }}})^{-1}]$,

$$\begin{aligned} \begin{array}{l} {\mathsf {E}}(||\underline{{\hat{\theta }}}_{w}-\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})= {\mathsf {E}}(||\underline{{\hat{\theta }}}-\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})+\\ +\mathrm{trace}\{[dd^{\mathrm{T}}-2{\bar{Q}}_{ii}-Q_{{\hat{i}}{\hat{i}}}]X\} \end{array} \end{aligned}$$

(28)

with $X=[{\bar{Q}}_{ii}+Q_{{\hat{i}}{\hat{i}}}]^{-1}D_{{\hat{i}}{\hat{i}}}[{\bar{Q}}_{ii}+Q_{{\hat{i}}{\hat{i}}}]^{-1}$, $Q_{{\hat{i}}{\hat{i}}}=A_{i}Q_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}$, $D_{{\hat{x}}{\hat{x}}}=Q_{{\hat{x}}{\hat{x}}}F(F^{\mathrm{T}}Q_{{\hat{x}}{\hat{x}}}F)^{-1}F^{\mathrm{T}}Q_{{\hat{x}}{\hat{x}}}=Q_{{\hat{x}}{\hat{x}}}-Q_{{\hat{x}}{\hat{x}}|F^{\mathrm{T}}x}$, $D_{{\hat{i}}{\hat{i}}}=A_{i}D_{{\hat{x}}{\hat{x}}}A_{i}^{\mathrm{T}}$. Hence ${\mathsf {E}}(||\underline{{\hat{\theta }}}_{w}-\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})\le {\mathsf {E}}(||\underline{{\hat{\theta }}}-\theta ||_{Q_{{\hat{\theta }}{\hat{\theta }}}}^{2})$ if and only if $\mathrm{trace}\{[dd^{\mathrm{T}}-2{\bar{Q}}_{ii}-Q_{{\hat{i}}{\hat{i}}}]X\} \le 0$, from which the result follows. $\square$

Proof of Lemma

The normal matrix of model (14) is given as

$$\begin{aligned} N=\left[ \begin{array}{lcl} k A^{\mathrm{T}}Q^{-1}A &{} \;e_{k}^{\mathrm{T}} \otimes A^{\mathrm{T}}Q^{-1}B \\ e_{k} \otimes B^{\mathrm{T}}Q^{-1}A &{} \;I_{k} \otimes B^{\mathrm{T}}Q^{-1}B \end{array} \right] \end{aligned}$$

with $Q=\mathrm{blockdiagonal}(Q_{\phi }, Q_{p}) \otimes R$ and for ‘ambiguity-float,’

$$\begin{aligned} A= \left[ \begin{array}{cc} e_{f} \otimes {\tilde{G}} &{} \;\Lambda \otimes I_{m-1}\\ e_{f} \otimes {\tilde{G}} &{} \;0 \end{array} \right] \;,\; B= \left[ \begin{array}{c} -\mu \otimes I_{m-1}\\ +\mu \otimes I_{m-1} \end{array} \right] . \end{aligned}$$

Hence, the reduced normal matrix for the ionospheric delays follows, with $P_{A}=A(A^{\mathrm{T}}Q^{-1}A)^{-1}A^{\mathrm{T}}Q^{-1}$ and $P_{A}^{\perp }=I-P_{A}$, as

$$\begin{aligned} \begin{array}{lcl} Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}} ^{-1} &{}=&{} [I_{k} \otimes B^{\mathrm{T}}Q^{-1}B-P_{k}\otimes B^{\mathrm{T}}Q^{-1}P_{A}B]\\ &{}=&{} [P_{k} \otimes B^{\mathrm{T}}Q^{-1}P_{A}^{\perp }B+P_{k}^{\perp }\otimes B^{\mathrm{T}}Q^{-1}B] \end{array} \end{aligned}$$

which upon inversion gives

$$\begin{aligned} Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}}= P_{k}\otimes [B^{\mathrm{T}}Q^{-1}P_{A}^{\perp }B]^{-1}+P_{k}^{\perp }\otimes [B^{\mathrm{T}}Q^{-1}B]^{-1} \end{aligned}$$

with

$$\begin{aligned} \begin{array}{lcl} [B^{\mathrm{T}}Q^{-1}P_{A}^{\perp }B]^{-1} &{}=&{} [\sigma _{{\hat{i}}}^{-2}R^{-1}P_{{\tilde{G}}}+\sigma _{{\hat{i}}|\rho }^{-2}R^{-1}P_{{\tilde{G}}}^{\perp }]^{-1}\\ &{}=&{} [\sigma _{{\hat{i}}}^{2}P_{{\tilde{G}}}R+\sigma _{{\hat{i}}|\rho }^{2}P_{{\tilde{G}}}^{\perp }R] \end{array} \end{aligned}$$

and

$$\begin{aligned}{}[B^{\mathrm{T}}Q^{-1}B]^{-1} = \sigma _{{\check{i}}|\rho }^{2}R. \end{aligned}$$

Note that in the above use was made of the following property: If P and $P^{\perp }$ are two complementary projectors (i.e., $P+P^{\perp }=I$), and S and T are invertible, then $[P\otimes S + P^{\perp } \otimes T]^{-1}=P \otimes S^{-1} + P^{\perp } \otimes T^{-1}$. The proof for the ‘baseline-constrained’ and ‘ambiguity-fixed’ cases, go along similar lines. $\square$

Proof of Theorem 2

With the choice ${\bar{Q}}_{ii}=\lambda Q_{{\hat{i}}{\hat{i}}}$, the MSE-expression of (28) becomes

$$\begin{aligned} \begin{array}{lcl} \mathrm{MSE}(\lambda ) &{} = &{} \mathrm{MSE}(\infty )+ \frac{\mathrm{trace}[(2\lambda +1)Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}}]}{(\lambda +1)^{2}} \; \times \quad \\ &{} &{} \qquad \left( \frac{d^{\mathrm{T}}[Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}}Q_{{\hat{i}}{\hat{i}}}^{-1}]d}{\mathrm{trace}[(2\lambda +1)Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}}]}-1\right) . \end{array} \end{aligned}$$

(29)

From setting the derivative to zero, $\tfrac{d\mathrm{MSE}}{d \lambda }=0$, and checking whether the second-order derivative is positive at this stationary point, the minimizer $\lambda _{\mathrm{min}}$ follows. With this minimizer, it then follows from (29) that $\mathrm{MSE}(\lambda ) \le \mathrm{MSE}(\infty )$ for $\lambda \ge \tfrac{1}{2}(\lambda _{\mathrm{min}}-1)$. $\square$

Proof Corollary 1a

We only prove the expression for $\lambda _{\mathrm{min}}^{(1)}$ as the proof for $\lambda _{\mathrm{min}}^{(2)}$ and $\lambda _{\mathrm{min}}^{(3)}$ goes along similar lines. For $F=I$, we have $D_{{\hat{i}}{\hat{i}}}=Q_{{\hat{i}}{\hat{i}}}$, since $Q_{{\hat{i}}{\hat{i}}|F^{\mathrm{T}}x}=0$. With the choice $A_{i}x \rightarrow (D_{k}^{\mathrm{T}} \otimes I_{m-1})i_{[k]}$, it follows that $Q_{{\hat{i}}{\hat{i}}}=A_{i}Q_{{\hat{i}}_{[k]}{\hat{i}}_{[k]}}A_{i}^{\mathrm{T}}=D_{k}^{\mathrm{T}}D_{k} \otimes \sigma _{{\check{i}}|\rho }^{2}R$ and therefore $\mathrm{trace}(Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}})=(k-1)(m-1)$. Substitution into (20) gives, with $P_{k}^{\perp }=D_{k}(D_{k}^{\mathrm{T}}D_{k})^{-1}D_{k}^{\mathrm{T}}=I_{k}-\tfrac{1}{k}e_{k}e_{k}^{\mathrm{T}}$,

$$\begin{aligned} \lambda _{\mathrm{min}}^{(1)} = \frac{i_{[k]}^{\mathrm{T}}[P_{k}^{\perp } \otimes \sigma _{{\check{i}}|\rho }^{-2}R^{-1}]i_{[k]}}{(k-1)(m-1)}. \end{aligned}$$

(30)

Writing $i_{[k]}=\mathrm{vec}({\mathcal {I}}_{k})$, with ${\mathcal {I}}_{k}=[i_{1}, \ldots , i_{k}]$, and making use of the property $\mathrm{trace}(A^{\mathrm{T}}BCD^{\mathrm{T}})=\mathrm{vec}(A)^{\mathrm{T}}(D \otimes B)\mathrm{vec}(C)$, gives

$$\begin{aligned} \begin{array}{lcl} i_{[k]}^{\mathrm{T}} (P_{k}^{\perp } \otimes R^{-1})i_{[k]} &{}=&{} \mathrm{trace}[{\mathcal {I}}_{k}R^{-1}{\mathcal {I}}_{k}P_{k}^{\perp }]\\ &{}=&{} \sum \limits _{t=1}^{k}(i_{t}-i_{{\bar{t}}})^{\mathrm{T}}R^{-1}(i_{t}-i_{{\bar{t}}}) \end{array} \end{aligned}$$

which, combined with (30), completes the proof. $\square$

Proof Corollary 1b

We only prove the expression for $\lambda _{\mathrm{min}}^{(b)}$ as the proof for $\lambda _{\mathrm{min}}^{(a)}$ goes along similar lines. As F is now chosen to only select the baseline and $A_{i}x \rightarrow I_{k} \otimes I_{m-1}$, we have $D_{{\hat{i}}{\hat{i}}}=Q_{{\hat{i}}{\hat{i}}}-Q_{{\hat{i}}{\hat{i}}|b}=P_{k} \otimes (\alpha -\beta )P_{{\tilde{G}}}R$. Therefore,

$$\begin{aligned} \begin{array}{lcl} Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}} &{}=&{} P_{k} \otimes \tfrac{1}{\alpha }(\alpha -\beta ) R^{-1}P_{G}R\\ \mathrm{trace}(Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}})&{}=&{} (1-\tfrac{\beta }{\alpha })\mathrm{trace}(P_{{\tilde{G}}})\\ Q_{{\hat{i}}{\hat{i}}}^{-1}D_{{\hat{i}}{\hat{i}}}Q_{{\hat{i}}{\hat{i}}}^{-1}&{}=&{} P_{k} \otimes \tfrac{1}{\alpha ^{2}}(\alpha -\beta )R^{-1}P_{{\tilde{G}}} \end{array} \end{aligned}$$

from which, after substitution into (20), the result follows. $\square$.