Information geometry-based action decision-making for target tracking by fixed-wing unmanned aerial vehicle

Abstract

In this article, we study the ground moving target tracking problem for a fixed-wing unmanned aerial vehicle equipped with a radar. This problem is formulated in a partially observable Markov process framework, which contains the following two parts: in the first part, the unmanned aerial vehicle utilizes the measurements from its radar and employs a Kalman filter to estimate the target’s real-time location; in the second part, the unmanned aerial vehicle optimizes its trajectory in a real-time manner so that the radar’s measurements can include more useful information. To solve the trajectory optimization problem, we proposed an information geometry-based partially observable Markov decision process method. Specifically, the cumulative amount of information in the observation is represented by Fisher information of information geometry, and acts as the criterion of the partially observable Markov decision process problem. Furthermore, to guarantee the real-time performance, an important trade-off between the optimality and computation cost is made by an approximate receding horizon approach. Finally, simulation results corroborate the accuracy and time-efficiency of our proposed method and also show our advantage in computation time compared to existing methods.

Introduction

Guiding unmanned aerial vehicle (UAV) to detect and track a suspicious ground target is an important requirement in many intelligence, surveillance, target acquisition, and reconnaissance (ISTAR) ¹ problems. Different from predefined trajectory tracking problems, tracking a moving target is more challenging, since it requires the UAV to respond promptly to the random movement of the target. Our article focuses on the problem of a UAV tracking a moving ground target in an uncertain environment. In many practical scenarios, the accurate observation and tracking of a target is achieved through adequate maneuver of the UAV equipped with a sensor. Therefore, the UAV needs to make a movement strategy based on the state of the target. However, the information of the target location obtained by the UAV sensor is usually incomplete or imperfect. In this regard, the task of target tracking mainly consists of two aspects ² : one is the estimation of the state of the target based on the measurements obtained by the sensor, and the other is the adjustment of the position/pose of the UAV based on the prediction of the target state to obtain better measurements.

For target state estimation, different strategies have been developed in the past decade. In order to estimate and predict the target state more accurately, the classical methods such as Kalman filter (KF), ^{3 –7} particle filter, ^{8 –11} and their modifications were widely used. For instance, an unscented KF was utilized to track an underwater submarine target/moving ship. ^6,12 Tang and Ozguner proposed the particle-filter-and-hospitability-map algorithm (PF-HMap) ⁸ to deal with the general target tracking maintenance problem with regional and intermittent measurements. In these abovementioned methods, the classical KF is one of the most widely used methods for estimation and tracking, due to its optimality, simplicity, and tractability.

For position/pose adjustment, one possible method is to track the estimated location of the target. A novel algorithm combining the tangent vector field guidance (TVFG) path-planning approach and the Lyapunov vector field guidance (LVFG) algorithm ¹³ was developed. Given the target position and current UAV dynamic state, this method is theoretically possible to obtain the shortest path with UAV operational constraints. In addition, based on the division of two kinds of possible path parttens, that is, the z type (sinusoidal type) and the whirling type, the ground target pursuit algorithm ^14,15 generates waypoints step by step and steers the UAV to the latest waypoint. Most of these planning methods are based on precomputed vector fields or alternative paths, which do not respond well to system uncertainties.

The other method for position/pose adjustment is based on the decision-making framework, which can better handle the system uncertainties. In this framework, the transition model of the tracking system is usually assumed to be Markovian. Then, the tracking policy selection problem can be formulated as a partially observed Markov decision process (POMDP), in which the state is only partially observed and one seeks to design a control policy which maps the state probability distributions to actions and maximizes the accuracy of target location. With its decision-making ability in an uncertain environment, ^16,17 the POMDP framework has been widely used in a variety of real-world scenarios. Prentice and Roy ¹⁸ addressed the problem of trajectory planning with imperfect state information, and model it as a linear-Guassian POMDP. Assuming that the UAV’s state and kinematics are known, Ragi and Chong ^3,4 employed the POMDP, combining with the motion constraints of the UAV, to design the guidance algorithm of UAV tracking a ground target. However, solving POMDP optimally is proven to be PSPACE-hard. ¹⁹ A large number of literature focuses on various heuristic or approximate solution techniques. ²⁰ Some literature use the receding horizon theory to approximate a POMDP. Sunberg et al. ²¹ proposed a receding horizon control approach to solve the information space dynamic programming of POMDPs. In the POMDP framework, the tracking decision problem was transformed into an optimization problem, and the choice of the optimization criterion plays a crucial role, ²² which directly determines the speed and performance of the solution. Different from the approach ^3,4 that used the trace of covariance matrix, we use the cumulative quantity of information from information geometry (IG) as the reward criterion of POMDP to evaluate the performance of the tracking strategy in essence and simplify the calculation.

IG is proposed by Rao, ²³ and axiomatized by Chentsov. ²⁴ It constructs a distance representation between statistical distributions. It deals with families of parameterized probability densities which carry a metric structure. This structure is derived via the well-known Fisher information metric. Recently, significant attention has been drawn in the area of IG. Costa et al. ²⁵ presented the Fisher information distance as a measure of difference between two probability distribution functions. The Fisher information distance is related to the information of the target estimation. Therefore, the target tracking problem is transformed into searching for the strategy that maximizes the accuracy of the target estimation. There have been many works on how IG methodology can be applied to the tracking process as it does in the signal processing. ^{26 –29} In their study, Akselrod et al. ²⁷ proposed a MDP model to formulate the collaborative sensor management for multi-target tracking decision processes, where the objective function is based on the Fisher information measure. However, in that study, the Fisher information is merely approximated by the posterior covariance matrix from the KF, which is not an exact derivation. Actually, it did not analyze the source and significance of the Fisher information in essence. In target tracking, the noise covariance depends on the performance of the sensor. That means the covariance matrix of noise matches to the Riemannian manifold of positive definite matrix. Wang et al. ²⁹ derived a complete form of relationship between the accuracy of the target location and the optimal sensor position, based on the maximum determinant of the Fisher information matrix (FIM) criterion, where the theory of IG to study the problem of bearing-only tracking was employed.

The above summarizes the research results about the construction of the tracking decision-making framework, the selection of the evaluation function, and the solution of the optimal decision. However, even though IG is a power tool in analyzing the statistics of the moving target, few studies provide rigorous analysis and design of the IG-based UAV ground target tracking problem. Previous work by Zhao et al. ³⁰ has shown the effectiveness of the IG method in the target tracking decision in two-dimensional (2-D) space by simulation results. In this article, we extend the problem of target tracking with an airborne radar to general three dimensions with new simulation results and more detailed analysis. The iterative form of FIM about the 3-D bearing-and-range radar, which takes into account the uncertainties of predicted target states, is derived. Then, the Fisher information distance on the Riemann manifold is regarded as the basis for POMDP, and the convergence of the POMDP approximate solution algorithm is proved theoretically. The detailed contributions of the proposed algorithm are as follows:

A novel IG-based POMDP frame is provided for the UAV to track a moving ground target in a 3-D uncertain environment. In this frame, the state of the target is observed by the 3-D bearing-and-range radar.

The optimization criterion based on the 3-D observation model in the POMDP is derived based on Fisher information in the view of the IG, which is the key to optimize the action polices in the UAV-to-target problem.

An approximate receding horizon control is developed to obtain an acceptable control strategy in the trade-off between the optimization and the computation cost. We also prove the convergence of approximation algorithm theoretically.

The rest of this article is organized as follows. In the “Problem formulation” section, the system model is given and the target tracking problem is formulated. The framework of the POMDP and its criterion of accumulative information are presented in “Target tracking decision-making based on IG” section, where the FIM is derived by iterative calculation via the predicted target state. In “Approximate receding horizon approach for POMDP” section, we introduce our proposed approximate receding horizon approach for the POMDP, and analyze the performance of the algorithm. “Simulation results” section presents the simulation and results, which is followed by the conclusion.

Problem formulation

In this article, the target tracking system is composed of a UAV and a target vehicle. The mission of the UAV is to observe and track the moving target on the ground. The UAV is equipped with a radar that can obtain the bearing and range measurements of the tracked object with limited precision and reliability.

Fixed-wing UAV model

The fixed-wing UAV dynamics augmented by the autopilot is a high dimensional, highly nonlinear, and extremely complex system. In our work, the UAV is supposed to fly at a constant altitude h, and the unicycle model is adopted to describe kinematics of fixed-wing UAV.

The UAV autopilot controls bank angle ϕ and forward acceleration a directly. The UAV state, which is a part of the world states, is ${\mathit{s}}^u={\left(x^u,y^u,v^u,{\theta }^u\right)}^{\text{T}}$ , where $\left(x^u,y^u\right)$ denotes the position of the UAV in the inertial coordinate system (CS) and $\left(v^u,{\theta }^u\right)$ represents the UAV’s speed and course angle. The UAV is assumed to be fully observable, and the discrete kinematic model of the UAV is as follows

$\{\begin{array}{l}{x}_{k+1}^u=x_k^u+v_k^{u}T\text{cos}{\theta }_k^u,\hfill \\ y_{k+1}^u=y_k^u+v_k^{u}T\text{sin}{\theta }_k^u,\hfill \\ v_{k+1}^u={\left[v_k^u+a_{k}T\right]}_{V_{\text{min}}}^{V_{\text{max}}},\hfill \\ {\theta }_{k+1}^u={\theta }_k^u+\left(gT\text{tan}({\varphi }_k)/v_k^u\right)\hfill \end{array}$

1

where ${\left[v\right]}_{V_{\text{min}}}^{V_{\text{max}}}=\text{max}\left\{V_{\text{min}},\text{min}(V_{\text{max}},v)\right\}$ , V_max and V_min are the upper and lower bounds of the ground speed of the fixed-wing UAV, g is the acceleration of gravity, and T is time step.

Target model

In this system, the target is on the ground and its exact state is not available. The state is given by ${\mathit{s}}^t={\left(x^t,y^t,{\dot{x}}^t,{\dot{y}}^t\right)}^{\text{T}}$ , including the location $\left(x^t,y^t\right)$ and velocity $\left({\dot{x}}^t,{\dot{y}}^t\right)$ . The mobility of the ground target is modeled as follows, which is commonly used in the literature

${\mathit{s}}_{k+1}^t=A{\mathit{s}}_k^t+{\mathbf{\nu }}_k,{\mathbf{\nu }}_k\sim \mathcal{N}(0,Q_k)$

2

where A is the state transition matrix with the following form

$A=\left[\begin{array}{cccc}1& 0& T& 0\\ 0& 1& 0& T\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right]$

3

and Q _k represents the process noise covariance matrix at time k.

The states of the target are measured by the active radar mounted on the UAV. The radar can measure the range and bearing information determined by relative state ${\mathit{s}}^o={\left(x,y,\dot{x},\dot{y}\right)}^{\text{T}}$ , which is represented by $\left(x,y\right)=\left(x^t-x^u,y^t-y^u,{\dot{x}}^t-{\dot{x}}^u,{\dot{y}}^t-{\dot{y}}^u\right)$ . Next, it is necessary to illustrate the observation model.

Observation model

In this work, the UAV only installs one radar as the measuring sensor. The airborne radar for target tracking provides measurement of the target in the sensor CS, which is a polar CS with range r and bearing φ. The measurement model with noise at time k is

${\mathit{o}}_k=\left[\begin{array}{c}{r}_k\\ {\phi }_k\end{array}\right]=\mathit{\eta }({\mathbf{s}}_k^o)+{\mathit{\omega }}_k$

4

where η (⋅) is the observation function. Then, $\mathit{\eta }({\mathbf{\text{s}}}_k^o)$ indicates the true position of the target in the polar CS of the sensor and ω _k is the measurement error. It is usually assumed that the measurement error follows a zero-mean normal distribution in the sensor CS, that is

${\mathit{o}}_k|{\mathit{s}}_k^{\mathit{o}}\mathrm{.}\sim N(\mathit{\eta }({\mathit{s}}_k^{\mathit{o}}),\mathbf{C}({\mathit{s}}_k^{\mathit{o}}))$

For the radar, the error-free target position in the sensor polar CS is represented as

$\begin{array}{c}\mathit{\eta }({\mathit{s}}_k^{\mathit{o}})=\left[\begin{array}{c}{\eta }_1({\mathit{s}}_k^{\mathit{o}})\\ {\eta }_2({\mathit{s}}_k^{\mathit{o}})\end{array}\right]\\ =\left[\begin{array}{c}\sqrt{x_k^2+y_k^2+h^2}\\ \text{arctan}\left(y_k/x_k\right)\end{array}\right]\end{array}$

5

As for an active sensor, the measurement error is dependent on the signal-to-noise ratio which is proportional to the fourth power of distance. Then the covariance matrix of the measurement is

$\mathbf{C}({\mathit{s}}_k^{\mathit{o}})=\left[\begin{array}{cc}{r}_k^4{\sigma }_r^2& 0\\ 0& {\sigma }_{\phi }^2\end{array}\right]$

6

where σ_r and σ_φ are the standard deviations of the range and bearing measurements, respectively. In equation (6), $r_k^4$ indicates the effect of the fourth power decay of the echo signal strength on the variance of the distance measurement. Strictly speaking, the power decay should be related to the actual distance rather than the distance measurement. Since the actual distance is not available at the UAV, we use distance measurement r_k instead.

This study considers a UAV decision-making problem, in which the goal is to design an algorithm to control a UAV for target tracking. Specifically, the UAV motion model is simplified as equation (1). It is mounted with a sensor that measures the relative position of the target. The observation model (equation (4)), combined with the assumed target model (equation (2)), is used to estimate and predict the states of the target. Based on these states the UAV makes action decision. The values of actions are limited to the maximum and minimum range. The objective is to obtain better observation for a more accurate estimation of the target state.

Tracking in Cartesian coordinate

In the target tracking problem, the target’s movement is described in Cartesian coordinate, while the measurement is available physically in the sensor polar CS. Thus, it is necessary to convert the measurements from the polar CS to the Cartesian coordinate. Specifically, in the Cartesian CS (selecting the UAV as the origin), the measurement model is converted to the following form

${{\mathit{o}}^{\prime }}_k=\mathbf{H}{\mathit{s}}_k^{\mathit{o}}+{{\mathit{\omega }}^{\prime }}_k$

7

where the measurement matrix is

$\mathbf{H}=\left[\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\end{array}\right]$

8

In equation (7), o ′ _k is the relative target position measured by the sensor, and ω ′ _k stands for the resulting measurement error.

Defining φ (⋅) as the transition function from polar CS to Cartesian CS, that is

$\left[\begin{array}{c}x\\ y\end{array}\right]=\mathit{\phi }\left(\left[\begin{array}{c}r\\ \phi \end{array}\right]\right)=\left[\begin{array}{c}r\text{cos}\phi \\ r\text{sin}\phi \end{array}\right]$

where (x, y) and (r, φ) are the coordinates in the polar CS and the Cartesian CS, respectively. Then, the relative position from the target to the UAV in Cartesian coordinate, defined as $\mathbf{H}{\mathit{s}}_k^{\mathit{o}}$ , has the following form

$\mathbf{H}{\mathit{s}}_k^{\mathit{o}}=\mathit{\phi }({\mathit{o}}_k-{\mathit{\omega }}_k)$

9

By Taylor series expansion of $\mathit{\phi }({\mathit{o}}_k-{\mathit{\omega }}_k)$ around the noisy measurement o _k , we obtain

$\mathbf{H}{\mathit{s}}_k^{\mathit{o}}=\mathit{\phi }({\mathit{o}}_k)-\mathbf{J}({\mathit{o}}_k){\mathit{\omega }}_k+\mathcal{O}({\mathit{\omega }}_k)$

10

where $\mathcal{O}({\mathit{\omega }}_k)$ stands for the higher order $(\ge 2)$ terms, and the Jacobian matrix $\mathbf{J}({\mathit{o}}_k)$ is evaluated at the noisy measurement o _k

$\mathbf{J}({\mathit{o}}_k)=\frac{\partial \mathit{\phi }}{\partial {\mathit{o}}_k}=\frac{h}{r_k}\left[\begin{array}{cc}\text{cos}({\phi }_k)& -r_k\text{sin}({\phi }_k)\\ \text{sin}({\phi }_k)& r_k\text{cos}({\phi }_k)\end{array}\right]$

11

Considering $\mathit{\phi }({\mathit{o}}_k)={{\mathit{o}}^{\prime }}_k$ and the equation (10), the exactly converted measurement model (equation (7)) can be written as

${{\mathit{o}}^{\prime }}_k=\mathbf{H}{\mathit{s}}_k^{\mathit{o}}+\underset{\mathit{\omega }{\text{'}}_k}{\underbrace{\mathbf{J}({\mathit{o}}_k){\mathit{\omega }}_k-\mathcal{O}({\mathit{\omega }}_k)}}$

12

The standard approach usually ignores $\mathcal{O}({\mathit{\omega }}_k)$ in the exact model (equation (12)), and treats ω ′ _k approximately as zero-mean with covariance

$\mathbf{R}({\mathit{s}}_k^{\mathit{o}})=\mathbf{J}({\mathit{o}}_k)\mathbf{C}({\mathit{s}}_k^{\mathit{o}})\mathbf{J}{({\mathit{o}}_k)}^{\text{T}}$

13

Target tracking decision-making based on IG

The decision-making problem of UAV under uncertain environment is modeled as a POMDP. In the POMDP framework, the decision-making program obtains measurements from the airborne radar and makes the action decisions for the UAV. For more measurement information of the target, we should predict the states of the target and make appropriate action decisions based on a certain criterion. In this article, we represent the reward criterion by the accumulative information, which is derived from the FIM in the IG. This section introduces the POMDP-based target tracking algorithm, especially the calculation of the decision-making criterion. In the second part of this section, we describe the traditional criterion and analyze its shortcomings. Then we propose the decision-making criterion based on the IG in the last part of this section. The key point is that we use the iterative calculation form of the FIM for the range-and-bearing radar to redefine the decision-making criterion in essence.

The framework of POMDP-based decision-making

POMDP is a stochastic process controlled by a decision-maker. In the target tracking problem, since the system state transition is a random process with the Markov property, and the states of the target cannot be obtained directly, we use the POMDP to select a sequence of actions for the UAV to reduce the uncertainty of target localization. In general, an infinite horizon POMDP is defined by a tuple $M=〈\mathit{S},\mathit{U},\mathit{T},\mathit{\Omega },\mathit{O},\mathit{R},b_0〉$ . The following defines the key components of the POMDP with respect to a target tracking problem.

States. The state set S is a discrete Borel space and the state at time k, including the state of the UAV ${\mathit{s}}_k^u$ and the state of the target ${\mathit{s}}_k^t$ , and the calculated relative state of the target to the UAV ${\mathit{s}}_k^{\mathit{o}}$ .

Actions. In this system, the action of the UAV is the control quantity ${\mathit{u}}_k\in \mathit{U}$ . Specifically, the action decision at time k is given by ${\mathit{u}}_k=[a_k,{\varphi }_k{]}^{\text{T}}$ .

State transition probabilities. T is a set of conditional transition probabilities between states. The UAV state transition function is the kinematic equation (equation (1)). The target state is unknown and approximated by a constant velocity model (equation (2)).

Observations and observation probabilities. Ω is the set of possible observation, and O is the set of conditional probabilities of the observation when an action is performed. In the tracking system, the UAV is assumed to be fully observable. We use a radar to measure the target state. The observation law in the sensor CS is equation (4) and the linearized conversion in the Cartesian CS is equation (7). We use the latter, when we estimate the current target state based on historical data. However, when we predict the accuracy of target location, we use the former, which is the research subject of this article.

Belief state. Belief state is the distribution of the real state. Specially, the state distribution at time k = 0 is b₀. Since the UAV states are fully observable, the belief state is $b_k({\mathit{s}}^{\mathit{u}})=\delta \left({\mathit{s}}^{\mathit{u}}-{\mathit{s}}_k^{\mathit{u}}\right)$ . At the same time, the probability distribution of target positions is assumed to be Gaussian, and can be expressed approximately as $b_k({\mathit{s}}^t)\sim N\left({\mathit{\mu }}_k,{\mathbf{P}}_k\right)$ , where μ _k and P _k are derived from the standard KF. If an optimal action policy exists in POMDP, there is an optimal action sequence that only depends on the belief state feedback. ³¹ Then we can make decision according to the belief state.

Reward function. The real-valued reward function R defines the reward of the action. It is used to compare different alternative action policies. In the target tracking problem, the core of evaluation is the accuracy of the target location. Some literatures describe it as the trace of covariance. In our work, for the purpose of better reflecting the essential characteristics of radar measurement, we propose the IG method to represent the reward function.

Remark 1

Partially observable Markov process. At each time period, the system is in some state s ∈ S . The agent takes an action u ∈ U, which causes the transition of the system to state s ′ with probability $\mathit{T}({\mathit{s}}^{\prime }|\mathit{s},\mathit{u})$ . At the same time, the agent receives an observation o ∈ Ω which depends on the new state of the system with probability $\mathit{O}(\mathit{o}|{\mathit{s}}^{\prime },\mathit{u})$ . Finally, the agent receives a reward equal to $R(\mathit{s},\mathit{u})$ . Then the process repeats.

The goal of the POMDP for the agent is to choose actions at each time step that maximize its long-run average expected reward

$J_{\infty }^{\pi }\left(\mathit{s}\right)=\underset{H\to \infty }{\text{lim inf}}\frac{1}{H}\text{E}\left[{\displaystyle \sum _{k=0}^{\mathbf{H}-1}}R({\mathit{s}}_k,{\mathit{u}}_k)|{\mathit{s}}_0=\mathit{s}\right]$

14

where H is the time horizon. In our target tracking system, the aim of decision-making is to find the optimal control policies for the UAV. This strategy is to make the UAV better able to observe and track the target, so as to maximize the accuracy of the target position estimation.

A policy is a sequence $\pi =\left\{{\pi }_k\right\}$ of probability π_k on the action space given the state (belief state) and action history satisfying the constraint. We denote by Δ the set of all policies. Once we have a POMDP model and a set Δ of admissible policies, we need to express a performance index or objective function. It is a function to assess the system’s performance when a given policy π ∈ Δ is used and the initial state (or belief state) of the system is s (or b ). Then we are able to select the optimal strategy. Unlike the actual state, the belief state is completely observable. It is necessary to represent the expected reward based on the belief state b

$R^{\prime }\left(b,\mathit{u}\right)={\displaystyle \int }R\left(\mathit{s},\mathit{u}\right)b(\mathit{s})\text{d}\mathit{s}$

15

This reward $R^{\prime }\left(b,\mathit{u}\right)$ then represents the expected reward function of POMDP. The objective is transformed into

$J_{\infty }^{\pi }\left(b\right)=\underset{H\to \infty }{\text{lim inf}}\frac{1}{H}\left[{\displaystyle \sum _{k=0}^{H-1}}{R}^{\prime }(b_k,{\mathit{u}}_k)|b_0=b\right]$

16

according to the nominal belief state optimization. ⁴ The objective function is approximated as follows

$J_{\infty }^{\pi }\left(b\right)\approx \underset{H\to \infty }{\text{lim inf}}\frac{1}{H}\left[{\displaystyle \sum _{k=0}^{H-1}}R({\widehat{\mathit{s}}}_k,{\mathit{u}}_k)|{\widehat{\mathit{s}}}_{\mathit{o}}={\displaystyle \int }\mathit{s}b(\mathit{s})\text{d}\mathit{s}\right]$

17

where ${\widehat{\mathit{s}}}_1,{\widehat{\mathit{s}}}_2\mathrm{,...}$ is a nominal belief state sequence. The nominal target belief state sequence can be obtained from the target state transfer law (equation (2)) with exactly zero-noise sequence.

Covariance-based conventional criterion function

In some literatures, ^3,4 the belief of the target mentioned above can be identified with the state of the tracker $({\mathit{\mu }}_k,{\mathit{P}}_k)$ , which are the posteriori state estimation and the posteriori error covariance matrix, and are given by

$\{\begin{array}{l}{\mu }_{k+1}=\mathbf{A}{\mu }_k\hfill \\ P_{k+1|k}=\mathbf{A}{P}_k{\mathbf{A}}^{\text{T}}+Q_k\hfill \\ P_{k+1}={\left[P_{k+1|k}^{-1}+H_{k+1}^{\text{T}}{R_{k+1}}^{-1}{H}_{k+1}\right]}^{-1}\hfill \end{array}$

18

The reward function is defined to represent the uncertainty of the target location, which is usually represented by the mean-squared error between the tracker and the target. Then, the objective function is as follows

$J_{\infty }^{\pi }\left(b\right)\approx \underset{H\to \infty }{\text{lim inf}}\frac{1}{H}\left[{\displaystyle \sum _{k=0}^{H-1}}(-\text{tr}{P}_{k+1})\right]$

19

The Kalman filtering equation is a linearized conversion of the actual system model. The measurement model in the sensor CS is equation (4), and the observation model in the Cartesian CS is in the form of equation (7). In contrast with the long-run reward from the trace of the Kalman filtering error covariance matrix (equation (19)), the FIM in the IG processes the original measurement data of the radar. It is derived from the measurement model in the natural sensor CS directly and has a clearly physical meaning, which reflects the volume of information from the measurement data. The greater the cumulative information is, the more accurate the measurement would be. Hence the Fisher information from the IG can better evaluate the accuracy of the predicted states. We use the Fisher information distance (which will be explained later) on the statistical manifold as the basis of the POMDP.

Cumulative information in IG

IG ²³ offers a comprehensive result about statistical models by regarding them as geometrical objects. From the perspective of IG, the set of belief states forms the statistical manifold of a particular geometric structure. The information is defined on the Fisher information divergence between the current belief state and the belief state after a measurement has been made. The cumulative information which the sensor may acquire from the target is characterized by the Fisher information distance. In a discrete measurement sampling scenario, the sum of the determinant of FIM is used to approximate the information distance. ³²

In this article, the determinant of the FIM is used to characterize the volume of information obtained by measurements. The performance index is the expected long-run average reward, therefore, it is a feasible way to predict the information each time through iterations. The following introduces the iteration process of the discrete FIM of the radar measurements.

We define G _k to express the FIM at time k, then the reward of the POMDP is $R\left({\mathit{s}}_k,{\mathit{u}}_k\right)=|{\mathit{G}}_{k+1}|$ , where $|\cdot |$ denote the determinant of the matrix. When we predict the volume of information in the coming period of time, the FIM at time k is defined ^5,23 by

${\mathit{G}}_k=\text{E}\left[({\nabla }_{\mathit{s}}\text{ln}p({\mathit{o}}_{1:k}|{\mathit{s}}_k)){({\nabla }_{\mathit{s}}\text{ln}p({\mathit{o}}_{1:k}|{\mathit{s}}_k))}^{\text{T}}\right]$

20

where $p\left({\mathit{o}}_{1:k}|{\mathit{s}}_k\right)$ is the batch measurement likelihood. In this tracking scenario, the batch measurement likelihood is defined as equation (21)

$\begin{array}{l}p\left({\mathit{o}}_{1:k}|{\mathit{s}}_k\right)={\displaystyle \prod _{m=1}^k}\frac{1}{\sqrt{2\pi |\mathit{C}({\mathit{s}}_k,m)|}}\\ \text{exp}\left(-\frac{1}{2}{({\mathit{o}}_m-\mathit{\eta }({\mathit{s}}_k,m))}^{\text{T}}{\mathit{C}}^{-1}({\mathit{s}}_k,m)({\mathit{o}}_m-\mathit{\eta }({\mathit{s}}_k,m))\right)\end{array}$

21

where $\mathit{\eta }({\mathit{s}}_k,m)$ denotes error-free target position of the sensor polar, evaluated at time m in terms of the state s _k at k. That is to say, $\mathit{\eta }({\mathit{s}}_k,m),m=1,2,\dots ,k$ is given by

$\begin{array}{c}\mathit{\eta }({\mathit{s}}_k,m)=\left[\begin{array}{c}{\eta }_1({\mathit{s}}_k,m)\\ {\eta }_2({\mathit{s}}_k,m)\end{array}\right]\\ =\left[\begin{array}{c}\sqrt{\Delta {x_{k,m}}^2+\Delta {y_{k,m}}^2+h^2}\\ \text{arctan}\left(\frac{\Delta x_{k,m}}{\Delta y_{k,m}}\right)\end{array}\right]\end{array}$

22

Besides, $\mathbf{C}({\mathit{s}}_k,m)$ is the measurement noise covariance evaluated at time m in terms of the state s _k at k, which is given by

$\mathbf{C}({\mathit{s}}_k,m)=\left[\begin{array}{cc}{r}_{k,m}^4{\sigma }_r^2& 0\\ 0& {\sigma }_{\phi }^2\end{array}\right]$

23

and

$r_{k,m}^2=\Delta {x_{k,m}}^2+\Delta {y_{k,m}}^2+h^2$

In the above equations, considering the linear motion target model described previously, we have

$\{\begin{array}{c}\Delta x_{k,m}=x_k^t-(k-m)T{\dot{x}}_k^t-x_m^{\mathit{u}}\\ \Delta y_{k,m}=y_k^t-(k-m)T{\dot{y}}_k^t-y_m^{\mathit{u}}\end{array}$

24

Then the FIM at time k can be calculated in a recursive form as follows

${\mathbf{G}}_k={\mathbf{G}}_{k-1}+\Delta {\mathbf{G}}_k$

25

where

$\begin{array}{l}{\left[\Delta {\mathbf{G}}_k\right]}_{ij}={\left[\frac{\partial \mathit{\eta }({\mathit{s}}_k,k)}{\partial {\mathit{s}}_i}\right]}^{\text{T}}{\mathbf{C}}^{-1}\left({\mathit{s}}_k,k\right)\left[\frac{\partial \mathit{\eta }\left({\mathit{s}}_k,k\right)}{\partial {\mathit{s}}_j}\right]\\ +\frac{1}{2}tr\left({\mathbf{C}}^{-1}\left({\mathit{s}}_k,k\right)\frac{\partial \mathbf{C}\left({\mathit{s}}_k,k\right)}{\partial {\mathit{s}}_i}{\mathbf{C}}^{-1}\left({\mathit{s}}_k,k\right)\frac{\partial \mathbf{C}\left({\mathit{s}}_k,k\right)}{\partial {\mathit{s}}_j}\right)\end{array}$

26

In equation (26), ${\left[\Delta {\mathbf{G}}_k\right]}_{ij}$ represents the i th row and the j th column of elements in the matrix. The specific derivation process is given in the Online Supplementary material.

Since the sensor collects measurements only at discrete points, the accumulative information should only consider those points when the sensor takes measurements. It is assumed that the decision-making points are consistent with the measurement points. Therefore, the sum of determinant of the FIM ²⁹ is used to approximate the accumulative information $({\mathit{s}}_1^{\mathit{u}},{\mathit{s}}_2^{\mathit{u}}\mathrm{,...})$ , that is

$\mathcal{D}({\mathit{s}}_1,{\mathit{s}}_2\mathrm{,...})\approx {\displaystyle \sum _{k=1}^{\infty }}|{\mathbf{G}}_k|$

27

The reward function based on the IG is

$J_{\infty }^{\pi }\left(b\right)\approx \underset{H\to \infty }{\text{lim inf}}\frac{1}{H}\left[{\displaystyle \sum _{k=0}^{\mathbf{H}-1}}|{\mathbf{G}}_{k+1}|\right]$

28

The larger the volume of information is, the more accurate the estimation of the target state will be and the better the tracking performance the algorithm will have. Our decision-making algorithm aims to find the optimal sequence of actions which maximizes the predicted cumulative information. From equation (20), we obtain that the FIM is determined by the relative position between the UAV and the target. However, the states of the target are unobservable. Thus, when we make decision, we use the predicted belief states to present the reward of the POMDP, that is, the predicted determinant of the FIM. That is to say, we use

${\mathit{\mu }}_{k+1}=\mathbf{A}{\mathit{\mu }}_k$

to predict the belief state of the target, and $\mathit{\mu }={\displaystyle \int }\mathit{s}b(\mathit{s})\text{d}\mathit{s}$ is used to replace the real state.

Approximate receding horizon approach for POMDP

Since it is intractable to solve the policy for maximizing the objective function (equation (17)) exactly, we use an fixed finite-horizon POMDP to create a policy to solve the infinite-horizon POMDP approximatively in this article, which is called “approximate receding horizon approach.”

Algorithm design

The idea of the approximate receding horizon approach is that at the current time t, we obtain the optimal policy sequence of the POMDP over a finite horizon $\left[t,t+H\right]$ , and at each time t, we predict the future states over the horizon length H based on the current state of the system. The finite horizon H indicates the number of predicted steps. Once an optimal sequence is found, the first action is applied as a control command to the UAV and the whole process repeats. We can take the process at time 0 as an example for analysis. That is, the current state is assumed to be s ₀. Then the problem boils down to find the sequence of actions over a time horizon H to maximize the expected cumulative reward.

The specific algorithm is as follows.

Algorithm 1.

The approximate receding horizon algorithm

Algorithm analysis

In this section, we analyze the performance of the approximate receding horizon approach for infinite-horizon average reward.

When we say the algorithm is stable, we mean that the difference between the value function and the optimal value function is bounded. Then the proposed algorithm is stable and the following theorem can be obtained.

Theorem 1

By defining the receding H-horizon control policy as ${\pi }_H\in \Delta$ with $H<\infty$ , we have

$\underset{b}{\text{sup}}|V_{\infty }^{{\pi }_H}-V_{\infty }^{*}|\le \frac{\Vert R^{\prime }\Vert }{1-\alpha }\cdot {\alpha }^n$

29

Proof

Recall that the state of the UAV is observable and controllable, and the state of the target obeys normal distribution.

The expectation of the system state at time k is

$\text{E}({\mathit{s}}_k)=({\mathit{s}}_k^{\mathit{u}},{\mathit{\mu }}_k)$

30

Then we have $\text{E}({\mathit{s}}_k)\in \mathit{S}$ , since the real state space S is a Borel space.

Define $\widehat{\mathit{s}}=\text{E}(\mathit{s})$ , $\mathit{K}=\left\{(\widehat{\mathit{s}},\mathit{u}(\widehat{\mathit{s}}))|\widehat{\mathit{s}}\in S,\mathit{u}\in U\right\}$ . Since the state space is continuous and compact, the system satisfies the ergodicity conditions. In other words, there exists a state expectation ${\widehat{\mathit{s}}}^{*}\in S$ such that the transition probability $p\left({\widehat{\mathit{s}}}^{*}|\mathit{k}\right)=1$ for all k ∈ K .

The ergodicity condition can be expressed in another way. ³³ That is, there exists a constant α < 1 such that

$\underset{k,k^{\prime }}{\text{sup}}\Vert p(\cdot |\mathit{k})-p(\cdot |{\mathit{k}}^{\prime })\Vert \le 2\alpha$

where the sup is over all k and k ′ in K , and ∥⋅∥ denotes the variation norm for signed measures.

Recall the one-step reward $R^{\prime }\left(b_k,{\mathit{u}}_k\right)$ of the POMDP based on the IG. The reward is the accumulative information when the system state transforms from s _k to s _k+1.

Thus we have

$R^{\prime }\left({\mathit{b}}_k,{\mathit{u}}_k\right)=|{\mathbf{G}}_{k+1}|$

31

where ${\mathbf{G}}_{k+1}$ is deduced from equation (25). The determinant of the FIM is:

$|{\mathbf{G}}_{k+1}|=\frac{8r_{k+1}^2{\sigma }_r^2+1}{r_{k+1}^6{\sigma }_r^2{\sigma }_{\phi }^2}$

32

The smallest distance between the target and the UAV is the height of the UAV which is h. From this it can be seen that the reward is bounded, that is

$|R^{\prime }|\le \frac{8h^2{\sigma }_r^2+1}{h^6{\sigma }_r^2{\sigma }_{\phi }^2}$

The above specification shows that the system satisfies the condition under which Hernández-Lerma ³³ show that

$0\le V_{\infty }^{*}-V_{\infty }^{{\pi }_H}\le \frac{\Vert R^{\prime }\Vert }{1-\alpha }\cdot {\alpha }^n$

33

We can see that the receding horizon approach provides a good approximation for the optimal infinite-horizon average reward and the error approaches zero geometrically with α.

The decision-making of UAV needs high real-time performance, so we tend to choose algorithm which is lower in the complexity of computation time. Next, the article analyzes the complexity of the computation time of the objective function in our algorithm (equation (28)), and compares it with that based on KF (equation (19)).

Using the reward derived from the IG to calculate the objective function, we should obtain G _k from iterative computations of equations (39) to (41) in the Online Supplementary material. In the same way, using the reward derived from the KF to compute the objective function, we need to obtain P _k from equation (18) and compute the measurement covariance R _k in the Cartesian CS. According to the above analysis, the required operations for the IG-based objective function and the KF-based objective function during time horizon H are listed in Table 1. In these operations, we ignore the calculation of state transitions of the UAV and the target, which are the same in both algorithms.

Table 1.

Required operations for the IG-based objective function and the KF-based objective function.

	Addition	Multiplication	Other operations
IG	30H−2	59H	0
KF	273H−1	635H	1H sin, 1H cos

IG: information geometry; KF: Kalman filter.

It is clearly shown in Table 1 that the required operations, either addition and multiplication or other operations of the IG-based objective function (equation (28)), are much less than those of the KF-based objective function (equation (19)), in other words, it can greatly improve computational efficiency to use information accumulation as the objective function of the optimal decision-making.

Remark 2

Several notes on operations. (1) The operations of matrix inversion in equation (18) are counted based on the adjoint matrix inversion method. (2) In the following simulation, the MATLAB command fmincon is used to minimize the objective function. (3) Due to the randomness of the tracking system, the difference of initial conditions, and the limitation of iteration steps, the computation time of optimal strategy is not linear with the computational complexity of the objective function and the tracking error also has a certain randomness.

Simulation results

This section presents the simulations to verify the performance of the proposed algorithm, and compares the simulation results with other method. In the comparative algorithm, the trace of the covariance in the KF is used as the objective function.

The simulation is implemented in MATLAB R2016b (maci64), where the action decision is obtained from the MATLAB command $f\mathit{\text{min}}c\mathit{o}n$ . The speed of the UAV is limited between 25 m/s and 35 m/s, and the fixed flight height is set to 100 m. To corroborate the effectiveness of our algorithm for different target movements and sensor parameters, we design four groups of typical simulations first.

In these four experiments, the time horizon H is set to 6, ³ which means that in every time step, the UAV plans 6 future time steps in advance, but only the first action is implemented. In each simulation group, we carry out 200 Monte Carlo simulation experiments. Each experiment performs 100 steps and takes 0.5 s per step. In the simulations, the root-mean-square error (RMSE) between the target and the tracker is employed to measure the location estimation accuracy for different algorithms (our algorithm and a KF-based algorithm; the KF-based algorithm means the performance criterion is based on the trace of the corresponding covariance matrix). Additionally, the computation time is taken into account for comparing the real-time performance of different algorithms.

In the first two groups of simulation experiments, the target is subjected to rectilinear motion, and the average linear velocity of the target is 15 m/s. The initial position of the UAV is (0, 50), and the initial position of the target is (0, 0). Figure 1 shows the trajectories of the UAV for target tracking.

Figure 1.

The trajectories of a UAV tracking a rectilinear motion target. UAV: unmanned aerial vehicle.

In the first simulation group, the sensor measurement standard deviation parameters are ${\sigma }_r{=10}^{-4}/\text{m}$ , ${\sigma }_{\phi }=5\times {10}^{-3}\pi \text{rad}$ . Figure 2(a) and (b) compares the two algorithms in RMSE and computation time, respectively. From Figure 2, we can see that our algorithm has a very close performance compared to the KF-based algorithm in the sense of location RMSE. Nevertheless, our algorithm is much faster than the KF-based algorithm, and the computation time of our algorithm is usually half than that in the KF-based algorithm.

Figure 2.

The performance of group I: (a) the RMSE of tracking and (b) the computation time. RMSE: root-mean-square error.

In the second simulation group, the sensor measurement standard deviation parameters are ${\sigma }_r=8\times {10}^{-4}/\text{m}$ , ${\sigma }_{\phi }{=10}^{-2}\pi \text{rad}$ . Figure 3 shows the performance of the tracking. From Figure 3(a), we can see that our algorithm gradually has a better RMSE performance; and after the 80th step, our algorithm outperforms the KF-based algorithm. From Figure 3(b), we can see that the computation time for our algorithm is usually half of that in the KF-based algorithm.

Figure 3.

The performance of group II: (a) the RMSE of tracking and (b) the computation time. RMSE: root-mean-square error.

For the first two simulation groups, when the sensor measurement noise is larger, our proposed method performs better. Moreover, the decision-making takes only half of the computation time which is required by the KF-based algorithm, which is significant in the practical application of UAVs.

In the next two simulation groups, the target moves in a circle, and the angular velocity of the target is 0.1 rad/s, and the radius of motion is 200 m. The UAV’s initial position is (0, 100), and the initial position of the target is (0, 200). Figure 4 shows the trajectories of the UAV and the target.

Figure 4.

The trajectories of a UAV tracking a circular motion target. UAV: unmanned aerial vehicle.

In the third group of simulations, the sensor measurement standard deviation parameters are set to ${\sigma }_r{=10}^{-4}\text{m}$ , ${\sigma }_{\phi }=5\times {10}^{-3}\pi \text{rad}$ . The following also illustrates the simulation results of tracking a target moving circularly in Figure 5(a) and (b).

Figure 5.

The performance of group III: (a) the RMSE of tracking and (b) the computation time. RMSE: root-mean-square error.

In the fourth group of simulations, the sensor measurement standard deviation parameters in group III are ${\sigma }_r=8\times {10}^{-4}\text{m}$ , ${\sigma }_{\phi }{=10}^{-2}\pi \text{rad}$ . Figure 6(a) and (b) compares the two algorithms in RMSE of tracking and average computation time.

Figure 6.

The performance of group IV: (a) the RMSE of tracking and (b) the computation time. RMSE: root-mean-square error.

As it shows in the last two scenarios, our method produces better tracking performance and is more time-efficient in solving.

The simulations indicate that the IG-based algorithm can obtain strategies faster in different sensor measurement precision and target movement forms. As shown in the previous section where the computational complexity of the objective function is analyzed, the decision-making algorithm based on the IG can greatly save the computation time. Besides, when the observation error increases, the tracking accuracy decreases.

In addition, we also present partial numerical results. For each set of experiments, we calculate the average distance error and average computation time. We summate every step distance error over the simulation runtime, and the mean of these errors (from each Monte Carlo run) is called the average tracking error. Similarly, for every step of the simulation runtime, we record the average of computation time, which is called the average computation time. The data statistics are listed in Table 2, and the better data are bolded. Obviously, our algorithm is more time-efficient in solving in each simulation. It is clear that the calculations of inverse matrix in KF cost more time than our method. It is also important to emphasize that the decision time for each step is shorter than the predefined time (0.5 s) per step. That is to say, the algorithm can meet the requirement of real-time calculation. However, in view that IG method takes shorter for decision-making, there is a greater degree of pre-improvement of the performance of decision-making by the short decision period. At the same time, we notice that the tracking error of the IG-based algorithm does not always perform better than the KF-based algorithm. The main reason is that, in this discrete measurement sampling scenario, the information distance is approximated by the sum of the determinant of FIM, which sometimes reduces the estimation accuracy for some scenarios. However, as shown in Figure 5(a) and (b) and Table 2, the IG-based algorithm performs better for high-order motion system, such as the circular motion target. Nevertheless, the IG-based method is more time-efficient than the KF-based method, which is more suitable for practical applications, especially when the computational resource is limited.

Table 2.

Statistical results of average computation time and average tracking error.

Group	Computation Time		Tracking error		Group	Computation time		Tracking error
	KF	IG	KF	IG		KF	IG	KF	IG
Group I: rectilinear motion, ${\sigma }_r{=10}^{-4}/\text{m}$ , ${\sigma }_{\phi }=5\times {10}^{-3}\pi \text{rad}$ .	0.2663	0.1397	0.9396	2.1303	Group III: circular motion, ${\sigma }_r{=10}^{-4}/\text{m}$ , ${\sigma }_{\phi }=5\times {10}^{-3}\pi \text{rad}$ .	0.2772	0.1421	5.1901	5.8166
	0.2502	0.1383	1.4695	1.5368		0.2536	0.1653	4.7843	5.4185
	0.2461	0.1149	1.5470	1.6527		0.2524	0.1295	5.2428	4.2371
	0.2377	0.1144	1.6872	1.8515		0.2386	0.1312	5.1264	4.8794
	0.2237	0.1175	1.0837	2.5165		0.2471	0.1326	4.5247	5.0728
	0.2258	0.1326	1.3376	1.7106		0.2642	0.1289	5.6331	5.4632
	0.2329	0.1194	1.2946	1.6813		0.2613	0.1171	5.0692	5.4448
	0.2264	0.1119	1.3023	2.0482		0.2497	0.1198	5.1437	4.5239
	0.2314	0.1030	0.9930	1.1756		0.2338	0.0910	5.3624	5.4249
	0.2310	0.1203	1.2539	1.2345		0.2273	0.1090	5.1987	5.2313
Group II: rectilinear motion, ${\sigma }_r=8\times {10}^{-4}/\text{m}$ , ${\sigma }_{\phi }{=10}^{-2}\pi \text{rad}$ .	0.2758	0.1956	12.1068	20.2923	Group IV: circular motion, ${\sigma }_r=8\times {10}^{-4}/\text{m}$ , ${\sigma }_{\phi }{=10}^{-2}\pi \text{rad}$ .	0.2805	0.1775	31.8603	25.7812
	0.2415	0.1720	16.1810	89.6004		0.2595	0.1934	30.2535	25.0373
	0.2483	0.1628	17.9907	24.6986		0.2580	0.1735	22.3440	29.4811
	0.2353	0.1500	18.2345	23.9607		0.2574	0.2027	34.4343	31.0254
	0.2318	0.1740	19.6474	23.1907		0.2360	0.1651	24.2302	40.2282
	0.2367	0.1693	11.5808	29.1510		0.2461	0.1739	26.1126	27.7403
	0.2255	0.1738	8.8416	25.7922		0.2497	0.1649	23.2740	22.4048
	0.2594	0.1649	12.5577	21.2700		0.2426	0.1734	28.6919	25.3280
	0.2889	0.1641	30.6833	28.5152		0.2452	0.1635	30.4491	17.5406
	0.2686	0.1654	19.0282	23.3245		0.2522	0.1677	25.4797	26.2637

KF: Kalman filter; IG: information geometry.

In order to study the effect of different prediction horizons on computation time and tracking performance, we also carry out two groups of experiments. In these two groups, the movement of the target is different and each group including 200 times Monte Carlo simulations at different time horizons. The average computation time and the location RMSE of tracking for these 200 experiments with different time horizons are presented in Figures 7 and 8.

Figure 7.

The UAV tracks the target of linear motion with different time horizons, (a) the RMSE of tracking and (b) the computation time. RMSE: root-mean-square error; UAV: unmanned aerial vehicle.

Figure 8.

The UAV tracks the target of circular motion with different time horizons (a) the RMSE of tracking and (b) the computation time. RMSE: root-mean-square error; UAV: unmanned aerial vehicle.

These figures reflect that the computation time reduces and the tracking performance gets worse with the time horizon decreasing. Therefore, when we select the time horizon, it need to make a trade-off between the computation load and the tracking performance.

Conclusion and future work

In this article, we have studied the moving ground target tracking problem of a fixed-wing UAV. More specifically, a POMDP-based action decision-making method has been proposed, which gives the optimal sequence to maximize the target information observed by the radar. In this method, we have introduced the FIM as the criterion of the proposed method with the aid of IG. Simulation results corroborate the effectiveness of our proposed method; and show that compared to the classical KF-based method, our method has higher time-efficiency in computations. In our future work, we would extend our method to multi-target tracking scenarios.

Supplemental material

Supplemental_material - Information geometry-based action decision-making for target tracking by fixed-wing unmanned aerial vehicle: From algorithm design to theory analysis

Supplemental_material for Information geometry-based action decision-making for target tracking by fixed-wing unmanned aerial vehicle: From algorithm design to theory analysis by Yunyun Zhao, Xiangke Wang, Yirui Cong, and Lincheng Shen in International Journal of Advanced Robotic Systems