**1. Introduction**

Global Navigation Satellite systems (GNSS) can provide accurate position and velocity information in outdoor environments, and its errors do not accumulate over time [1]. The disadvantages are that it can only provide less accuracy attitude information, the output frequency is low (1–20 Hz), and it is vulnerable to environmental interference. In contrast, the Inertial Navigation System (INS) is less dependent on the environment, and relies entirely on the angular velocity and acceleration information that is measured by the Inertial Measurement Unit (IMU), which can provide high-frequency navigation information [2,3]. But the position error will disperse over time due to the integral acquisition of positional information, resulting large errors in navigation results. Therefore, combining the respective advantages of GNSS and INS to obtain the navigation results with high accuracy, high interference immunity, and high frequency is a hot topic of research in the field of navigation at present [4,5].

The Kalman filter (KF) and its upgrade variants are the most widely utilized algorithms for INS and GNSS information fusion [6,7]. The traditional Kalman filter algorithm can only be applied to linear systems, but most of the information in real navigation systems are nonlinear. Bucy et al. [8] proposed the extended Kalman filter (EKF), which linearizes the nonlinear function around the current estimate, and truncates the first-order linearization of the Taylor expansion of the nonlinear function. The remaining higher-order terms are ignored, and their performance depends on the degree of local nonlinearity. The unscented Kalman filter (UKF) was proposed to further improve the performance under nonlinear systems by making the nonlinear system equations applicable to linear assumptions through lossless transformations [9,10]. By approximating the posterior probability density of the state with a series of deterministic samples, the problem of the EKF

**Citation:** Zhao, S.; Zhou, Y.; Huang, T. A Novel Method for AI-Assisted INS/GNSS Navigation System Based on CNN-GRU and CKF during GNSS Outage. *Remote Sens.* **2022**, *14*, 4494. https://doi.org/10.3390/rs14184494

Academic Editors: Yuwei Chen, Changhui Jiang, Qian Meng, Bing Xu, Wang Gao, Panlong Wu, Lianwu Guan and Zeyu Li

Received: 31 July 2022 Accepted: 6 September 2022 Published: 9 September 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

accuracy dispersion under a highly nonlinear system is avoided. But the UKF has low accuracy in the high-dimensional case, the cubature Kalman filter (CKF) that is based on the spherical radial volume criterion is applied to data fusion, which can effectively approximate the Gaussian density function with higher accuracy, convenient parameter selection, and good convergence effect [11–13]. In order to improve the fusion accuracy in complex measurement environments, robust Kalman algorithms have also started to attract the attention of researchers [14,15]. To solve the problem of error model that is caused by measurement anomalies, Chen et al. [16] proposed a new cardinal maximum correlation entropy Kalman filter, which uses the robust maximum correlation entropy criterion (MCC) as the optimality criterion to solve the state estimation problem under outlier interference by maximizing the correlation entropy between states and measurements. Yun et al. [17] proposed a variational Bayesian-based state estimation algorithm to improve the CKF accuracy under dynamic model mismatch and outlier disturbance.

When the GNSS signals are unavailable, KF operates in predictive mode and corrects INS measurements according to the system error model. At this time, the accuracy of data fusion that relies only on the KF is not effective and navigation performance deteriorates rapidly. To improve the integrated navigation accuracy during GNSS outage, machine learning has started to be applied to integrated navigation systems. Ning et al. [18] proposed an optimal radial basis function (RBF)-based neural network that can improve the overall positioning accuracy during short-term GNSS signal outages. Hang et al. [19] proposed a new hybrid intelligence algorithm combining a discrete gray predictor (DGP) and a multilayer perceptron (MLP) neural network that provides pseudo-GPS positions during GNSS failures and uses GNSS position information from the last few moments to predict positions for future moments. Compared with traditional artificial neural networks, recurrent neural networks are more advantageous in combined navigation systems and can make full use of historical information [20–22]. Liu et al. [23] proposed a multi-channel long-short term memory (LSTM) network to predict the increments of vehicle position, which reduces the navigation error in case of GNSS outages by an order of magnitude. In practical applications, a large amount of historical data before the GNSS outage needs to be trained when the GNSS outage occurs, so the training efficiency of neural networks also has high requirements. Tang et al. [24] proposed a hybrid algorithm that was based on the gated recurrent unit (GRU) and adaptive Kalman filter (AKF), and the experimental results showed that GRU outperformed LSTM in terms of prediction accuracy and training efficiency. Zhi et al. [25] proposed a convolutional neural network-long short-term memory (CNN-LSTM) model, which uses convolutional neural network (CNN) to quickly extract the features of the input and LSTM network to output the pseudo-GPS signal, further improving the training efficiency. However, most of the current articles use the offline simulation, assuming that the GNSS failure time is known and do not consider the time that is required to train the model online. Al Bitar et al. [26] proposed a novel real-time training strategy for regular training on the past one minute data, with the disadvantage that only short historical data are used and the accuracy is poor when the time of GNSS outage is long.

To overcome the shortcomings of the traditional methods, our paper proposes a new AI-assisted method for the integrated high-precision INS/GNSS navigation system. The method consists of two parts: first, CKF is used to provide more accurate neural network training samples. Then, by building a CNN-GRU network to predict the position increments during GNSS outage, the CNN is utilized to quickly extract the multi-dimensional sequence features, and GRU is used to model the time series. In addition, a new real-time training strategy is proposed for practical application scenarios, where the duration of the GNSS outage time and the motion state information of the vehicle are taken into account in the training strategy. The experiments verify that the proposed algorithm has the advantages of high prediction accuracy and high training efficiency.

The rest of the paper is organized as follows: Section 2 introduces the INS error propagation model and the integrated navigation model that is based on CKF, Section 3 introduces the proposed CNN-GRU network, Section 4 performs the road test and result analysis, and the conclusion is presented in Section 5.

#### **2. Mathematical Integrated System Model**

The INS and GNSS are loosely coupled as shown in Figure 1. The INS and GNSS can complete position and velocity independently. Then, the position and velocity errors of INS are estimated using CKF, and these errors are used to correct the navigation output of INS and achieve the correction of gyro and accelerometer drift, thus reducing the impact of INS errors. In this section, we derive the error propagation model of INS, and then introduce the integrated navigation model that is based on CKF.

**Figure 1.** Flowchart of the loosely coupled integrated system.

#### *2.1. The Error Propagation Model*

By integrating the angular velocity that is measured by gyro, the attitude direction of the inertial component can be obtained. By using the directional cosine matrix of the attitude direction, the specific force components of the accelerometer observations along each axis of the carrier system can be converted to the navigation coordinate system, the velocity and position can be calculated.

In the process of attitude calculation, the navigation coordinate system that is obtained is regarded as the real navigation coordinate system. In practice, due to the interference of various factors, the calculated navigation coordinate system will have deviations compared to the real navigation coordinate system, and the attitude error equation can be expressed as:

$$
\dot{\phi} = \phi \times \omega\_{in}^{\text{n}} + \delta \omega\_{in}^{\text{n}} - \delta \omega\_{ib'}^{\text{n}} \tag{1}
$$

where *δω<sup>n</sup> in* = *δω<sup>n</sup> ie* + *δω<sup>n</sup> en* is the angular velocity error in the navigation coordinate system, *δω<sup>n</sup> ie* is the angular velocity error of the Earth's rotation, *δω<sup>n</sup> en* is the rotation error of the navigation system, *δω<sup>n</sup> ib* = *<sup>C</sup><sup>n</sup> <sup>b</sup> δω<sup>b</sup> ib* = *<sup>C</sup><sup>n</sup> b* (*δKG* + *δG*)*ω<sup>b</sup> ib* + *<sup>ε</sup><sup>b</sup>* is the gyro measurement error, *ε<sup>b</sup>* is the gyro drift error, and *δKG* and *δG* are the gyro scale factor error and nonorthogonal installation error, which can be expressed as:

$$
\begin{aligned}
\delta \mathcal{K}\_{\mathcal{G}} &= \text{diag}\left( \begin{bmatrix}
\ \delta \mathcal{K}\_{\mathcal{G}x} & \delta \mathcal{K}\_{\mathcal{G}y} & \delta \mathcal{K}\_{\mathcal{G}z} \\
0 & \delta \mathcal{G}\_z & -\delta \mathcal{G}\_y \\
\delta \mathcal{G}\_y & -\delta \mathcal{G}\_x & 0
\end{bmatrix} \right) \tag{2}
$$

Both calibration residuals and installation errors can be considered as constants, while random drift can be expressed as a cumulative model of random constants and a first-order Markov model. In the presence of attitude error and specific force measurement error, the velocity differential error can be expressed as:

$$\delta \dot{\boldsymbol{v}}^{\boldsymbol{n}} = f\_{sf}^{\boldsymbol{n}} \times \boldsymbol{\phi} + \left(2\delta\omega\_{\rm ir}^{\boldsymbol{n}} + \delta\omega\_{\rm en}^{\boldsymbol{n}}\right) \times \boldsymbol{v}^{\boldsymbol{n}} - \left(2\omega\_{\rm ir}^{\boldsymbol{n}} + \omega\_{\rm en}^{\boldsymbol{n}}\right) \times \delta \boldsymbol{v}^{\boldsymbol{n}} + \delta f\_{sf}^{\boldsymbol{n}} + \delta \mathbf{g}^{\boldsymbol{n}} \tag{3}$$

where *δ f <sup>n</sup> s f* is the accelerometer measurement error, which can be expressed as the cumulative model of accelerometer zero deviation and white noise.

$$\delta \mathbf{g}^n = \begin{bmatrix} 0 \\ -\beta\_3 (\sin 2L \cdot \delta \mathbf{h} + 2h \cos 2L \cdot \delta L) \\ -[\mathbf{g}\_\varepsilon \sin 2L(\beta - 4\beta\_1 \cos 2L) \delta L - \beta\_2 \delta h] \end{bmatrix} \tag{4}$$

where *<sup>β</sup>* is the gravity flattening, *<sup>β</sup>*<sup>1</sup> = 2.32718 × <sup>10</sup>−5*s*<sup>−</sup>2, *<sup>β</sup>*<sup>2</sup> = 3.08 × <sup>10</sup>−6*s*<sup>−</sup>2, *<sup>β</sup>*<sup>3</sup> = 8.08 × <sup>10</sup>−9*s*<sup>−</sup>2.

The position error is [27]:

$$\begin{cases} \delta \dot{L} = \frac{1}{\overline{R}\_M + h} \delta v\_N - \frac{v\_N}{\left(\overline{R}\_M + h\right)^2} \delta h\\ \delta \dot{\lambda} = \frac{\sec L}{\overline{R}\_N + h} \delta v\_E + \frac{v\_F \sec L \tan L}{\overline{R}\_N + h} \delta L - \frac{v\_F \sec L}{\left(\overline{R}\_N + h\right)^2} \delta h\\ \delta \dot{h} = \delta v\_{II} \end{cases} \tag{5}$$

#### *2.2. The Integrated Navigation Model Based on CKF*

In the integrated navigation system, the CKF estimation is adopted to estimate the system state vector which is a 15 array vector:

$$dX = \begin{bmatrix} \mathcal{Q}\_{\mathcal{E}} \ \mathcal{q}\_{\mathcal{U}} \ \mathcal{q}\_{\mathcal{U}} \ \delta V\_{\mathcal{E}} \ \delta V\_{\mathcal{U}} \ \delta L \ \delta \lambda \ \delta h \ \nabla\_{\mathcal{X}} \ \nabla\_{\mathcal{Y}} \ \nabla\_{\mathcal{Z}} \ \varepsilon\_{\mathcal{X}} \ \varepsilon\_{\mathcal{Y}} \ \varepsilon\_{\mathcal{Z}} \end{bmatrix}^{T} \tag{6}$$

where *ϕe*, *ϕn*, and *ϕ<sup>u</sup>* denote the attitude error angles of the INS in the east, north, and zenith directions. *δVe*, *δVn*, and *δV*<sup>u</sup> denote the velocity error of INS in the three directions. *δL*, *δλ*, and *δh* denote the latitude error, longitude error, and altitude error. ∇*x*, ∇*y*, and ∇*<sup>z</sup>* are the accelerometer bias vectors, and *εx*, *εy*, and *ε<sup>z</sup>* are the gyro bias vectors.

The equation of state and measurement equation of the system can be expressed as:

$$\begin{cases} \begin{aligned} \dot{X} &= FX + GW \\ Z &= HX + V \end{aligned} \tag{7}$$

where *F* is the system state transfer matrix, which consists of the INS system error. *G* is the system noise matrix. *Z* is the measurement vector, which is set to the position error vector and velocity error vector of INS and GNSS. *H* is the observation matrix. *W* is the process noise vector. *V* denotes the measured noise vector.

The CKF is based on the 3rd order volume rule and is suitable for high-dimensional filtering problems. The computation consists of two steps: the time update and measurement update phase, as shown in Algorithm 1.

#### **Algorithm 1**: Cubature Kalman Filter.

**Require: <sup>x</sup>**<sup>ˆ</sup> *<sup>k</sup>*−1|*k*−1, **<sup>P</sup>***k*−1|*k*−1, **<sup>Q</sup>***k*−1, **<sup>R</sup>***<sup>k</sup>* **for** i = 1, . . . , N **do x***i <sup>k</sup>*−<sup>1</sup> <sup>=</sup> *Sk*−1|*k*−1*ξ<sup>i</sup>* <sup>+</sup> **<sup>x</sup>**<sup>ˆ</sup> *<sup>k</sup>*−1|*k*−<sup>1</sup> **w***i <sup>k</sup>*−<sup>1</sup> <sup>=</sup> 1/*<sup>N</sup>* **end for Prediction phase: ^ <sup>x</sup>***k*|*k*−<sup>1</sup> <sup>=</sup> *<sup>N</sup>* ∑ *i*=1 **w***i <sup>k</sup>*−1<sup>f</sup> **x***i k*−1 **<sup>P</sup>***k*|*k*−<sup>1</sup> <sup>=</sup> *<sup>N</sup>* ∑ *i*=1 **w***i k*−1 f **x***i k*−1 − **^ x***k*|*k*−<sup>1</sup> f **x***i k*−1 − **^ x***k*|*k*−<sup>1</sup> *T* + **<sup>Q</sup>***k*−<sup>1</sup> **for i** = 1, . . . , N **do x***i <sup>k</sup>*|*k*−<sup>1</sup> <sup>=</sup> *Sk*|*k*−1*ξ<sup>i</sup>* <sup>+</sup> **<sup>x</sup>**<sup>ˆ</sup> *<sup>k</sup>*|*k*−<sup>1</sup> **w***i <sup>k</sup>*|*k*−<sup>1</sup> <sup>=</sup> 1/*<sup>N</sup>* **end for Update phase: <sup>z</sup>**<sup>ˆ</sup> *<sup>k</sup>*|*k*−<sup>1</sup> <sup>=</sup> *<sup>N</sup>* ∑ *i*=1 **w**i *<sup>k</sup>*|*k*−<sup>1</sup> <sup>h</sup> **x***i k*|*k*−1 **P***zz <sup>k</sup>*|*k*−<sup>1</sup> <sup>=</sup> *<sup>N</sup>* ∑ *i*=1 **w***i k*|*k*−1 h **x***i k*|*k*−1 − **<sup>z</sup>**<sup>ˆ</sup> *<sup>k</sup>*|*k*−<sup>1</sup> <sup>h</sup> **x***i k*|*k*−1 − **<sup>z</sup>**<sup>ˆ</sup> *<sup>k</sup>*|*k*−<sup>1</sup> *T* + **R***<sup>k</sup>* **P***x*<sup>z</sup> *<sup>k</sup>*|*k*−<sup>1</sup> <sup>=</sup> *<sup>N</sup>* ∑ *i*=1 **w***i k*|*k*−1 **x***i <sup>k</sup>*|*k*−<sup>1</sup> <sup>−</sup> **^ x***k*|*k*−<sup>1</sup> h **x***i k*|*k*−1 − **<sup>z</sup>**<sup>ˆ</sup> *<sup>k</sup>*|*k*−<sup>1</sup> *T*
