Robust Underwater Direction-of-Arrival Tracking Based on AI-Aided Variational Bayesian Extended Kalman Filter

Abstract

The AI-aided variational Bayesian extended Kalman filter (AI-VBEKF)-based robust direction-of-arrival (DOA) technique is proposed to make reliable estimations of the bearing angle of an uncooperative underwater target with uncertain environment noise. Considering that the large error of the guess of the initial mean square error matrix (MSEM) will lead to inaccurate DOA tracking results, an attention-based deep convolutional neural network is first proposed to make reliable estimations of the initial MSEM. Then, by utilizing the AI-VBEKF estimating scheme, the uncertain measurement noise caused by the unknown underwater environment along with the bearing angle of the target can be estimated simultaneously to provide reliable results at every DOA tracking step. The proposed technique is demonstrated and verified by both of the simulations and the real sea trial data from the South China Sea in July 2021, and both the robustness and accuracy are proven superior to the traditional DOA-estimating methods.

1. Introduction

Estimating the bearing angle of an uncooperative underwater target is vitally important in underwater acoustic engineering [1,2,3,4]. Considering that the unknown underwater environment sometimes causes uncertain measurement noise that will degrade the estimation accuracy of the interested target’s bearing angle, the robust bearing angle estimating technique is necessary and has attracted a number of researchers’ attention [5,6,7,8]. In addition, the techniques of passive bearing angle estimation can be based on two different ways, namely the DOA estimation and the DOA tracking methods. The DOA estimation merely depends on the measurements from a certain period. Besides the measurements, the DOA tracking techniques also utilize the prior motion information of an interested target. Although the DOA estimation has already made many milestones in various bearing angle estimation missions [9,10,11], its ignorance of the target’s kinematics makes the DOA estimation techniques always have superior estimating results when the target is nearly static and the environment is stable. When the target is maneuvering or the measurement noise shifts because the underwater environment is uncertain, the DOA estimation methods can hardly have satisfied results. In addition, in order to make reliable and accurate estimations of an uncooperative underwater target’s bearing angle, the DOA estimation techniques usually need a wide observation window, which will lead to heavy on-board computational load and which will occupy a lot of on-board limited storage.

Considering the advantages and drawbacks of the DOA estimation methods, DOA tracking techniques have been proposed. Instead of merely utilizing the measurement information, DOA tracking methods take both the target’s motion model and the current measurements into consideration. Based on the Bayesian estimation principle, the DOA tracking techniques can reach higher estimating accuracy and can occupy less on-board computing or storage resources, especially when the underwater target is maneuvering [12,13,14,15]. In addition, the unpredicted moving objects or environmental changes in the ocean can sometimes make the measurements uncertain, which will degrade the bearing angle estimation accuracy or even make the estimating procedure diverge. Under this circumstance, the robust DOA tracking technique is necessary to guarantee the precision of the whole estimation process. Among all the robust DOA tracking techniques, two categories can be divided depending on the basic mathematical principles, namely the data-inspiring robust tracking techniques [16] and the Bayesian inferring techniques [17,18]. These two categories of robust DOA tracking techniques are both derived based on the widely used Kalman filter technique [19,20,21]. The former one utilizes good estimations from a long tracking period to adjust the weight of the bad measurements when the measurement noise becomes uncertain. As a result, the weights of the poor measurements can become very small to make the current bearing angle estimations mainly depend on the historical good estimations. These kinds of data-inspiring robust DOA tracking techniques have proven their efficiency in a number of research studies [22,23,24,25,26] nevertheless, their convergence is still an open question. When the hyper-parameters of the inspired robust DOA tracking are set inappropriately, the estimations usually have poor accuracy and sometimes diverge.

To overcome the drawbacks of the data-inspiring robust tracking techniques, the Bayesian inferring techniques, or the so-called variational Bayesian (VB) robust tracking techniques, are proposed and have been utilized in a various of research studies [27,28,29]. Unlike the data-inspiring robust tracking methods, the VB robust tracking techniques are strictly derived mathematically. By assuming the prior probability density functions (PDFs) of the tracking parameters and the interested time-varying parameters (i.e., covariance matrix of the uncertain measurements noise) along with optimizing the Kullback–Leibler Divergence (KLD) among the posterior PDFs, the VB technique can analytically estimate the tracking parameters and other interested parameters simultaneously. Different from the data-inspiring robust tracking methods, the VB robust tracking techniques have strict mathematical foundations and can guarantee convergence and accuracy if the prior PDFs are set properly. Under this circumstance, the VB robust tracking techniques do not need good estimations from a long period to output reliable results when measurements become uncertain, and they are more theoretically complete.

In addition, the error between the guess of the initial state and its real value will influence the final DOA tracking accuracy since the DOA tracking system has high nonlinearities. If the error is large enough, the tracking results will even be diverged [30]. Under this circumstance, before operating a DOA tracking algorithm, the initial MSEM should be set properly. However, in the real DOA tracking scenario, it is almost impossible to obtain an accurate error covariance matrix since the real initial states cannot be accurately obtained. In addition, the uncertain underwater environment always affects the state estimating results and makes the initial MSEM much harder to accurately predetermine. To deal with this problem, some researchers suggested utilizing a large guess of the initial MSEM to make the tracking system have fast converging speed and to deal with the large error between the initial guesses of the states and their real values. However, since the underwater target DOA tracking system has high nonlinearities, setting a large guess of the initial MSEM sometimes makes the whole tracking procedure become unstable or even diverge. As a consequence, how to set a proper guess of the initial MSEM with respect to the error of the initial guess of the state and its real value is still an open question. Usually, the initial MSEM is set by engineering experience, which is similar to an expert system that tries to have common solutions. Although some adaptive tracking techniques have been designed to weaken the influences by the initial error [31], the unproperly set initial MSEM still affects the overall tracking performance. Under these circumstances, a robust and accurate technique for properly setting the guess of the initial MSEM is necessary. Usually, the traditional Bayesian estimating techniques rely heavily on the mathematical model of the systems or parameters to be estimated, but the guess of the initial MSEM does not have a determined and accurate mathematical model with respect to the system state, which makes the traditional methods fail. However, the value of the guess of the initial MSEM has a relationship with the theoretical measurements calculated by the state by the measurement model, initial guess of the state, and its real value. Under this circumstance, the data characteristics other than the mathematical characteristics of the guess of the initial MSEM can be analyzed and utilized. As a result, although the traditional techniques have few general solutions for the problem of properly setting the guess of the initial MSEM, the latest data-driven AI techniques allow the problem discussed above to be solved.

Except for depending heavily on the mathematical derivations, the AI techniques can only rely on the statistical characteristics to infer the inner nonlinear relationship between the inputs and the outputs from a certain system. As a result, the AI techniques are more suitable for the problems that are hard to be mathematically or analytically depicted but have plenty of data resources. After the back propagated (BP) training method was proposed by Rumelhart et al. [32] and after the deep learning scheme that brought the deep artificial neural network (DANN) into light by Hinton et al. [33], the DANN has shown its powerful abilities in many fields, such as automation, parameter estimation, target tracking, etc. Among all of the DANNs, the deep convolutional neural networks (DCNNs) and the variations of the DCNNs (i.e., Resnet [34], Senet [35]) are the most notable ones, with robust and accurate performances, especially for regression and classification problems that are hard to be analytically modeled. By adopting the convolutional and pooling ideas into one DANN, the DCNN can be modeled much deeper and can easily be trained. Based on the DCNN technique, Lecun et al. [36] proposed a DCNN named Lennet-5 to deal with the digital image classification task. Ref. [37] considered the training efficiency and the gradient vanishing problem caused by the DCNNs and proposed the shortcut concept to make the DCNNs have a much deeper structure so that the DCNNs can deal with more complicated classification and regression issues. Ref. [38] considered that the different features may have different weights during the training process of the designed DCNNs and proposed an attention-based scheme to make the designed DCNNs more robust and accurate. Besides the principle developments of the DCNNs, the applications of them have also attracted a number of researchers from various fields, especially from the field of underwater target localization, detection, classification, and tracking [37,38,39,40,41,42,43]. Niu et al. [37,38,39,40] have published a series of research studies to make the underwater target localization task more accurate and robust to the uncertain underwater environment. Wang et al. [41,42] utilized the DCNNs technique in estimating the unknown parameters in the underwater acoustical channel. Based on the DCNNs technique, Ref. [40] developed intelligent target classification algorithms that can obtain higher accuracy and classification speed than the traditional methods. Ref. [43] proposed a DCNN-based seabed parameters inversion method. However, to the best of the authors’ knowledge, very few research studies have considered using the advanced DCNN technique for the underwater target DOA tracking problem. As a result, the potential and strength of the DCNNs in the underwater target DOA tracking scenario, especially in dealing with the proper setting of the guess of the initial MSEM, have not been proven yet.

Based on the above analysis, an AI-aided variational Bayesian extended Kalman filter (AI-VBEKF)-based robust direction-of-arrival (DOA) technique is proposed to make reliable estimations of the bearing angle of an uncooperative underwater target with uncertain environment noise. The main contributions of this study are summarized as follows.

Firstly, a uniform circular array (UCA) is considered to provide measurements of the underwater target. By adopting the UCA as the measurement system, the port and starboard ambiguity problem is overcome. In addition, the uniform aperture at all bearing angles is processed.

Secondly, considering the effects caused by the unproperly set guess of the initial MSEM, an attention-based DCNN is designed to make reliable initial guesses of it and to make the later DOA tracking process steady and accurate.

Thirdly, considering that the unknown underwater disturbance can sometimes make the measurement noise uncertain, the AI-VBEKF is designed to robustly estimate the bearing angle of the underwater target with a shifting covariance matrix of the uncertain measurement noise.

Finally, based on the sea trail data from the South China Sea in July 2021, the proposed AI-VBEKF is verified. The robust and accurate estimation results proved the superior characteristics of the proposed DOA tracking method.

The remainder of this paper is organized as follows: Section 2 shows the kinematic model and the measurement model of the DOA tracking problem. In Section 3, the VB-EKF using a UCA is derived first. In addition, the attention-based DCNN is proposed to make a good estimation of the guess of the initial MSEM. Then, based on the VB-KEF and proposed DCNN, the whole frame of the AI-VBKEF is proposed. In Section 4, the simulation and experiment verification results are shown. Finally, the conclusions are drawn in Section 5.

2. Problem Formulation

2.1. Kinematic Model of the DOA Tracking Process

Since the underwater target is usually performing at a constant velocity to keep invisibility and is far from the hydrophone arrays, the motion model of the target’s bearing angle can be modeled as the constant velocity (CV) model.

Let ${\theta }_k$ denote the bearing angle of the target (${\theta }_k$ is the angle between the target and the positive x-axis with the positive counterclockwise direction) and ${\dot{\theta }}_k$ denote the change rate of ${\theta }_k$. Then, the CV can be expressed as:

$$x_k=F_{k|k-1}{x}_{k-1}+G_{k-1}{w}_{k-1}$$

(1)

where $x_k={\left({\theta }_k,{\dot{\theta }}_k\right)}^{\mathrm{T}}$ is the system state at time step k, ${\left(\cdot \right)}^{\mathrm{T}}$ denotes the matrix transposition, and $w_{k-1}$ denotes the zero mean Gaussian process noise with the covariance matrix $Q_{k-1}$ caused by the unknown underwater environment. $Q_{k-1}$ is named the process noise covariance matrix (PNCM). $F_{k|k-1}$ and $G_k$ are the state transition matrix for CV motion and the noise driving matrix at time step k, which are given by:

$$F_{k|k-1}=\left[\begin{array}{cc}1& T\\ 0& 1\end{array}\right]$$

(2)

$$G_{k-1}=\left[\begin{array}{c}{T}^2/2\\ T\end{array}\right]$$

(3)

where T is the tracking interval.

2.2. Measurement Model Based on the Received Signal of UCA

Considering that the bearing angle of the target varies in a wide range in long-time tracking, the UCA is used since it is free from the port and starboard ambiguity problem. The configuration of the UCA and the underwater target is shown in Figure 1.

Consider a narrowband acoustic signal is emitted by a target. The frequency, amplitude and the initial phase are f, a, and ${\varphi }_0$. Then, we express the target signal at time step k as:

$$s\left(k\right)=\alpha \exp \left(j\left(2\pi fkT+{\varphi }_0\right)\right)$$

(4)

The emitted signal given by Equation (4) is received by the UCA with P elements at time step k. The P elements are uniformly placed around the circle with the radius r and the center at the origin of the cartesian coordinates. It is assumed that the range of the target satisfies the far-field condition to regard the arrived signal as a plane wave, and the sound velocity is c. Then, the time delay of the signal relative to the center of the UCA at the p-th element ($p=1,2,\dots ,P$) is $\cos ({\theta }_k-2\pi p/P)r/c$. In addition, considering the signal given by Equation (4), the shift of the signal caused by the time delay at the p-th element is expressed as:

Consider that a P-element UCA with radius of r received the emitted signal at time step k. The arrived signal is assumed to be the plane wave, and the sound velocity is denoted by c. Then, the array manifold at the p-th element is expressed as:

$$A_p\left({\theta }_k\right)=\exp \left(j2\pi f\frac{r}{c}\cos \left({\theta }_k-2\pi \frac{p}{P}\right)\right)$$

(5)

Thus, the received signal at the p-th element is expressed as:

$$z_p\left(k\right)=\mathrm{real}\left(A_p\left({\theta }_k\right)s\left(k\right)\right)+v_k$$

(6)

where $\mathrm{real}(\cdot )$ denotes taking real part of a complex number, ${\theta }_k$ denotes the bearing angle of the target at time step k, and $v_k$ denotes the zero mean Gaussian measurement noise with the variance ${\sigma }_{r,k}^2$ at time step k from the measurement system and the unknown ocean environment.

We define the received signal vector of the UCA as $z_k={[z_1(k),z_2(k),\dots ,z_P(k)]}^{\mathrm{T}}$, the output of the UCA at time step k is:

$$z_k=\mathrm{real}\left(A\left({\theta }_k\right)s\left(k\right)\right)+v_k$$

(7)

where $A({\theta }_k)$ denotes the array manifold vector, which is given by

$$A({\theta }_k)={\left[A_1({\theta }_k),A_2({\theta }_k),\dots ,A_P({\theta }_k)\right]}^{\mathrm{T}}$$

(8)

M snapshots of the received signal are taken in a measurement. The snapshots are taken from a small period within every tracking interval T to guarantee that the target remains stationary. Then, we can express the measurement vector at tracking step k as:

We take M snapshots of the received signal from a small period within every tracking interval T in one measurement to guarantee that the target remains stationary. Then, the measurement vector at time step k is given as follows:

$$Z_k=\left[\begin{array}{l}z\left(k\right)\\ z\left(k+\tau \right)\\ z\left(k+2\tau \right)\\ ⋮\\ z\left(k+\left(M-1\right)\tau \right)\end{array}\right]=h_k\left(X_k\right)+V_k=\left[\begin{array}{l}\mathrm{real}\left(A\left({\theta }_k\right)s\left(k\right)\right)\\ \mathrm{real}\left(A\left({\theta }_k\right)s\left(k+\tau \right)\right)\\ \mathrm{real}\left(A\left({\theta }_k\right)s\left(k+2\tau \right)\right)\\ ⋮\\ \mathrm{real}\left(A\left({\theta }_k\right)s\left(k+\left(M-1\right)\tau \right)\right)\end{array}\right]+V_k$$

(9)

where τ denotes the interval between adjacent snapshots, $M\tau \ll T$, $h_k(x_k)$ denotes the nonlinear measurement function, and $V_k$ denotes the measurement noise.

The measurement noise of each element is assumed to be independent and identically distributed. Then, $V_k$ at time step k is defined as a zero mean Gaussian noise with covariance matrix $R_k=I_{MP}{\sigma }_{r,k}^2$, where ${\sigma }_{r,k}^2$ is the variance of the measurement noise of each element given in Equation (6), and $I_{MP}$ is the MP-order identity matrix. The measurement noise covariance matrix (MNCM) $R_k$ is usually assumed to be fixed in the underwater tracking problem. However, in the real scenario of underwater tracking, the unknown ocean environment can result in uncertain measurement noise, which makes the MNCM become fluctuant and hard to deal with by using traditional tracking methods. This problem seriously degrades the tracking performance.

3. Methods

3.1. EKF for DOA Tracking

Since the measurement model depicted by Equation (9) has high nonlinearities, the extended Kalman filter (EKF) scheme is used to derived a DOA tracking algorithm in this section.

From Equation (1), the one-step prediction of the state estimate ${\widehat{x}}_{k|k-1}$ is depicted as:

$${\widehat{x}}_{k|k-1}=F_{k|k-1}{\widehat{x}}_{k-1}$$

(10)

where $F_{k|k-1}$ is the state transition matrix given by Equation (2), and ${\widehat{x}}_{k-1}$ is the state estimate at tracking step k − 1. The one-step prediction of the MSEM $P_{k|k-1}$ is expressed as:

$$P_{k|k-1}=F_{k|k-1}{P}_{k-1}{F}_{k|k-1}{}^T+G_k{Q}_k{G}_k{}^T$$

(11)

where $P_{k-1}$ is the MSEM at step k − 1, $G_k$ is the noise driving matrix given by Equation (3), and $Q_k$ is the covariance matrix of the process noise.

Then, considering the UCA-based measurement model given in Equation (9), the Kalman filter gain $K_k$ is expressed as:

$$K_k=P_{k|k-1}{H}_k{}^T{\left(H_k{P}_{k|k-1}{H}_k^T+R_k\right)}^{-1}$$

(12)

where $H_k$ is the first order Taylor expansion of the nonlinear measurement model. $R_k$ is the MNCM. According to Equation (9), $H_k$ can be calculated as:

$$H_k={\left[\begin{array}{ccccccccccccc}{h^{\prime }}_{1,1}& {h^{\prime }}_{2,1}& \cdots & {h^{\prime }}_{P,1}& {h^{\prime }}_{1,2}& {h^{\prime }}_{2,2}& \cdots & {h^{\prime }}_{P,2}& \cdots & {h^{\prime }}_{1,M}& {h^{\prime }}_{2,M}& \cdots & {h^{\prime }}_{P,M}\\ 0& 0& \cdots & 0& 0& 0& \cdots & 0& \cdots & 0& 0& \cdots & 0\end{array}\right]}^T$$

(13)

The elements of the matrix ${h^{\prime }}_{p,m}$($p=1,2,\dots ,P$, $m=1,2,\dots ,M$) are given by [44]:

$${h^{\prime }}_{p,m}=-2\pi f\frac{r}{c}\sin \left({\widehat{\theta }}_{k|k-1}-2\pi \frac{p}{P}\right)\mathrm{real}\left(j{A}_p\left({\widehat{\theta }}_{k|k-1}\right)\widehat{s}\left(k+\left(m-1\right)\tau \right)\right)$$

(14)

where $\widehat{s}\left(k+\left(m-1\right)\tau \right)$ denotes the estimate of the target signal, which is given by [45]:

$${\widehat{s}}_{k|k}^{(i)}(k)={\left(A^{\mathrm{H}}\left({\widehat{\theta }}_{k|k}^{(i)}\right)A\left({\widehat{\theta }}_{k|k}^{(i)}\right)\right)}^{-1}{A}^{\mathrm{H}}\left({\widehat{\theta }}_{k|k}^{(i)}\right)\mathrm{Hilbert}\left(z\left(k\right)\right)$$

(15)

where ${\widehat{\theta }}_{k|k}^{(i)}$ is the predicted bearing angle, i.e., the first term of ${\widehat{x}}_{k|k-1}$, $A\left({\widehat{\theta }}_{k|k}^{(i)}\right)$ is obtained by using Equations Equations (5) and (8) ${\left(\cdot \right)}^{\mathrm{H}}$ and $\mathrm{Hilbert}\left(\cdot \right)$ denote the conjugate transposition and the Hilbert transform, respectively.

The state estimate ${\widehat{x}}_k$ at tracking step k is expressed as:

$${\widehat{x}}_k={\widehat{x}}_{k|k-1}+K_k\left(Z_k-h\left({\widehat{x}}_{k|k-1}\right)\right)$$

(16)

where $Z_k$ is the measurement of the UCA at step k. The one-step prediction of MSEM $P_{k|k-1}$ is modified by $K_k$, i.e.,

$$P_k=\left(I-K_k{H}_k\right)P_{k|k-1}$$

(17)

where $P_k$ denotes the MSEM at tracking step k. From Equations (10)–(17), the current bearing angle estimation can be calculated by the inputs of ${\widehat{x}}_0$ and $P_0$.

From the above analysis, the DOA tracking process not only depends on the measurements but also utilizes the target’s prior motion information. Therefore, the DOA tracking techniques can be robust to the target’s motion; especially, the kinematic model of the target is accurate. However, from Equation (12), it is obvious that the measurement noise $R_k$ will affect the Kalman gain, which will hugely influence the final tracking precision. In addition, from Equation (10), it is obvious that the initial values of the state and MSEM will influence the final tracking accuracy. As a result, in the scenario of accurate and robust underwater DOA tracking, both the uncertain measurement noise and proper setting of the ${\widehat{x}}_0$ and $P_0$ need to be considered.

3.2. VB-EKF for Robust DOA Tracking

Regarding the measurement model with high nonlinearity depicted by Equation (7), the extended Kalman filter (EKF) technique is used for DOA tracking in this section. Furthermore, considering the fluctuations of the MNCM caused by the unknown underwater environment, the variational Bayesian approach is utilized to improve the tracking performance by estimating the MNCM. Thus, the VB-EKF for DOA tracking using the UCA is proposed, and the details are given as follows.

3.2.1. Choice of Prior Distribution

In the framework of the standard EKF [5], the one-step predicted probability density distribution (PDF) $p(x_k|z_{1:k-1})$ and the likelihood PDF $p(z_k|x_k)$ are assumed to be subject to Gaussian distributions as follows:

$$p\left(x_k|z_{1:k-1},P_{k|k-1}\right)=\mathrm{N}\left(x_k;{\widehat{x}}_{k|k-1},P_{k|k-1}\right)$$

(18)

$$p\left(z_k|x_k,R_k\right)=\mathrm{N}\left(z_k;h\left(x_k\right),R_k\right)$$

(19)

where $\mathrm{N}(\cdot ;\mu ,\Sigma )$ denotes the PDF of the Gaussian distribution with mean μ and covariance matrix Σ, and $h(x_k)$ is the nonlinear measurement function given by Equation (7). ${\widehat{x}}_{k|k-1}$ and $P_{k|k-1}$ denote the one-step prediction of state and MSEM, respectively, which are given by Equations (10) and (11).

In order to infer $x_k$ along with $R_k$, a conjugate prior distribution needs to be selected for the fluctuant MNCM $R_k$ since a conjugate distribution can guarantee the same functional forms of the prior distribution and the posterior distribution. In the Bayesian theory, inverse Wishart distribution is usually used as the conjugate prior to the covariance matrix of a Gaussian distribution with known mean [6]. Since $R_k$ is the covariance matrix of Gaussian distribution, the prior distribution $p(R_k|z_{1:k-1})$ is selected as an inverse Wishart distribution given by:

$$p\left(R_k|z_{1:k-1}\right)=\mathrm{IW}\left(R_k;{\widehat{u}}_{k|k-1},{\widehat{U}}_{k|k-1}\right)$$

(20)

where $\mathrm{IW}(\cdot ;\lambda ,\Psi )$ denotes the PDF of the inverse Wishart distribution with degree of freedom (dof) λ and inverse scale matrix Ψ [3], ${\widehat{u}}_{k|k-1}$ and ${\widehat{U}}_{k|k-1}$ are the dof and the inverse scale matrix of $p(R_k|z_{1:k-1})$, respectively.

The posterior distribution $p(R_{k-1}|z_{1:k-1})$ is also subject to an inverse Wishart distribution as follows:

$$p\left(R_{k-1}|z_{1:k-1}\right)=\mathrm{IW}\left(R_{k-1};{\widehat{u}}_{k-1|k-1},{\widehat{U}}_{k-1|k-1}\right)$$

(21)

To guarantee that $p(R_k|z_{1:k-1})$ is the inverse Wishart distribution given by Equation (20), the previous approximate posteriors is spread through a forgetting factor $\rho \in \left(0,1\right]$, which indicates the extent of the time-fluctuations of the MNCM. Then, the prior dof ${\widehat{u}}_{k|k-1}$ and the prior inverse scale matrix ${\widehat{U}}_{k|k-1}$ are given as follows:

$${\widehat{u}}_{k\mid k-1}=\rho \left({\widehat{u}}_{k-1\mid k-1}-n-1\right)+n+1$$

(22)

$${\widehat{U}}_{k\mid k-1}=\rho {\widehat{U}}_{k-1\mid k-1}$$

(23)

where n denotes the order of the MNCM $R_k$.

3.2.2. Variational Approximations of Posterior PDFs

According to the variational Bayesian approximation, the joint posterior PDF of the state $x_k$ and the MNCM $R_k$ is approximated to

$$p\left(x_k,R_k\mid z_{1:k}\right)\approx q\left(x_k\right)q\left(R_k\right)$$

(24)

where $q(x_k)$ and $q(R_k)$ are the approximate posterior PDF of $x_k$ and $R_k$, respectively. The variational Bayesian approximation is formed by minimizing the Kullback–Leibler divergence (KLD) between the true joint distribution $p(x_k,R_k\mid z_{1:k})$ and the approximate distribution $q(x_k)q(R_k)$, i.e.,

$$\left\{q\left(x_k\right),q\left(R_k\right)\right\}=\arg \min \mathrm{KLD}\left(q\left(x_k\right)q\left(R_k\right)\parallel p\left(x_k,R_k\mid z_{1:k}\right)\right)$$

(25)

where $\mathrm{KLD}(q(x)\parallel p(x))$ denotes the KLD between $q(x)$ and $p(x)$, and

$$\mathrm{KLD}\left(q\left(x\right)\parallel p\left(x\right)\right)={\displaystyle \int q(x)\log \frac{q(x)}{p(x)}dx}$$

(26)

The optimal solution of Equation (25) satisfies the following Equations [4]:

$$\log q(x_k)={\mathrm{E}}_{R_k}\left[\log p\left(x_k,R_k,z_{1:k}\right)\right]+c_x$$

(27)

$$\log q(R_k)={\mathrm{E}}_{x_k}\left[\log p\left(x_k,R_k,z_{1:k}\right)\right]+c_R$$

(28)

where ${\mathrm{E}}_{x_k}[\cdot ]$ and ${\mathrm{E}}_{R_k}[\cdot ]$ denote the expectation with regard to $x_k$ and $R_k$, respectively, and $c_x$ and $c_R$ denote the constants with respect to $x_k$ and $R_k$, respectively. Since the variational parameters of $q(x_k)$ and $q(R_k)$ are coupled, a fix-point iteration process is applied to solve Equations (27) and (28), i.e., the approximate posterior PFD $q(x_k)$ is updated to $q^{\left(i+1\right)}(x_k)$ at the $i+1$-th iteration using the posterior PDF $q^{\left(i\right)}(R_k)$, and $q(R_k)$ is updated to $q^{\left(i+1\right)}(R_k)$ using the posterior $q^{\left(i\right)}(x_k)$.

According to Equations (18)–(20), the joint PDF is expressed as

$$\begin{array}{l}p\left(x_k,R_k,z_{1:k}\right)\\ =p\left(z_k|x_k,R_k\right)p\left(x_k|z_{1:k-1}\right)p\left(R_k|z_{1:k-1}\right)p\left(z_{1:k-1}\right)\\ =\mathrm{N}\left(z_k;h\left(x_k\right),R_k\right)\mathrm{N}\left(x_k;{\widehat{x}}_{k|k-1},P_{k|k-1}\right)\times \mathrm{IW}\left(R_k;{\widehat{u}}_{k|k-1},{\widehat{U}}_{k|k-1}\right)p\left(z_{1:k-1}\right)\end{array}$$

(29)

(1)

Update of $x_k$

The posterior $q^{(i+1)}(x_k|z_{1:k-1})$ is updated according to Equations (12), (16) and (17) as:

$$q^{(i+1)}\left(x_k|z_{1:k-1}\right)=\mathrm{N}\left(x_k;{\widehat{x}}_{k|k}^{(i+1)},{\widehat{P}}_{k|k}^{(i+1)}\right)$$

(30)

where the mean vector ${\widehat{x}}_{k|k}^{(i+1)}$ and the covariance matrix ${\widehat{P}}_{k|k}^{(i+1)}$ are given as follows:

$${\mathit{K}}_k^{(i+1)}={\mathit{P}}_{k|k-1}{\left({\mathit{H}}_k^{(i)}\right)}^{\mathrm{T}}{\left({\mathit{H}}_k^{(i)}{\mathit{P}}_{k|k-1}{\left({\mathit{H}}_k^{(i)}\right)}^{\mathrm{T}}+{\widehat{\mathit{R}}}_k^{(i)}\right)}^{-1}$$

(31)

$${\widehat{\mathit{x}}}_{k|k}^{(i+1)}={\widehat{\mathit{x}}}_{k|k-1}+{\mathit{K}}_k^{(i+1)}\left(z_k-h\left({\widehat{\mathit{x}}}_{k|k-1}\right)\right)$$

(32)

$${\mathit{P}}_{k|k}^{(i+1)}={\mathit{P}}_{k|k-1}-{\mathit{K}}_k^{(i+1)}{\mathit{H}}_k^{(i)}{\mathit{P}}_{k|k-1}$$

(33)

where $H_k^{(i)}$ denotes the Jacobian matrix of the measurement function. $H_k^{(i)}$ is obtained by substituting ${\widehat{\mathit{x}}}_{k|k}^{(i)}$ into $H_k$ given by Equation (13).

Based on the new estimated state ${\widehat{x}}_{k|k}^{(i)}$ and MSEM $P_{k|k}^{(i)}$, a more accurate approximation of $h(x_k)$ can be obtain by performing linearization with ${\widehat{x}}_{k|k}^{(i)}$ [6], i.e.,

$$h\left(x_k\right)=h\left({\widehat{x}}_{k|k}^{(i)}\right)+H_k^{(i)}\left(x_k-{\widehat{x}}_{k|k}^{(i)}\right)$$

(34)

The estimate of the signal ${\widehat{s}}_{k|k}^{(i)}(k)$ is also used instead of the real one to calculate the measurement function $h({\widehat{x}}_{k|k}^{(i)})$ by using Equation (7). Thus, Equation (32) is upgraded to

$${\widehat{\mathit{x}}}_{k|k}^{(i+1)}={\widehat{\mathit{x}}}_{k|k-1}+{\mathit{K}}_k^{(i+1)}\left(z_k-h\left({\widehat{\mathit{x}}}_{k|k}^{(i)}\right)-{\mathit{H}}_k^{(i)}\left({\widehat{\mathit{x}}}_{k|k-1}-{\widehat{\mathit{x}}}_{k|k}^{(i)}\right)\right)$$

(35)

(2)

Update of $R_k$

According to Equations (18)–(20), $\log q^{\left(i\right)}(R_k)$ is given by

$$\begin{array}{l}\log q^{(i+1)}\left(R_k\right)\\ =-0.5\left(n+{\widehat{u}}_{k|k-1}+2\right)\log \left|R_k\right|-0.5\mathrm{tr}\left({\widehat{U}}_{k|k-1}{R}_k^{-1}\right)-0.5{\left(z_k-h\left(x_k\right)\right)}^{\mathrm{T}}{R}_k^{-1}\left(z_k-h\left(x_k\right)\right)+c_R\\ =-0.5\left(n+{\widehat{u}}_{k|k-1}+2\right)\log \left|R_k\right|-0.5\mathrm{tr}\left(\left(B_k^{(i)}+{\widehat{U}}_{k|k-1}\right)R_k^{-1}\right)+c_R\end{array}$$

(36)

where

$$B_k^{(i)}={\mathrm{E}}^{\left(i\right)}\left[\left(z_k-h\left(x_k\right)\right){\left(z_k-h\left(x_k\right)\right)}^{\mathrm{T}}\right]$$

(37)

Similar to Equation (34), $h\left(x_k\right)$ is linearized as

$$h\left(x_k\right)=h\left({\widehat{x}}_{k|k}^{(i+1)}\right)+H_k^{(i+1)}\left(x_k-{\widehat{x}}_{k|k}^{(i+1)}\right)$$

(38)

where $H_k^{(i+1)}$ denotes the Jacobian matrix of the measurement function $h(x_k)$ at ${\widehat{x}}_{k|k}^{(i+1)}$, and $H_k^{(i+1)}$ is obtained by substituting ${\widehat{x}}_{k|k}^{(i+1)}$ into Equation (13). Substituting Equation (38) into Equation (37), we obtain:

$$\begin{array}{l}{B}_k^{(i)}\\ ={\mathrm{E}}^{\left(i\right)}\left[\left(z_k-h\left(x_k\right)\right){\left(z_k-h\left(x_k\right)\right)}^{\mathrm{T}}\right]\\ ={\mathrm{E}}^{\left(i\right)}\left[\left(z_k-h\left({\widehat{x}}_{k|k}^{(i)}\right)-H_k^{(i)}\left(x_k-{\widehat{x}}_{k|k}^{(i)}\right)\right)\right.\left.{\left(z_k-h\left({\widehat{x}}_{k|k}^{(i)}\right)-H_k^{(i)}\left(x_k-{\widehat{x}}_{k|k}^{(i)}\right)\right)}^{\mathrm{T}}\right]\\ =\left(z_k-h\left({\widehat{x}}_{k|k}^{(i)}\right)\right){\left(z_k-h\left({\widehat{x}}_{k|k}^{(i)}\right)\right)}^{\mathrm{T}}+H_k^{(i)}{\mathrm{E}}^{\left(i\right)}\left[\left(x_k-{\widehat{x}}_{k|k}^{(i)}\right){\left(x_k-{\widehat{x}}_{k|k}^{(i)}\right)}^{\mathrm{T}}\right]H_k^{\mathrm{T}}\\ =\left(z_k-h\left({\widehat{x}}_{k|k}^{(i)}\right)\right){\left(z_k-h\left({\widehat{x}}_{k|k}^{(i)}\right)\right)}^{\mathrm{T}}+H_k^{(i)}{P}_{k|k}^{(i)}{H}_k^{\mathrm{T}}\end{array}$$

(39)

From Equation (36), $q^{(i+1)}(R_k)$ is updated as

$$q^{(i+1)}\left(R_k\right)=\mathrm{IW}\left(R_k;{\widehat{u}}_k^{(i+1)},{\widehat{U}}_k^{(i+1)}\right)$$

(40)

where the dof ${\widehat{u}}_k^{(i+1)}$ and the inverse scale matrix ${\widehat{U}}_k^{(i+1)}$ are given as follows:

$${\widehat{u}}_k^{(i+1)}={\widehat{u}}_{k\mid k-1}+1$$

(41)

$${\widehat{U}}_k^{(i+1)}=B_k^{(i)}+{\widehat{U}}_{k\mid k-1}$$

(42)

Then, according to Equation (29), $\log q^{\left(i\right)}(x_k)$ is given by:

$$\begin{array}{c}\log q^{(i+1)}\left(x_k\right)=-0.5{\left(z_k-h\left(x_k\right)\right)}^{\mathrm{T}}{\mathrm{E}}^{(i+1)}\left[R_k^{-1}\right]\left(z_k-h\left(x_k\right)\right)\\ -0.5{\left(x_k-{\widehat{x}}_{k|k-1}\right)}^{\mathrm{T}}{P}_{k|k-1}^{-1}\left(x_k-{\widehat{x}}_{k|k-1}\right)+c_x\end{array}$$

(43)

where ${\mathrm{E}}^{(i+1)}[R_k^{-1}]$ is given by:

$${\mathrm{E}}^{(i+1)}\left[R_k^{-1}\right]=\left({\widehat{u}}_k^{(i+1)}-m-1\right){\left({\widehat{U}}_k^{(i+1)}\right)}^{-1}$$

(44)

The modified one step predicted PDF $p^{(i+1)}\left(z_k|x_k\right)$ at the I + 1th iteration is defined as:

$$p^{(i+1)}\left(z_k|x_k\right)=\mathrm{N}\left(z_k;H_k{x}_k,{\widehat{R}}_k^{(i+1)}\right)$$

(45)

where the modified MNCM ${\widehat{R}}_k^{(i+1)}$ are formulated as:

$${\widehat{R}}_k^{(i+1)}={\left\{{\mathrm{E}}^{(i+1)}\left[R_k^{-1}\right]\right\}}^{-1}={\widehat{U}}_k^{(i+1)}/\left({\widehat{u}}_k^{(i+1)}-m-1\right)$$

(46)

Finally, after N fixed-point iterations, the variational approximations of the posterior PDFs are given as follows:

$$q\left(x_k\right)\approx q^{(N)}\left(x_k\right)=\mathrm{N}\left(x_k;{\widehat{x}}_{k|k}^{(N)},P_{k|k}^{(N)}\right)=\mathrm{N}\left(x_k;{\widehat{x}}_{k|k},P_{k|k}\right)$$

(47)

$$q\left(R_k\right)\approx q^{(N)}\left(R_k\right)=\mathrm{IW}\left(R_k;{\widehat{u}}_k^{(N)},{\widehat{U}}_k^{(N)}\right)=\mathrm{IW}\left(R_k;{\widehat{u}}_{k\mid k},{\widehat{U}}_{k\mid k}\right)$$

(48)

Combining the EKF depicted by Equations (10)–(17) and the VB estimating technique, the state can be robustly estimated under uncertain measurement noise.

3.3. AI-VBEKF for Robust DOA Tracking

Considering the high nonlinearities of the measurement model by Equation (7), the initial values will affect the tracking results. Since the EKF is based on the Taylor series expansion to linearize the nonlinear system model, the truncation error will become unacceptable if the error between the initial guess of the state and its real value is quite large. In addition, if the initial MSEM cannot be set properly, the DOA tracking results will be in low accuracy or even diverged. In the real underwater DOA tracking scenario, the initial DOA of the target can be obtained by traditional DOA estimation methods (i.e., CBF, MVDR, MUSIC, et.al.), but the initial MSEM cannot be determined since the initial error between the real state and the estimated one is unknown. If the initial MSEM is set too small, the DOA tracking algorithm will converge very slowly or cannot converge at all. On the contrary, if the initial MSEM is set too large, the DOA tracking process will be soon diverged since the Kalman gain will become uncontrolled during the measurement update steps for nonlinear tracking systems. As a result, a proper preset initial MSEM not only affects the final DOA tracking accuracy but also determines the convergence of the total tracking procedure. However, in the real DOA tracking scenario, the initial MSEM can only be set by engineering experience. To solve this problem and deal with the uncertain measurement noise, an AI-aided robust DOA tracking algorithm is proposed in this subsection.

In order to minimalize the final DOA tracking error caused by the inaccurate initial guess of the MSEM, an AI-aided robust DOA tracking algorithm is proposed by this subsection. Firstly, considering that the different initial state guesses will lead to different theoretical measurement calculations by the measurement model, the error between the real measurements and the theoretical ones contains the information of the error of the initial guess of the MSEM between its true value. Thus, based on the difference between the theoretical measurement calculations and the real measurements, the initial MSEM can be estimated via deep learning techniques. By utilizing the covariance matrix of the difference between the theoretical measurement calculations and the real measurements from all bearing angles with different initial errors as inputs, an attention-based deep convolutional neural network is proposed by this subsection to output reliable initial MSEM. Then, utilizing the estimated MSEM as the input of the VB-EKF, the AI-VBEKF for a robust DOA tracking algorithm is finally carried out in this subsection.

3.3.1. Input Data Processing

Supposing the initial state is selected as:

$$x_0={\left({\theta }_0,{\dot{\theta }}_0\right)}^{\mathrm{T}}$$

(49)

Using the measurement model depicted by Equation (9), the error between the theoretical measurement and the real measurement can be presented as:

$${\tilde{Z}}_k=Z_k-h_0\left(x_0\right)=\left[\begin{array}{l}z\left(0\right)\\ z\left(\tau \right)\\ z\left(2\tau \right)\\ ⋮\\ z\left(\left(M-1\right)\tau \right)\end{array}\right]-\left[\begin{array}{l}\mathrm{real}\left(A\left({\theta }_0\right)s\left(0\right)\right)\\ \mathrm{real}\left(A\left({\theta }_0\right)s\left(\tau \right)\right)\\ \mathrm{real}\left(A\left({\theta }_0\right)s\left(2\tau \right)\right)\\ ⋮\\ \mathrm{real}\left(A\left({\theta }_0\right)s\left(\left(M-1\right)\tau \right)\right)\end{array}\right]$$

(50)

Using Equation (50), the quadratic form of the error between the theoretical measurement and the real measurement can be presented as:

$$e_V={\tilde{Z}}_k{\tilde{Z}}_{{}_k}^T$$

(51)

where $e_V$ is a matrix with the dimension of $MP\times MP$.

From Equation (7), it can be found that the real measurement has a strong relationship with the real DOA angle of the target and its angular velocity. When the value of the real state changes, the measurement will also change. In addition, the estimated state influences the value of the theoretical measurement computed by the measurement model. Thus, the computed $e_V$ has the information of the difference of the real values of the state and its estimated values. From the view of data-driven techniques, $e_V$ can be utilized to output this error and can make good estimation of the guess of the initial MSEM.

In addition, from Equation (51), it is obvious that the size of $e_V$ depends on the number of the snapshots in one measurement time interval. In order to make reliable DOA estimations and DOA tracking, the value of M is usually set high during one measurement updating step (i.e., the same as the sampling rate of the sonar system). As a result, the size of $e_V$ is usually large, and the nonlinear relationship between the real guess of the initial MSEM and its real values is hard to be dug out by shallow BP neural networks. In order to make sufficient utilization of the data resources obtained by the sonar system and make reliable estimations of the guess of the initial MSEM, an attention-based DCNN is proposed in the following subsection.

3.3.2. Attention-Based Deep Convolutional Neural Network

As mentioned in Section 3.3.1, since the traditional BP neural networks have limited ability in dealing with a large matrix and the computed $e_V$ has the form of a symmetric matrix such as an image, an attention-based DCNN is proposed by this subsection. Among the latest DCNNs, the attention-based DCNN is the most popular one and has been utilized in a number of target classification, localization, and tracking missions. Unlike the basic DCNN that cannot be developed very deep and the Resnet that cannot give different features different weights, which make the overall performance degrade, the attention-based DCNN not only utilizes the advantages of the Resnet to make the DCNN much deeper than the traditional ones, but also develops the attention scheme that values different features with different weights to make the DCNN have better tracking performance. The basic structure of the attention-based DCNN is shown in Figure 2 [35].

From Figure 2 [35], the main difference between the attention-based DCNN and the other DCNNs is the “squeeze and excitation” block. For any given data with the dimension $C\prime \times W\prime \times H\prime$, firstly a transformation mapping $F_{tr}$ is operated to generate the feature maps U, which can be depicted as:

$$F_{tr}:X\to U,X\in {ℝ}^{H^{\prime }\times W^{\prime }\times C^{\prime }},U\in {ℝ}^{H\times W\times C}$$

(52)

From Equation (52), it is obvious that the $F_{tr}$ operation has no difference with other DCNNS (i.e., $F_{tr}$ is the traditional DCNN or Resnet). As a result, the attention-based DCNNs can be fitted into any existing DCNNs to enhance their performance.

After the $F_{tr}$ mapping operation, the “squeeze and excitation” block can operate. Firstly, the feature maps U is operated via a squeeze process that aims to compress the information of the feature maps into a one-dimensional vector. By operating in this way, the compressed feature of every channel of the feature maps U is generated. The squeeze operation can be depicted as the following equation:

$$z_c=F_{sq}\left(u_c\right)=\frac{1}{H\times W}\sum _H^{i=1}\sum _W^{j=1}{u}_c(i,j)$$

(53)

where $u_c$ is the $c$-th channel of the feature maps U, and $F_{sq}$ is the squeeze operation.

From Equation (53), it can be found that the squeeze operation operates on every feature map with the dimension of $H\times W$. After the squeeze operation, every feature map of U turns to a real number so that the information is maximally compressed. Then, the excitation operation is processed to calculate the weights of every feature map. The excitation operation aims to fully capture the channel-wise dependencies of the input feature maps U. By utilizing full connection layers and the ReLU activation function, the nonlinear mapping process made by the excitation operation can be depicted as:

$$s=F_{ex}(z_c,W)=\sigma (g(z_c,W))=\sigma \left(W_2\delta \left(W_1{z}_c\right)\right)$$

(54)

where $s$ is the output of the excitation operation, $F_{ex}$ is the excitation operation, $W_1$ is the first fully connected (FC) layer that can be regarded as a nonlinear operation to extract the inner information of the input $z_c$, $\delta \left(\mathbf{·}\right)$ is the ReLU activating function, $W_2$ is the second FC layer, and $\sigma \left(\mathbf{·}\right)$ is a typical activating function that can be selected as a sigmoid function in general.

From Equation (54), after the excitation operation, the inner relationships between every feature map in U can be represented by $s$. Thus, $s$ is the core of the attention-based DCNN, for it depicts the weights of every feature inside the whole feature maps U shown by Equation (52). As long as $s$ is obtained, the weights vector of every feature map is known. As a result, the final output of the attention-based DCNN can be represented as:

$${\tilde{x}}_c=F_{\mathrm{scale}}\left(u_c,s_c\right)=s_c\cdot u_c$$

(55)

where $F_{\mathrm{scale}}$ is the final scaling operation that refers to the channel-wise multiplication between the $c$-th element of the output of the excitation operation and the $c$-th feature map of U. It is obvious from Equation (52) and Equation (55) that the original feature maps U is weighted by the “squeeze and excitation” block. This weighting mechanism is similar to paying attention to different features such that the whole deep learning technique is called the attention-based DCNN.

From Equation (52) and the above statements, the first mapping operation $F_{tr}$ can be any type of DCNN technique. As the Resnet has proven its superior performance in underwater target localization missions [38], the $F_{tr}$ is selected as the residual model [38] and the attention-based DCNN blocks combined with the residual model can be shown as the following block as shown in Figure 3:

Equations (52)–(55) represent the whole theoretical derivations of the attention-based DCNN. In addition, since this DCNN is mainly based on the “squeeze and excitation” block, the mentioned attention-based DCNN is also called the Senet. Here, substitute Equation (51) into Equation (52), and change the initial input $X$ to $e_V$. The Senet can then be directly utilized to analyze the inner nonlinear relationship between the real values of the state and its estimated values.

3.3.3. Design of the AI-VBEKF

From Equation (9), if the measurement noise varies for the same initial error between the initial guess of the state and its real value, $e_V$ computed by Equation (51) will change by different $R_k$. Supposing the number of varying covariance matrix of $R_k$ is $H\prime$ and the measurement $Z_k$ has snapshots with the number of M, $e_V$ will have the dimensions of $H\prime \times M\times M$. Then, substituting Equation (51) into Equation (52) and changing the initial input $X$ to $e_V$ with $H\prime$ covariance matrix of the measurement noise and M snapshots, Equation (52) can be represented as:

$$F_{tr}:e_V\to U,e_V\in {ℝ}^{H^{\prime }\times M\times M},U\in {ℝ}^{H\times W\times C}$$

(56)

In addition, if the transformation mapping $F_{tr}$ is chosen as the residual model, and the following steps are the same as the aforementioned ones (Equations (53)–(55)), the attention-based DCNN combined with the residual model can be proposed. By using $e_V$ with different $R_k$ as inputs and the differences between the initial guess of the state and its real value as outputs, the attention-based DCNN for the initial guess of the error of the covariance matrix of the states can be proposed as in

Table 1. The attention-based DCNN for estimation of the guess of the initial MSEM.

Layer Name	[Kernel Size, Filters] $\times$ No. of Blocks
SE-Comv1	$[5\times 5,64]\times 1$
SE-Conv2	$\left[\begin{array}{cc}1\times 1& 64\\ 3\times 3& 64\\ 1\times 1& 256\\ fullconnection& \left[16,256\right]\end{array}\right]\times 3$
SE-Conv3	$\left[\begin{array}{cc}1\times 1& 128\\ 3\times 3& 128\\ 1\times 1& 512\\ fullconnection& \left[32,512\right]\end{array}\right]\times 4$
SE-Conv4	$\left[\begin{array}{cc}1\times 1& 256\\ 3\times 3& 256\\ 1\times 1& 1024\\ fullconnection& \left[64,1024\right]\end{array}\right]\times 6$
SE-Conv5	$\left[\begin{array}{cc}1\times 1& 512\\ 3\times 3& 512\\ 1\times 1& 2048\\ fullconnection& \left[128,2048\right]\end{array}\right]\times 3$
SE-Global	Global average pool, fully-connected layer with ReLU

.

Then, by combining the proposed attention-based DCNN, the EKF DOA tracking scheme and the VB robust estimating technique, the AI-VBEKF, which has the abilities of estimating the guess of the initial MSEM and performing robust DOA tracking, can be depicted as Algorithm 1.

Algorithm 1: AI-VBEKF

Input$x_0$. Calculate${\tilde{Z}}_k$ as Equation (50). Calculate $e_V$ as Equation (51). Initial MSEM ${\mathit{P}}_0$ estimation by attention-based DCNN. Inputs: ${\widehat{\mathit{x}}}_{k-1|k-1}$, ${\mathit{P}}_{k-1|k-1}$, ${\mathit{z}}_k$, ${\widehat{u}}_{k-1|k-1}$, ${\widehat{\mathit{U}}}_{k-1|k-1}$, $Q_{k-1}$, ρ, N. Time update ${\widehat{\mathit{x}}}_{k|k-1}={\mathit{F}}_{k|k-1}{\widehat{\mathit{x}}}_{k-1|k-1}$. ${\mathit{P}}_{k|k-1}={\mathit{F}}_{k|k-1}{\mathit{P}}_{k-1|k-1}{\mathit{F}}_{k|k-1}^{\mathrm{T}}+{\mathit{G}}_{k-1}{\mathit{Q}}_{k-1}{\mathit{G}}_{k-1}^{\mathrm{T}}$. ${\widehat{u}}_{k\mid k-1}=\rho ({\widehat{u}}_{k-1|k-1}-n-1)+n+1$,${\widehat{\mathit{U}}}_{k|k-1}=\rho {\widehat{\mathit{U}}}_{k-1|k-1}$. Iteration measurement update Initialization: ${\widehat{\mathit{x}}}_{k|k}^{(0)}={\widehat{\mathit{x}}}_{k|k-1}$, ${\widehat{\mathit{R}}}_k^{(0)}={\widehat{\mathit{U}}}_{k|k-1}/{\widehat{u}}_{k|k-1}$. For $i=0:N-1$ Update $q^{(i+1)}({\mathit{x}}_k)=\mathrm{N}({\mathit{x}}_k;{\mathit{x}}_{k|k}^{(i+1)},{\mathit{P}}_{k|k}^{(i+1)})$ $h({\widehat{\mathit{x}}}_{k|k}^{(i)})$ is calculated by Equation (7), ${\mathit{H}}_k^{(i)}$ is calculated by Equation (13). ${\mathit{K}}_k^{(i+1)}={\mathit{P}}_{k|k-1}{({\mathit{H}}_k^{(i)})}^{\mathrm{T}}{({\mathit{H}}_k^{(i)}{\mathit{P}}_{k|k-1}{({\mathit{H}}_k^{(i)})}^{\mathrm{T}}+{\widehat{\mathit{R}}}_k^{(i)})}^{-1}$, ${\widehat{\mathit{x}}}_{k|k}^{(i+1)}={\widehat{\mathit{x}}}_{k|k-1}+{\mathit{K}}_k^{(i+1)}(z_k-h({\widehat{\mathit{x}}}_{k|k}^{(i)})-{\mathit{H}}_k^{(i)}({\widehat{\mathit{x}}}_{k|k-1}-{\widehat{\mathit{x}}}_{k|k}^{(i)}))$, ${\mathit{P}}_{k|k}^{(i+1)}={\mathit{P}}_{k|k-1}-{\mathit{K}}_k^{(i+1)}{\mathit{H}}_k^{(i)}{\mathit{P}}_{k|k-1}$. Update $q^{(i+1)}({\mathit{R}}_k)=\mathrm{IW}({\mathit{R}}_k;{\widehat{u}}_{k|k}^{(i+1)},{\widehat{\mathit{U}}}_{k|k}^{(i+1)})$ $h({\widehat{x}}_{k|k}^{(i+1)})$ is calculated by Equation (7), ${\mathit{H}}_k^{(i+1)}$ is calculated by Equation (13). ${\mathit{B}}_k^{(i+1)}=(z_k-h({\widehat{\mathit{x}}}_{k|k}^{(i+1)})){(z_k-h({\widehat{\mathit{x}}}_{k|k}^{(i+1)}))}^{\mathrm{T}}+{\mathit{H}}_k^{(i+1)}{\mathit{P}}_{k|k}^{(i+1)}{({\mathit{H}}_k^{(i+1)})}^{\mathrm{T}}$, ${\widehat{\mathit{U}}}_{k|k}^{(i+1)}={\widehat{\mathit{U}}}_{k|k-1}+{\mathit{B}}_k^{(i+1)}$, ${\widehat{u}}_{k|k}^{(i+1)}={\widehat{u}}_{k|k-1}+1$, ${\widehat{\mathit{R}}}_k^{(i+1)}={\widehat{\mathit{U}}}_{k|k}^{(i+1)}/{\widehat{u}}_{k|k}^{(i+1)}$. End for ${\widehat{\mathit{x}}}_{k|k}={\widehat{\mathit{x}}}_{k|k}^{(N)}$, ${\mathit{P}}_{k|k}={\mathit{P}}_{k|k}^{(N)}$, ${\widehat{u}}_{k|k}={\widehat{u}}_{k|k}^{(N)}$, ${\widehat{\mathit{U}}}_{k|k}={\widehat{\mathit{U}}}_{k|k}^{(N)}$. Outputs: ${\widehat{\mathit{x}}}_{k|k}$, ${\mathit{P}}_{k|k}$, ${\widehat{u}}_{k|k}$, ${\widehat{\mathit{U}}}_{k|k}$.

3.3.4. Performance Metrics

In the underwater target DOA tracking scenario, the angular velocity of the target is quite small, for the target usually moves at a relatively far distance from the observer. To obtain the accurate initial MNCM of the target, only the difference between the real initial angle and its estimation will be determined by the attention-based DCNN. According to the goal of the estimation of the initial guess of the error of the covariance matrix of the states, the metrics for quantifying the performance of the proposed attention-based DCNN can be described as:

$$E_{acc,\theta }=\frac{1}{N}{\displaystyle \sum _{i=1}^N\frac{\left|{\overline{\theta }}_0-{\theta }_0\right|}{{\theta }_0}}$$

(57)

where ${\overline{\theta }}_0$ is the estimated initial guess of the DOA angle, and ${\theta }_0$ is the real initial values of the DOA angle and its angular velocity. $N$ is the total number of the training data.

4. Results and Discussions

4.1. Simulations

4.1.1. Simulation Scenario and Data Set Generation of the Attention-Based DCNN

An underwater target is assumed to perform in CV motion mode. The initial state $x_0$ is set to [270, −0.03], and process noise covariance matrix $Q_{k-1}$ is set to $2\times {10}^{-4}°/{\mathrm{s}}^2$. The total simulation time is 5000 s with a 1 s interval. The frequency f, amplitude a and the initial phase ${\varphi }_0$ of the target signal are set to 170 Hz, 1 and 0 rad, respectively. The underwater sound velocity c is assumed to be 1500 m/s. A UCA with 12 elements is used to obtain the measurements. The sample frequency of the array signal is 1000 Hz.

The MNCM $R_k$ is set to $I_{MP}{\sigma }_r^2$, where ${\sigma }_r^2$ is set as 20 to simulate a low SNR tracking scenario. The uncertainties of the MNCM appear from 3600 to 4000 s to test the robustness of the VB-EKF to the uncertain measurement noise caused by the unknown underwater environment. Then, $R_k$ is expressed as:

$$R_k=\left\{\begin{array}{ll}{I}_{MP}{\sigma }_r^2& 0\le k<3600,4000<k\le 5000\\ 10I_{MP}{\sigma }_r^2& 3600\le k\le 4000\end{array}\right.$$

(58)

The proposed AI-VBEKF for DOA tracking given by Algorithm 1 is tested by the simulated measurement data. The EKF, the SH-EKF for DOA tracking [44] are also tested for comparison. In addition, to illustrate the superior performance of the proposed AI-VBEKF, $P_k$ is initialized to different initial guesses as $P_0=\mathrm{diag}\left([25,{10}^{-5}]\right)$, $P_0=\mathrm{diag}\left([50,{10}^{-5}]\right)$, and $P_0=\mathrm{diag}\left([100,{10}^{-5}]\right)$. The process noise covariance matrix $Q_k$ for EKF, SH-EKF and VB-EKF is set to the true value. The MNCM of EKF is set to a constant $I_{MP}$. The initial estimate of the MNCM and the forgetting factor b are set to $I_{MP}$ and 0.9 for the SH-EKF. The dof ${\widehat{u}}_0$, the inverse scale matrix ${\widehat{\mathit{U}}}_0$ of the inverse Wishart distribution and the forgetting factor ρ for VB-EKF are set to 16, $16I_{MP}$, and 0.9995, respectively.

To numerically compare the tracking precision of the above methods, the bearing angle estimation error (BEE) is defined as a performance metrics, which is given as follows:

$$BEE\left(k\right)=\left|{\widehat{\theta }}_k-{\theta }_k^{}\right|$$

(59)

where ${\widehat{\theta }}_k$ is the bearing angle estimate at step k, and ${\theta }_k^{}$ is the real bearing angle.

Furthermore, 500 Monte Carlo simulations are performed to statistically evaluate the performance of the VB-EKF. The root mean square error (RMSE) for DOA tracking results at time step k is defined as

$${\mathrm{RMSE}}_{\theta }=\sqrt{\frac{1}{N}{{\displaystyle \sum _{n=1}^N\left({\widehat{\theta }}_{k,n}-{\theta }_k\right)}}^2}$$

(60)

where N is the number of the Monte Carlo trials, ${\widehat{\theta }}_{k,n}$ is the estimate of bearing angle in the n-th trail. The average RMSE of the bearing angle estimate in total time is defined as:

$${\stackrel{\_}{RMSE}}_{\theta }=\sqrt{\frac{1}{KN}{{\displaystyle \sum _{k=1}^K{\displaystyle \sum _{n=1}^N\left({\widehat{\theta }}_{k,n}-{\theta }_k\right)}}}^2}$$

(61)

where K is the total number of tracking steps. A smaller RMSE represents a higher precision of the DOA tracking.

4.1.2. Data Set Generation and Training of the Attention-Based DCNN

The training data for the attention-based DCNN for initial MSEM estimations are generated as the following steps:

Firstly, ${\theta }_0$ of the initial state ${\mathbf{x}}_0$ is set as ${\theta }_0\in \left[0,360\right]$ deg with intervals of 1 deg.

In addition, the offset of the initial guess of ${\theta }_0$ is set as: $\Delta {\theta }_0\in \left[0,20\right]$ deg with the interval of 0.1 deg.

Then, the MNCM $R_k$ is set as:$R_k=\alpha I_{MP}{\sigma }_{r,k}^2$, where the parameter $\alpha$ is set as $\alpha \in \left[1,2\right]$ with the interval of 0.1 to generate different measurement noise.

Finally, assuming that a UCA has 12 elements and the number of snapshots in one measurement is set as 12, by utilizing the different generated ${\theta }_0$ along with $R_k$ and using Equations (4)–(9), (50) and (51), the data set for training and testing can be generated.

From the above data-generating methods, the dimension of every input data of the attention-based DCNN is $144\times 144$, and the total number of data are 1,440,000.

By utilizing the generated training data and the proposed attention-based DCNN depicted by Table 1, the proposed attention-based DCNN can be properly trained, and the averaged cross entropy loss for the training process is shown in Figure 4.

4.1.3. Simulation Results

The DOA tracking results of the EKF, SH-EKF and AI-VBEKF are shown in Figure 5, and Table 1 shows the average RMSEs. Since the CBF is only a DOA estimation technique, we merely compared different DOA tracking algorithms in the simulation section to show the robustness and accuracy of the proposed AI-VBEKF. When the initial MSEM is set properly ($P_0=\mathrm{diag}\left([25,{10}^{-5}]\right)$), all of the DOA tracking techniques can output acceptable DOA tracking results. The trajectory obtained by the EKF, SH-EKF, and AI-VBEKF fluctuated in the beginning. After about 1000 s, the algorithms converged, and then, the tracking stabilized. The reason is that the error in the initial state estimate and MSEM affect the performance at the beginning of the tracking. The error was gradually eliminated by the correction of the new measurements after a certain period. When the uncertainties of the measurement noise appear, the DOA estimation accuracy of the EKF degrades. However, by utilizing the prior information of the kinematics model of the underwater target, the DOA tracking results by the EKF are still within a certain threshold (see in Figure 5b).

Since the SH-EKF estimates the MNCM in real time by using the Sage–Husa algorithm to provide accurate Kalman filter gain, the precision of SH-EKF is higher than EKF. By giving the prior distribution of MNCM, the AI-VBEKF uses the VB iteration method to jointly estimate the MNCM and the target state, which provides a more accurate estimation of the MNCM. Thus, the AI-VBEKF provides the most high-precision tracking trajectory. Therefore, by utilizing the online measurement noise covariance estimating techniques, the SH-EKF and AI-VBEKF can achieve much higher DOA tracking accuracy.

From the RMSE results ($P_0=\mathrm{diag}\left([25,{10}^{-5}]\right)$) shown by

Table 2. The average RMSEs of EKF, SH-EKF and VB-EKF.

	${\mathit{P}}_0$	EKF	SH-EKF	VB-AEKF
$\overline{RMS{E}_{\theta }}$ (°)	$P_0=\mathrm{diag}\left([25,{10}^{-5}]\right)$	7.4	2.8	2.7
	$P_0=\mathrm{diag}\left([50,{10}^{-5}]\right)$	13.2	7.4	2.7
	$P_0=\mathrm{diag}\left([100,{10}^{-5}]\right)$	61.7	7.3	2.7

, similar to the conclusion above, the EKF, SH-EKF and AI-VBEKF converged after a period, and SH-EKF and AI-VBEKF can provide high-precision tracking even with the uncertain measurement noise. The AI-VBEKF can achieve the highest DOA tracking results so that the superior performance of the proposed AI-VBEKF is verified again.

Figure 5c,f and Table 2 ($P_0=\mathrm{diag}\left([50,{10}^{-5}]\right)$,$P_0=\mathrm{diag}\left([100,{10}^{-5}]\right)$) give the DOA tracking results under different initial MSEMs. From the DOA tracking results, it is obvious that the properly set initial MSEM can strongly affect the overall DOA tracking results. By utilizing the attention-based DCNN to estimate the initial MSEM before operating the DOA tracking procedure, the AI-VBEKF can achieve good tracking accuracy and robustness. However, without properly setting the initial MSEM, the DOA tracking results of the SH-EKF and EKF degrade, especially for the case that the initial MSEM is largely deviated from its proper values ($P_0=\mathrm{diag}\left([100,{10}^{-5}]\right)$).

4.2. Experimental Scenario

4.2.1. Experimental Scenario

We used the data from an acoustic experiment in the South China Sea in July 2021 to verify the performance of the proposed method. A ship with an acoustic signal emission system kept moving and emitting an acoustic signal at 170 Hz for 7500 s with a source level of 145 dB during the experimental procedure.

Figure 6a shows the underwater buoy system used in the experiment, which was placed 20 m under the sea surface. A 12-element uniform circular hydrophone array with a radius of 1 m was fixed on this underwater buoy along with a digital acquisition system. As shown in Figure 6b, a compass system (the rectangle) was fixed on the underwater buoy in order to measure the bearing of hydrophone 1 to offset the rotation angle of the buoy φ (measured clockwise wrt north). The digital acquisition system sampled the array signal at a frequency of 8192 Hz during the acoustic signal emission. A CTD profiler was used to measure the sound velocity to be constant at 1544 m/s at a depth from 0 to 20 m.

In order to provide the true bearing angle of the target with regard to the circular hydrophone array for reference, the GPSs are used to record the locations of the acoustic signal emission system and the underwater buoy system. Then, the bearing angle trajectory of the target was calculated clockwise with regard to north by using the GPS data.

4.2.2. Verification of the AI-VBEKF for DOA Tracking

The CBF-based DOA estimation [46], the EKF, SH-EKF [44] and the proposed AI-VBEKF are tested by using the experimental data. The interval of adjacent steps was set to 1 s. The bearing angle of the CBF-based DOA estimation scans from 0° to 360° every 1°. The initial state estimates of the EKF, SH-EKF and the AI-VBEKF are both set to ${\widehat{X}}_0=[{\widehat{\theta }}_{CBF}^0,0]$, where ${\widehat{\theta }}_{CBF}^0$ is the estimate of bearing angle obtained by using CBF-based DOA estimation at k = 0. In addition, to illustrate the superior performance of the proposed AI-VBEKF, $P_k$ is initialized to different initial guesses as $P_0=\mathrm{diag}\left([25,{10}^{-5}]\right)$, $P_0=\mathrm{diag}\left([50,{10}^{-5}]\right)$, and $P_0=\mathrm{diag}\left([100,{10}^{-5}]\right)$. The covariance matrix of the process noise $Q_k$ is set as ${10}^{-4}$. The MNCM was set as $I_{MP}$ for the EKF. The estimation of the MNCM ${\widehat{R}}_k$ is initialized to $I_{MP}$, and the forgetting factor b is set to 0.995 for the SH-EKF. The parameters ${\widehat{u}}_0$, ${\widehat{U}}_0$ and ρ of the AI-VBEKF are set to 16, $16I_{MP}$ and 0.997, respectively.

To numerically compare the tracking precision of the above methods, the bearing angle estimation error (BEE) and the average BEE (ABEE) are defined as the performance metrics, which are given as follows:

$$BEE\left(k\right)=\left|{\widehat{\theta }}_k-{\theta }_k^{GPS}\right|$$

(62)

$$ABEE\left(k\right)=\sqrt{\frac{1}{K}{{\displaystyle \sum _{k=1}^K\left({\widehat{\theta }}_k-{\theta }_k^{GPS}\right)}}^2}$$

(63)

where ${\widehat{\theta }}_k$ is the bearing angle estimate at step k, and ${\theta }_k^{GPS}$ is the true bearing angle obtained by the GPS data.

The trajectories obtained by the AI-VBEKF under different $P_0$ are shown in Figure 7 and

Table 3. ABEEs of different DOA tracking methods under different initial MSEMs.

Data			ABEE (°)
		CBF	EKF	SH-EKF	AI-VBEKF
$P_0=\mathrm{diag}\left([25,{10}^{-5}]\right)$	Raw data	13.2	8.2	7.1	6.3
	With noise	20.1	12.2	7.3	6.5
$P_0=\mathrm{diag}\left([50,{10}^{-5}]\right)$	Raw data	13.2	11.3	11.6	6.8
	With noise	20.1	12.9	11.8	7.2
$P_0=\mathrm{diag}\left([100,{10}^{-5}]\right)$	Raw data	13.2	18.3	14.2	6.4
	With noise	20.1	18.5	14.9	6.6

, which show the BEEs and the ABEEs, respectively. Figure 7 and Table 3 also show the results of the CBF, EKF and the SH-EKF for comparison. Specifically, considering the complicated actual underwater environment, the CBF used 10 measurements to obtain bearing angles every 10 s for a higher precision. In order to further test the robustness of the proposed AI-VBEKF to uncertain measurement noise, a period of Gaussian noise with the covariance of $2I_{MP}$ was added to the raw measurement data from 2600 to 3000 s, and then the methods were tested again. The results are shown in Figure 7 and Table 3.

Figure 7 shows that the trajectories of all methods fluctuated around the true trajectory under the effect of the underwater environmental noise. Figure 7 shows that the trajectories obtained by all methods using raw experimental data fluctuated around the true trajectory due to the effect of the underwater environmental noise. In addition, it can be found that the initial guess of $P_0$ influences the final tracking accuracy. From Figure 7a,b, when $P_0$ is properly set, all the DOA tracking and DOA estimation methods can output reliable results. The trajectory obtained by the CBF-based DOA estimation is the most fluctuant and provided a large error. The precision of the EKF is higher than that of the CBF-based DOA estimation since the EKF considers the kinematic model of the underwater target.

The trajectory obtained by the robust DOA tracking algorithms, namely the SH-EKF and the AI-VBEKF, are closest to the true trajectory. The reason is that the robust DOA tracking algorithms provide more accurate Kalman filter gains than EKF due to the high-precision estimation of the MNCM, which can significantly improve the precision of DOA tracking. In addition, it can be found that the proposed AI-VBEKF has higher DOA tracking accuracy compared the SH-EKF since the VB technique is more analytical, such that it makes a precise mathematical model of the unknown MNCM and makes it more accurately estimated. In addition, the AI-VBEKF utilizes the attention-based DCNN to estimate the initial MSEM, which makes the whole DOA tracking process more stable and accurate. However, the data-inspiring SH-EKF outputs the DOA tracking results with lower accuracy since the precision of the experimental data is not accurate enough such that the tracking window of good estimations is not sufficient. Table 3 shows the high precision and robustness of the AI-VBEKF.

Figure 7c–f and Table 3 illustrate the DOA tracking results under different $P_0$. It can be found that the tracking results of the EKF and SH-EKF degrade when $P_0$ is not set properly. As analyzed in Section 3.3, the high nonlinearities and complexed underwater environment make the unproperly set $P_0$ seriously influence the DOA tracking result. However, since the proposed AI-VBEKF does not rely on the preset $P_0$, its DOA tracking accuracy only has slight changes with different initial MSEMs. As a result, the proposed AI-VBEKF is robust to the initial MSEMs and has the highest tracking accuracy.

Figure 8a–f and Table 3 give the results obtained by using the data with added noise. Similar conclusions can be made as the above. Furthermore, the AI-VBEKF kept robust tracking even during the existence of the added noise under different initial MSEMs, while the performance of the SH-EKF, EKF and the CBF-based DOA estimation obviously degraded when the noise was added or the initial MSEMs were not properly set.

5. Conclusions

An AI-VBEKF robust DOA tracking method was proposed and verified. The theoretical derivations were comprehensively given, and the performances were proven by simulations and the real sea trial data. It can be found in the experiment results that the proposed attention-based DCNN can make good estimation of the initial MSEM after being well trained by the training data with different initial angles, initial errors, and covariance matrix of the measurement noise. Based on the experimental results, the proposed AI-VBEKF can not only have more accurate estimation results, but it is also robust to the uncertain underwater disturbance and the initial guesses of the MSEMs. From the overall consideration of accuracy and robustness, the proposed AI-VBEKF can be considered as an alternative underwater DOA tracking algorithm, especially when the underwater environment is uncertain and the initial MSEM cannot be determined properly.