Active Acoustic Source Tracking Exploiting Particle Filtering and Monte Carlo Tree Search

TL;DR

This paper presents a novel approach for active acoustic source tracking in robotics by combining particle filtering with Monte Carlo Tree Search to plan paths that improve source localization over long-term sequences.

Contribution

It introduces a flexible path planning algorithm that integrates MCTS and particle filtering, optimizing long-term tracking performance in robotic acoustic source localization.

Findings

Enhanced tracking accuracy demonstrated in experiments

Long-term planning outperforms greedy strategies

Effective integration of MCTS with particle filtering

Abstract

In this paper, we address the task of active acoustic source tracking as part of robotic path planning. It denotes the planning of sequences of robotic movements to enhance tracking results of acoustic sources, e.g., talking humans, by fusing observations from multiple positions. Essentially, two strategies are possible: short-term planning, which results in greedy behavior, and long-term planning, which considers a sequence of possible future movements of the robot and the source. Here, we focus on the second method as it might improve tracking performance compared to greedy behavior and propose a flexible path planning algorithm which exploits Monte Carlo Tree Search (MCTS) and particle filtering based on a reward motivated by information-theoretic considerations.

Tables

Table 1. TABLE I: Parameters of the state-space model

		${\sigma }_{q,x}={\sigma }_{q,y}\left[\text{m}\right]$	${\sigma }_{q,v}\left[\frac{\text{m}}{\text{s}}\right]$	${\sigma }_{q,\theta }\left[°\frac{1}{\text{s}}\right]$
Sensor	($q=r$)	$0$	$0$	$5$
Source	($q=s$)	$0$	$0.025$	$10$

Equations

$${\theta }_t^q$$

$$v_t^q$$

$$\left(\begin{array}{c}{x}_t^q\\ y_t^q\end{array}\right)$$

$$\displaystyle p(z_{t}|\boldsymbol{x}_{t})=\int_{z_{t}-\frac{\rho}{2}}^{z_{t}+\frac{\rho}{2}}\mathcal{N}\Big{(}\check{z};$$

$$\left({\Gamma }_{t+1}\right)=\text{SIR}({\Gamma }_t).$$

$$R($$

$$+\sum _{i=1}^I\log \left(p(z_t|{\mathbf{x}}_t^{(i)})\left(\sum _{\nu =1}^{I}p({\mathbf{x}}_t^{(i)}|{\mathbf{x}}_{t-1}^{(\nu )},u_t){\omega }_{t-1}^{(\nu )}\right)\right){\omega }_t^{(i)}$$

$$G({\mathbf{h}}_{t+K})=\sum _{k=1}^K{\gamma }^{k-1}R({\mathbf{h}}_{t+k})$$

$$V({\mathbf{h}}_t,u_{t+1})=\mathbb{E}[G({\mathbf{h}}_{t+K})|u_{t+1},{\mathbf{h}}_t],$$

$$\text{vspace}\ast -.2cm{u}_{t+1}^{\text{opt}}=\underset{u\in \mathcal{U}}{\mathrm{argmax}}V({\mathbf{h}}_t,u).$$

$$u_{t+1}^{\text{UCB}}=\underset{u\in \mathcal{U}}{\mathrm{argmax}}\left(\hat{V}({\mathbf{h}}_t,u)+c\sqrt{\frac{\ln N({\mathbf{h}}_t)}{N_u({\mathbf{h}}_t,u)+ϵ}}\right).$$

Full text

I Introduction

One of the main challenges in modern robotics is intuitive human-robot interaction [1] for which reliable information about the location of the human is essential. If microphones are mounted on the robot, e.g., for the purpose of voice communication, acoustic data can be used for inferring localization information. Typically, the relative Angle of Arrival (AoA) between the robot and the source is estimated by, e.g., MUSIC [2] or SRP-PHAT [3]. Inferring distance information from microphone data is generally difficult in reverberant rooms and often relies on labeled training data, e.g., [4].

As an option to improve localization accuracy one can estimate multiple AoAs at different sensor positions in parallel and fuse them to obtain a position estimate. This can be achieved for example by exploiting several spatially distributed microphone arrays, i.e., an acoustic sensor network [5]. Alternatively, AoAs can also be measured sequentially, e.g., by a microphone array-equipped robot that is moving in a room [6]. The latter approach greatly complicates the task as the robot-source configuration will generally vary over time. Thus, sequential fusing algorithms have to be employed. On the bright side, however, sequential measurements can exploit the controllability of the sensor topology of the robot relative to the target. Thereby, either the movement of the sensor array as a whole, or relative movements of the sensors within the array itself, e.g., by moving head or limbs of the robot, can be exploited such that target tracking can be supported. This idea is referred to as active sensing. For this problem, heuristic strategies have been suggested, e.g., moving a robotic head [7] or its limbs [8] furnished with microphones in the direction of the estimated source. Besides heuristic policies, different algorithms have been proposed to solve this task by maximizing a specific objective function, e.g., [6], [9], [10] or [11]. These approaches can be categorized as greedy or long-term path planning algorithms. Greedy one-step algorithms, as proposed, e.g., in [11] and [12], model the objective to depend only on one future hypothesis. This short-term planning is prone to suboptimal tracking performance as it ignores measurements while moving the sensor and subsequent ones. Thus, long-term motion planning algorithms have been proposed in, e.g., [6] for active localization, i.e., static sources, where the objective function is evaluated over several hypothetical future robot-source configurations. Conditional independence of the objective function given prior movements is assumed to deal with the exponentially increasing amount of possible future hypotheses [6]. In [9], an extension based on MCTS has been proposed which does not rely on the previously assumed conditional independence. Moreover, to circumvent the limitation of discrete source positions [6], the state posterior is modelled by a mixture of Gaussians [9].

In this paper, we introduce a path planning algorithm for active tracking which is also based on MCTS but employs a sequential Monte Carlo method [13], for updating the state posterior, i.e., belief. As sequential Monte Carlo methods approximate the ideal continuous Bayes filter [14], we are capable of overcoming the limitations of Kalman filter based planning algorithms, e.g., the approximation of the non-linearity by a Taylor series [14], and address a broader class of non-linear movement models which are commonly encountered in acoustic source tracking applications [15]. Representing the state posterior by weighted particles also allows to model arbitrary PDFs, e.g., multi-modal PDFs resulting from front-back ambiguity [16]. As planning objective we employ the expected information gain of applying a specific action, based on an estimate of the differential entropy of the belief [17]. Due to the very general description of the state posterior by weighted particles, we are also capable of exploiting higher-order statistics included in the planning objective. Furthermore, compared to [6] and [9], uncertainty of robotic movements is taken into account during planning. This is decisive for real-world applications as the execution of the same control command in various scenarios leads to different movements depending on, e.g., the composition of the ground. Additionally, while [6] and [9] have simulated only static sources, we are focusing on moving sources as a common scenario, e.g., when involving talking humans. Finally, the efficacy of the proposed algorithm and the effect of different planning depths is shown by simulations.

II Acoustic Source Tracking

In the following section, the proposed system is described by a state-space model. Based on this, an appropriate tracking algorithm is discussed.

II-A State-Space Model

The system state at time $t$, captured by the state vector $\boldsymbol{x}_{t}=\begin{pmatrix}(\boldsymbol{x}_{t}^{r})^{\text{T}}&(\boldsymbol{x}_{t}^{s})^{\text{T}}\end{pmatrix}^{\text{T}}$, comprises the robotic - i.e., sensor array - state $\boldsymbol{x}_{t}^{r}$ and the source state $\boldsymbol{x}_{t}^{s}$. For simplicity, the array and the source are assumed to move in the same 2D plane. A state describing either the array $q=r$ or the source $q=s$ is given by $\boldsymbol{x}^{q}_{t}=\begin{pmatrix}x^{q}_{t}&y^{q}_{t}&\theta^{q}_{t}&v^{q}_{t}\end{pmatrix}^{\text{T}}$ including the position in Cartesian coordinates $\begin{pmatrix}x^{q}_{t}&y^{q}_{t}\end{pmatrix}^{\text{T}}$, orientation $\theta^{q}_{t}$ and speed $v^{q}_{t}$. At each time step the system can be influenced by a specific action $u_{q,t}$ controlling the angular speed. Following [15], the system dynamics is modelled by the nonlinear constant-velocity movement model

[TABLE]

Hereby, $w_{t}^{(\cdot)}\sim\mathcal{N}(0,\sigma^{2}_{(\cdot)})$ denotes additive Gaussian noise, parametrized by the variances $\sigma_{q,x}^{2}$, $\sigma_{q,y}^{2}$, $\sigma_{q,v}^{2}$ and $\sigma_{q,\theta}^{2}$, to account for model uncertainties. $\text{wrap}_{[a,b]}(\cdot)$ is a function which wraps its argument to the range $[a,b]$ and $\Delta T$ the sampling interval between successive observations. Note that only the robot can be controlled and thus $u_{s,t}=0\degree\frac{1}{\text{s}}$. The state equations (1) - (3) do not just account for uncertain source movements but also model the array movement to be noisy. Inference about the states is drawn from observations as follows: While we assume the array state to be known, its movement is stochastic. Information about the latent source states is obtained by observing the sound field as sensed by the robot’s microphone array. Hence, under the assumption that the source is continuously active and located in the far-field, an estimate of the AoA ${z}_{t}\in\mathcal{Z}$, e.g., a finite uniform grid of AoA hypotheses with a given angular resolution $\rho$, between the robot and the source can be computed from the microphone measurements at each time step. In this paper, the conditional measurement probability mass function (PMF) $p({z}_{t}|\boldsymbol{x}_{t})$ of observing ${z}_{t}$ while being in state $\boldsymbol{x}_{t}$ is assumed to be a quantized Gaussian

[TABLE]

with mean $\mu_{\text{m}}(\cdot)$ and variance $\sigma^{2}(\cdot)$ being typically non-linear functions of the state. These parameters have to be either modelled or learned. We suggest to learn them in a training phase. In [6], it was assumed that these parameters mainly depend on the distance $d(\boldsymbol{x})$ and the true AoA $\phi_{\text{true}}(\boldsymbol{x})$ between the robot and the source. However, if the robot is equipped with a uniform linear array (ULA), it can additionally be assumed that sources which are located symmetrically with respect to the array axis inherit similar parameters. These assumptions greatly simplify the task of learning the observation model as its mean and variance do not depend on the absolute states of the robot and the source anymore.

II-B Bayesian Tracking

As the states of the system are partially hidden, we introduce the state posterior, i.e., belief, $p(\boldsymbol{x}_{t}|\boldsymbol{h}_{t})$, representing the knowledge about the states, given the vector $\boldsymbol{h}_{t}=\begin{pmatrix}u_{r,1}&z_{1}&\dots&u_{r,t}&z_{t}\end{pmatrix}^{\text{T}}$ of prior actions and observations [14]. This leads to the question how to fuse a new action-observation vector $\begin{pmatrix}u_{r,t+1}&{z}_{t+1}\end{pmatrix}^{\text{T}}$ into a given belief. Straightforward extensions of the Kalman filter update, i.e., extended or unscented Kalman filter [18], might not model accurately the belief due to the possibly multi-modal state posterior and highly nonlinear state transition equations. Thus, we employ a sequential Monte Carlo method [13]. Hereby, the belief is represented by a set $\Gamma_{t}=\{(\boldsymbol{x}_{t}^{(1)},\omega_{t}^{(1)}),\dots,(\boldsymbol{x}_{t}^{(I)},\omega_{t}^{(I)})\}$ of $I$ weighted particles whose elements are tuples of state samples $\boldsymbol{x}_{t}^{(i)}$ and weights $\omega_{t}^{(i)}$. In this paper a Sampling Importance Resampling (SIR) particle filter with systematic resampling is used to update the belief from $t$ to $t+1$ [13, 19, 20] which is abstractly described by

[TABLE]

III Path Planning

Paths are considered in this paper as consecutively executed control commands chosen from a finite set of $J$ discrete commands $u_{r,t}=u_{t}\in\mathcal{U}$. Thus, path planning can be interpreted as a sequential decision making problem, i.e., selecting the optimal action at each time step.

III-A Rewards and Returns

For comparing different paths and actions the concept of rewards and returns, as described in [21], is adopted. Rewards $R(\boldsymbol{h}_{t})$ can be interpreted as instantaneous feedback on the performance of a specific control command $u_{t}$. Ideally one wants to minimize an error norm between the true state and the estimated state. However, the true state is not known. Thus, a common idea is to minimize the belief uncertainty which has to be quantified. [12] proposed the negative differential entropy of the belief as intrinsic measure of uncertainty for active localization. We adopt this idea by utilizing the recursive weighted particle-based entropy estimator proposed in [17]

[TABLE]

corresponding to the history $\boldsymbol{h}_{t}^{\text{T}}=\begin{pmatrix}\boldsymbol{h}_{t-1}^{\text{T}}&u_{t}&z_{t}\end{pmatrix}$.

As rewards can only account for a single action and observation at a time, the concept of returns

[TABLE]

as cumulative discounted, i.e., weighted, future rewards is introduced [21]. Hereby, $\gamma\in[0,1]$ defines the discount factor which controls the behavior of striving for early rewards within the planning horizon of $K$ time intervals. As actions and observations are stochastic, the expected future return

[TABLE]

starting from history $\boldsymbol{h}_{t}$ and selecting action $u_{t+1}$, is introduced [21]. The action maximizing the expected future return is selected for path planning

[TABLE]

III-B Monte Carlo Tree Search Planning

Now we address the problem of how to efficiently compute the expected future returns when applying the actions $u_{t+1}\in\mathcal{U}$: One option is to use Monte Carlo simulation, i.e., simulating $L$ returns corresponding to hypothetical sequences of future actions and observations $\boldsymbol{h}_{t+K}=\begin{pmatrix}\boldsymbol{h}_{t}^{\text{T}}~{}\dots~{}u_{t+K}~{}{z}_{t+K}\end{pmatrix}^{\text{T}}$ followed by arithmetic averaging over all trials. However, this approach suffers from the fact that the number of possible future hypotheses is growing exponentially with increasing planning horizon $K$. To overcome this problem, the MCTS technique [22] has been developed which builds up a partial tree $\mathcal{T}_{\boldsymbol{h}_{t}}$, i.e., a set of simulated future hypotheses, as exemplified in Figure 1 for the action set $\mathcal{U}=\{u^{1},u^{2}\}$ and observation set $\mathcal{Z}=\{z^{1},z^{2},z^{3}\}$, and selects appropriate actions to concentrate on promising paths, i.e., individual sequences with high returns. The exploitation-exploration dilemma, encountered during the sequential decision-making on which paths to investigate, is hereby tackled by the Upper Confidence Bound (UCB) policy [23] formulated as

[TABLE]

Hereby, $N(\boldsymbol{h}_{t})$ denotes the number of previous visits of the simulated hypothesis $\boldsymbol{h}_{t}$, $N_{u}(\boldsymbol{h}_{t},u)$ is the number of previously simulated situations where action $u\in\mathcal{U}$ had been selected and $\epsilon$ is a small number to prevent division by [math]. The ratio of exploration and exploitation can be controlled via the hyperparameter $c>0$ [22].

Algorithm 1, 2 and 3 jointly describe the proposed path planning method, using MCTS and particle filtering. The main Algorithm 1 samples $L$ times from the current belief (line (l.) 2, Alg. 1) and invokes Algorithm 2, i.e., function sim, for updating the estimated future returns $\hat{V}(\boldsymbol{h}_{t},\cdot)$ by simulating possible future hypotheses $\boldsymbol{h}_{t+K}$ (l. 3, Alg. 1), including the computation of the corresponding returns $G(\boldsymbol{h}_{t+K})$ and arithmetic averaging, i.e., (l. 15, Alg. 2). These simulated future hypotheses are depicted in Figure 1 as various paths starting from history $\boldsymbol{h}_{t}$. As only the first new action-observation sequence is added to the tree $\mathcal{T}_{\boldsymbol{h}_{t}}$ (l. 11 Alg. 2), the respective hypotheses have varying lengths. After planning, the robot moves by executing the optimal control command $u_{t+1}^{\text{opt}}$ (l. 5 Alg. 1) and measures the real AoA $\tilde{z}_{t+1}\in\mathcal{Z}$ (l. 6 Alg. 1). If the true hypothesis is not part of the tree, the state posterior is updated (l. 8 Alg. 1). Subsequently, the tree is pruned (l. 9 Alg. 1), i.e., all hypotheses which do not include $\boldsymbol{h}_{t+1}$ are deleted, and the planning starts again.

sim (Alg. 2) is a recursive algorithm which builds up the partial tree $\mathcal{T}_{\boldsymbol{h}_{t}}$ consisting of the respective future movements and observations and computes the estimated expected future returns. If the maximum planning depth $K$ is not yet reached, a future hypothesis is simulated in Alg. 2 by initially selecting an action $u^{\text{UCB}}_{k+1}$ according to the UCB policy (l. 3 Alg. 2). The selection of actions is visualized in Figure 1 by either choosing $u_{t+k}^{1}$ or $u_{t+k}^{2}$. Subsequently, a future state $\boldsymbol{x}^{(i)}_{k+1}$ is sampled (l. 4 Alg. 2) from the PDF described by the state transition equations (1) - (3) conditioned on the sampled particle and selected action. Then, an observation ${z}_{k+1}$ is sampled from the measurement process (l. 5 Alg. 2). Sampling from the measurement process manifests itself in Figure 1 by branching after the simulated actions, i.e., from the empty to the filled circles. Based on the selected action $u^{\text{UCB}}_{k+1}$ and sampled measurement ${z}_{k+1}$, it is checked whether the simulated future hypothesis $\boldsymbol{h}_{k+1}$ is part of the tree $\mathcal{T}_{\boldsymbol{h}_{t}}$. If it is, the algorithm computes the return, i.e., recursive computation of (7), by first computing the reward (6) and afterwards recursively calling itself with the new particle $\boldsymbol{x}_{k+1}^{(i)}$ and incremented depth $k$ (l. 8 Alg. 2). If $\boldsymbol{h}_{k+1}$ has not been simulated before, the state posterior is updated (l. 10 Alg. 2) and $\boldsymbol{h}_{k+1}$ is added to the tree $\mathcal{T}_{\boldsymbol{h}_{t}}$ (l. 11 Alg. 2). Subsequently, the following unexplored part of the tree is investigated by a default algorithm def which is applied until reaching the planning horizon $K$ (l. 12 Alg. 2) [22]. After simulating one future return $\tilde{G}$, $N(\boldsymbol{h}_{t})$ and $N_{\text{UCB}}(\boldsymbol{h}_{t})$ are incremented (l. 13 Alg. 2, l. 14 Alg. 2) and the expected future return $\hat{V}(\boldsymbol{h}_{k},u_{k+1}^{\text{UCB}})$ is updated by computing the mean of the previously simulated returns and the current return $\tilde{G}$ (l. 15 Alg. 2).

The default Algorithm 3 , i.e., function def, consists of computing an action $u_{k+1}^{\text{DEF}}$ according to a default policy (l. 3 Alg. 3), followed by sampling the transition (l. 4 Alg. 3) and observation PDF to create possible future hypotheses (l. 5 Alg. 3) and recursively calling itself (l. 7 Alg. 3) with incremented depth. However, compared to sim the simulated hypothesis is not added to the tree $\mathcal{T}_{\boldsymbol{h}_{t}}$. The default policy cannot rely on experience from previous simulations for action selection [22]. Thus, in [9] uniform random sampling from the action set $\mathcal{U}$ has been used. We propose to incorporate the prior knowledge that moving closer to the target enhances localization results. However, in order to not fully rely on this naive assumption, our proposed default policy randomly uses either a uniform random behavior or an informed action with equal probability. The informed decision is made by sampling from the PDF defined by the array transition equations (1) - (3) with $q=r$ for each possible action and selecting the one moving closest to the estimated source position. Subsequently, the default algorithm recursively calls itself until reaching the planning horizon $K$ and gives back the simulated return $\tilde{G}$.

IV Results

In this section, we evaluate the proposed planning algorithm in a simulated environment. Sensor and source are modelled to move at a constant speed of $v_{t}^{r}=v_{t}^{s}=0.3~{}\frac{\text{m}}{\text{s}}$. The action set $\mathcal{U}$ of the sensor consists of $4$ different actions, namely staying in the current state, or changing the angular velocity $u_{t}$ to $-45~{}\degree\frac{1}{\text{s}}$, $0~{}\degree\frac{1}{\text{s}}$ or $45~{}\degree\frac{1}{\text{s}}$. The respective standard deviations of the sensor and source movement models, i.e., (1) - (3), are summarized in Table I. The control of the angular sensor movements is modelled to be noisy and thus taken into account for path planning. In contrast, the source movement model, i.e., $q=s$, allows for changes of velocity and orientation. The sampling interval $\Delta T$ is set to $1$ s. The particle filter state-space model was chosen according to the simulated model with additional source position uncertainties $\sigma_{s,x}=\sigma_{s,y}=0.1\text{ m}$. The number of particles was set to $I=1000$.

The sensor consists of a ULA with $6$ microphones and a spacing of $4.2$ cm. AoAs are estimated by SRP-PHAT [3] from microphone signals sampled at a sampling frequency of $8$ kHz and observation length of $4096$ samples per measurement. The observation parameters $\mu_{\text{m}}(\cdot,\cdot)$ and $\sigma^{2}_{\text{m}}(\cdot,\cdot)$ are estimated in a training phase by simulating various AoA measurements, i.e., different sensor and source states, for each position of a uniform distance-AoA grid with spacing $0.5$ m and $5\degree$. Positions in between are linearly interpolated. Several rooms with dimensions $7\text{m}\times 5\text{m}\times 3\text{m}$, and varying reverberation time of $0.4$ s, $0.6$ s and $0.8$ s and SNR $10$ dB, $20$ dB and $30$ dB, are simulated [24, 25].

The path planning parameters are set to $\gamma=0.9$, cf. (7), and $c=1.75$, cf. (10) to strive for long-term planning and to rely on previously simulated returns. As suggested in [22], the rewards are normalized to lie in the interval $[-1,0]$. The maximum number of simulated hypotheses $L$ is set to $500$. For reliable results, $300$ source paths and robot starting positions are created with uniformly distributed random initial positions. The initial source belief is modelled by particles whose positions and orientations are distributed uniformly in the enclosure. The first two actions are set to $u_{0}=u_{1}=0\degree\frac{1}{\text{s}}$ to start path planning from an informed belief. The Euclidean distance between the true and the estimated source position, obtained by the weighted mean of the respective particles, is computed for each time step and averaged over all $2700$ combinations of paths and rooms, i.e., reverberation times and SNRs. Figure 2 depicts mean and standard deviation of this error over time, starting from $t=2$, i.e., after the two fixed actions at the beginning, for varying planning depths, i.e., $K=1$, $5$ and $7$. Additionally, as reference movement behaviour we simulated random action selection, similar to [6, 9], and patrolling, e.g., [26], the room by always moving to the most distant unexplored corner.

From Figure 2 it becomes obvious that path planning, either by patrolling or the proposed algorithm, significantly increases the accuracy and reliability of position estimates, compared to selecting an action at random at each time step. Additionally, greedy one-step planning ($K=1$) outperforms the the patrolling reference strategy. The effect of planning depths can be seen by comparing $K=1$ and $K=5$. Both perform similar at the beginning, however drift apart in the long-term. This is consistent with the idea that some actions might pay off later. However, a saturation effect can be seen. Planning $K=7$ time steps ahead does not improve the results further which can be explained either by the increased uncertainty of future predictions or an insufficient number of simulated paths. As the same effect has also been observed for $L=1500$, we suppose this results from the increased uncertainty.

V Summary

In this paper, we proposed a path planning algorithm exploiting MCTS [9] and particle filtering. Based on the belief representation by weighted particles, an estimate of the expected differential entropy has been used as planning objective. Additionally, a new default policy for MCTS has been introduced. The efficacy of the algorithm and the effect of different planning horizons were demonstrated by simulated data of dynamic speech sources for a variety of challenging acoustic scenarios.