Physeter catodon localization by sparse coding

arXiv:1306.3058v1 [cs.LG] 13 Jun 2013

Sébastien PARIS
DYNI team, LSIS CNRS UMR 7296, Aix-Marseille University

sebastien.paris@lsis.org

Yann DOH
DYNI team, LSIS CNRS UMR 7296, Université Sud Toulon-Var

yanndoh.m2@gmail.com

Hervé GLOTIN
DYNI team, LSIS CNRS UMR 7296, Université Sud Toulon-Var

glotin@univ-tln.fr

Xanadu HALKIAS
DYNI team, LSIS CNRS UMR 7296, Université Sud Toulon-Var

halkias@univ-tln.fr

Joseph RAZIK
DYNI team, LSIS CNRS UMR 7296, Université Sud Toulon-Var

razik@univ-tln.fr

Abstract
This paper presents a spermwhale’ localization architecture using jointly a bag-offeatures (BoF) approach and machine learning framework. BoF methods are known, especially in computer vision, to produce from
a collection of local features a global representation invariant to principal signal transformations. Our idea is to regress supervisely from these local features two rough estimates of the distance and azimuth thanks
to some datasets where both acoustic events
and ground-truth position are now available. Furthermore, these estimates can feed
a particle filter system in order to obtain a
precise spermwhale’ position even in monohydrophone configuration. Anti-collision system and whale watching are considered applications of this work.

1. Introduction
Most of efficient cetacean localisation systems are
based on the Time Delay Of Arrival (TDOA) estimation from detected1 animal’s click/whistles signals
1

As click/whistles detector, matching filter is often prefered
Proceedings of the 30 th International Conference on Machine Learning, Atlanta, Georgia, USA, 2013. JMLR:
W&CP volume 28. Copyright 2013 by the author(s).

(Nosal & Frazer, 2006; Bénard & Glotin, 2009). Longbase hydrophones’array is involving several fixed, efficient but expensive hydrophones (Giraudet & Glotin,
2006) while short-base version is requiring a precise array’s self-localization to deliver accurate results. Recently (see (Glotin et al., 2011)), based on Leroy’s attenuation model versus frequencies (Leroy, 1965), a
range estimator have been proposed. This approach
is working on the detected most powerful pulse inside the click signal and is delivering a rough range’
estimate robust to head orientation variation of the
animal. Our purpose is to use i) these hydrophone’
array measurements recorded in diversified sea conditions and ii) the associated ground-truth trajectories of
spermwhale (obtained by precise TDAO and/or Dtag
systems) to regress both position and azimuth of the
animal from a third-party hydrophone2 (typically onboard, standalone and cheap model).
We claim, as in computer-vision field, that BoF approach can be successfully applied to extract a global
and invariant representation of click’s signals. Basically, the pipeline of BoF approach is composed of
three parts: i) a local features extractor, ii) a local feature encoder (given a dictionary pre-trained on
data) and iii) a pooler aggregating local representations into a more robust global one. Several choice for
encoding local patches have been developed in recent
years: from hard-assignment to the closest dictionary
basis (trained for example by Kmeans algorithm) to
2
We assume that the velocity vector is colinear with the
head’s angle.

Physeter catodon localization by sparse coding

a sparse local patch reconstruction (involving for example Orthognal Maching Pursuit (OMP) or LASSO
algorithms).

2. Global feature extraction by spare
coding
2.1. Local patch extraction
Let’s denote by C , {C j }, j = 1, . . . , H the collection
of detected clicks associated with the j th hydrophone
of the array composed by H hydrophones. Each matrix C j is defined by C j , {cji }, i = 1, . . . , N j where
cji ∈ Rn is the ith click of the j th hydrophone. For
our Bahamas2 dataset (Giraudet & Glotin, 2006), we
choose typically n = 2000 samples surrounding the
detected click. The total number of available clicks is
H
P
equal to N =
Nj.
i=1

As local features, we extract simply some local signal patches of p ≤ n samples (typically p = 128)
and denoted by z ji,l ∈ Rp . Furthermore all z ji,l are
`2 normalized. For each cji , a total of L local patches
n
e samZ ji , {z ji,l }, l = 1, . . . , L equally spaced of d L
ples are retrieved (see Fig. 1). All local patches associated with the j th hydrophone is denoted by Z j ,
{Z ji }, i = 1, . . . , N j while Z , {Z j } is denoting all
the local patches matrix for all hydrophones. A final
post-processing consists in uncorrelate local features
by PCA training and projection with p0 ≤ p dimensions.
2.2. Local feature encoding by sparse coding
In order to obtain a global robust representation of
c ⊂ C, each associated local patch z ⊂ Z are first
linearly encoded via the vector α ∈ Rk such as z ≈
Dα where D , [d1 , . . . , dk ] ∈ Rp×k is a pre-trained
dictionary matrix whose column vectors respect the
constraint dTj dj = 1. In a first attempt to solve this
linear problem, α can be the solution of the Ordinary
Least Square (OLS) problem:


1
2
lOLS (α|z; D) , min
kz − Dαk2 .
(1)
α∈Rk 2
OLS formulation can be extended to include regularization term avoiding overfitting. We obtain the ridge
regression (RID) formulation:


1
kz − Dαk22 + βkαk22 . (2)
lRID (α|z; D) , min
α∈Rk 2
This problem have an analytic solution α = (D T D +
βI k )−1 D T z. Thanks to semi-positivity of D T D +

βI k , we can use a cholesky factor on this matrix to
solve efficiently this linear system. In order to decrease
reconstruction error and to have a sparse solution, this
problem can be reformuled as a constrained Quadratic
Problem (QP):
1
lSC (α|z; D) , min kz − Dαk22 s.t. kαk1 = 1. (3)
α∈Rk 2
To solve this problem, we can use a QP solver involving
high combinatorial computation to find the solution.
Under RIP assumptions (Tibshirani, 1994), a greedy
approach can be used efficiently to solve and eq. 3 and
this latter can be rewritten as:
1
lSC (α|z; D) , min kz − Dαk22 + λkαk1 ,
k
α∈R 2

(4)

where λ is a regularization parameter which controls
the level of sparsity. This problem is also known as
basis pursuit (Chen et al., 1998) or the Lasso (Tibshirani, 1994). To solve this problem, we can use the
popular Least angle regression (LARS) algorithm.
2.3. Pooling local codes
The objective of pooling (Boureau et al.; Feng et al.)
is to transform the joint feature representation into
a new, more usable one that preserves important information while discarding irrelevant detail. For each
click signal, we usually compute L codes denoted V ,
{αi }, i = 1, . . . , L. Let define v j ∈ RL , j = 1, . . . , k as
the j th row vector of V . It is essential to use feature
pooling to map the response vector v j into a statistic value f (v j ) from some spatial pooling operation f .
We use v j , the response vector, to summarize the joint
distribution of the j th compounds of local features over
the region of interest (ROI). We will consider the `µ norm pooling and defined by:

fn (v; µ) =

L
X

! µ1
|vm |µ

s.t. µ 6= 0.

(5)

m=1

The parameter µ determines the selection policy for
locations. When µ = 1, `µ -norm pooling is equivalent to sum-pooling and aggregates the responses over
the entire region uniformly. When µ increases, `µ norm pooling approaches max-pooling. We can note
the value of µ tunes the pooling operation to transit
from sum-pooling to max-pooling.
2.4. Pooling codes over a temporal pyramid
In computer vision, Spatial Pyramid Matching (SPM)
is a technic (introduced by (Lazebnik et al.)) which
improves classification accuracy by performing a more

Physeter catodon localization by sparse coding
0.25

0.2

20

0.15

40

0.1

60

0.05

80
0

100
−0.05

120
−0.1
0

500

1000

1500

2000

100

200

300

400

500

600

700

800

900

1000

Figure 1. Left: Example of detected click with n = 2000. Right: extracted local features with p = 128, L = 1000 (one
local feature per column).

robust local analysis. We will adopt the same strategy
in order to pool sparse codes over a temporal pyramid
(TP) dividing each click signal into ROI of different
sizes and locations. Our TP is defined by the matrix
Λ of size (P × 3) (Paris et al.):

risk RM :
RM (V , D) ,

M
1 X1
kz i − Dαi k22 + λkαi k1
M i=1 2

(7)

s.t. dTj dj = 1.
Λ = [a, b, Ω],

(6)

where a, b, Ω are 3 (P × 1) vectors representing subdivision ratio, overlapping ratio and weights respectively. P designs the number of layers in the pyramid.
Each row of Λ represents a temporal layer of the pyramid, i.e. indicates how do divide the entire signal into
sub-regions possibly overlapping. For the ith layer, the
i
click signal is divided into Di = b 1−a
bi +1c ROIs where
ai , bi are the ith elements of vector a, b respectively.
P
P
For the entiere TP, we obtain a total of D =
Di
i=1

ROIs. Each click signal c (n × 1) is divided into temporal ROI Ri,j , i = 1, . . . , P , j = 1, . . . , Di of size
(bai .nc × 1). All ROIs of the ith layer have the same
weight Ωi . For the ith layer,
by bbi .nc
 ROIs are shifted

1 1 1
samples. A TP with Λ = 1 1
is designing a
1
2
4
2-layers pyramid with D = 1 + 4 ROIs, the entiere signal for the first layer and 4 half-windows of n2 samples
with 25% of overlapping for the second layer. At the
end of pooling stage over Λ, the global feature x ∈ Rd ,
d = D.k is defined by the weighted concatenation (by
factor Ωi ) of L pooled codes associated with c.
2.5. Dictionary learning
To encode each local features by sparse coding (see
eq. 4), a dictionary D is trained offline with an important collection of M ≤ N.L local features as input. One would minimize the regularized empirical

Unfortunatly, this problem is not jointly convex but
can be optimized by alternating method:
RM (V |D̂) ,

M
1 X1
kz i − D̂αi k22 + λkαi k1 ,
M i=1 2

(8)

which can be solved in parallel by LASSO/LARS and
then:
M
1 X1
kz i − D α̂i k22 s.t. dTj dj = 1.
M i=1 2
(9)
Eq. 9 have an analytic solution involving a large matrix (k × k) inversion and a large memory occupation
for storing the matrix V (k × M ). Since M is potentially very large (up to 1 million), an online method to
update dictionary learning is prefered (Mairal et al.).
Figure 2 depicts 3 dictionary basis vectors learned via
sparse coding. As depicted, some elements reprensents
more impulsive responses while some more harmonic
responses.

RM (D|V̂ ) ,

3. Range and azimuth logistic
regression from global features
After the pooling stage, we extracted unsupervisly N
global features X , {xi } ∈ Rd×N . We propose to
regress via logistic regression both range r and azimuth az (in x − y plan, when animal reach surface
to breath) from the animal trajectory groundtruth denoted y. For the current train/test splitsets of the

Physeter catodon localization by sparse coding
4

x 10

0.1

−1

0

−1.1

−0.1
H

−1.2

−0.2
0

20

40

60

80

100

120

8

−1.3

0.1

−1.4

m

0

7

−1.5

H

−0.1
−1.6

9

H

−0.2
0

20

40

60

80

100

120

−1.7

0.2

H10

−1.8

0
−1.9

−0.2
−0.4

−2
0.7

0

20

40

60

80

100

H11
0.8

0.9

1

1.1
m

120

Figure 2. Example of trained dictionary basis with sparse
coding.

S
S
data, such as X = X train X test , y = y train y test
and N = Ntrain + Ntest , ∀ {xi , yi } ∈ X train × y train ,
we minimize:
(
)
NX
train
1 T
−yi w T
x
i
θ
b θ = arg min
w
log(1 + e
w w +C
) ,
wθ 2 θ θ
i=1
(10)
where yi denotes ri and azi for θ = r and θ = az respectively. Eq. 10 can be efficiently solved for example
with Liblinear software (Fan et al., 2008). In the test
part, range and azimuth for any xi ∈ X test are reconb Tr xi and by a
b Taz xi
tructed linearly by rbi = w
cz i = w
respectively.

1.2

1.3

1.4

1.5
4

x 10

Figure 3. The 2D trajectory (in xy plan) of the single
sperm whale observed during 25 min and corresponding
hydrophones positions.

l
l
l
represent the ground truth,
where yi,j
, ybi,j
and Ntest
the estimate and the number of test samples for the lth
hydrophone respectively. The global ARMSE is then
H
P
1
calculated by ARM SE = H
ARM SE(l).
l=1

4.2. `µ -norm pooling case study
For prilimary results, we investigate the influence of
the µ parameter during the pooling stage. We fix the
number of dictionary basis to k = 128 and the temporal pyramid equal to Λ1 = [1, 1, 1], i.e. we pool sparse
codes on whole the temporal click signal. A value of

4. Experimental results
ARMSE on range

4.1. bahamas2 dataset

1900

To extract local features, we chose n = 2000, p = 128
and L = 1000 (tuned by model selection). For both
the dictionary learning and the local features encoding, we chose λ = 0.2 and fixed 15 iterations to train
dictionary on a subset of M = 400.000 local features drawn uniformaly. We performed K = 10 crossvalidation where training sets reprensented 70% of the
total of extracted global features, the rest for the testing sets. Logistic regression parameter C is tuned by
model selection. We compute the average root mean
square error (ARMSE) of range/azimuth
estimates per
s
hydrophone: ARM SE(l) =

1
K

K
P

l
NP
test

i=1

j=1

l −y
l )2
(yi,j
bi,j

Λ1

1800
1700

ARMSE (m)

This dataset (Giraudet & Glotin, 2006) contains a total of N = 6134 detected clicks for H = 5 different
hydrophones (named H 7 , H 8 , H 9 , H 10 and H 11 and
with N 7 = 1205, N 8 = 1238, N 9 = 1241, N 10 = 1261
and N 11 = 1189 respectively).

1600
1500
1400
1300
1200
1100
0.25

0.5

1

1.5

2

2.5
µ

3

4

6

8

10

20

Figure 4. ARM SE vs. µ for range estimation.

µ = {3, 4} seems to be a good choice for this pooling
procedure. For µ ≥ 20, results are similar to those obtained by max-pooling. For azimuth, we observe also
the same range of µ values.

Physeter catodon localization by sparse coding

4.3. Range and azimuth regression results
Here, we fixed the value of µ = 3 and we varied the
number of dictionary basis k from 128 to 4096 elements. We also investigated the influence of the temporal pyramid and we give results fortwo particu
lary choices: Λ1 = [1, 1, 1] and Λ2 =

1

1

1
3

1
3

1
1

.

For Λ2 , the sparse are first pooled over all the signal
then pooled over 3 non-overlapping windows for a total of 1 + 3 = 4 ROIs. In order to compare results
of our presented method, we also give results for an
hand-craft feature (Glotin et al., 2011) specialized for
spermwhales and based on the spectrum of the most
energetic pulse dtected inside the click. This specialized feature, denoted Spectrum feature, is a 128 points
vector.
ARMSE on range
Λ1

1200

Λ2
Spectrum Feat

ARMSE (m)

1100

1000

900

temporal pyramid for pooling permits also to improve
slightly results.

5. Conclusions and perspectives
We introduced in the paper, for spermwhale localization, a BoF approach via sparse coding delivering
rough estimates of range and azimuth of the animal,
specificaly towarded for mono-hydrophone configuration. Our proposed method works directly on the
click signal without any prior pulses detection/analysis
while being robust to signal transformation issue by
the propagation. Coupled with non-linear filtering
such as particle filtering (Arulampalam et al., 2002),
accurate animal position estimation could be perform
even in mono-hydrophone configuration. Applications
for anti-collision system and whale whatching are targeted with this work.
As perspective, we plan to investigate other local features such as spectral features, MFCC (Davis & Mermelstein, 1980; Rabiner & Juang, 1993), Scattering
transform features (Andén & Mallat). These latter
can be considered as a hand-craft first layer of a deep
learning architecture with 2 layers.

800

References
700

600
128

256

512

1024
k

2048

3072

4096

Figure 5. ARM SE vs. k for range estimation with µ = 3.

Andén, Joakim and Mallat, Stéphane. Multiscale scattering for audio classification. In ISMIR, 11.
Arulampalam, M. Sanjeev, Maskell, Simon, and Gordon, Neil. A tutorial on particle filters for online
nonlinear/non-gaussian bayesian tracking. IEEE
Trans. SP, 50:174–188, 2002.

ARMSE on azimuth

Bénard, Frédéric and Glotin, Hervé. Whales localization using a large array : performance relative to cramer-rao bounds and confidence regions.
In e-Business and Telecommunications, pp. 294–
306. Springer - Verlag, Berlin Heidelberg, september
2009.

120
115
110

ARMSE (deg)

105
100
95
90

Boureau, Y-Lan, Ponce, Jean, and Lecun, Yann.
A theoretical analysis of feature pooling in visual
recognition. In ICML’ 10.

85
Λ1

80

Λ2

75

Spectrum Feat
128

256

512

1024
k

2048

3072

4096

Figure 6. ARM SE vs. k for azimuth estimation with µ =
3.

For both range and azimuth estimate, from k = 2048,
our method outperforms results of the Spectrum feature and particulary for azimuth estimate. Using a

Chen, Scott Shaobing, Donoho, David L., Michael, and
Saunders, A. Atomic decomposition by basis pursuit. SIAM Journal on Scientific Computing, 20:
33–61, 1998.
Davis, S. and Mermelstein, P. Comparison of parametric representations for monosyllabic word recognition in continuously spoken sentences. IEEE Trans.
ASSP, 28:357–366, 1980.

Physeter catodon localization by sparse coding

Fan, Rong-En, Chang, Kai-Wei, Hsieh, Cho-Jui,
Wang, Xiang-Rui, and Lin, Chih-Jen. LIBLINEAR:
A library for large linear classification. JMLR, 2008.
Feng, Jiashi, Ni, Bingbing, Tian, Qi, and Yan,
Shuicheng. Geometric lp -norm feature pooling for
image classification. In CVPR ’11.
Giraudet, Pascale and Glotin, Hervé. Real-time 3d
tracking of whales by precise and echo-robust tdoas
of clicks extracted from 5 bottom-mounted hydrophones records of the autec. Applied Acoustics,
67:1106–1117, 2006.
Glotin, H., Doh, Y., Abeille, R., and Monnin, A. Physeter distance estimation using sub-band leroy transmission loss model. In 5th Internationnal Workshop
on Detection, Classification, Localization and Density Estimation of Marine Mammals using Passive
Acoustics, 2011.
Lazebnik, Svetlana, Schmid, Cordelia, and Ponce,
Jean. Beyond bags of features: Spatial pyramid
matching for recognizing natural scene categories.
In CVPR ’06, pp. 2169–2178.
Leroy, C. Sound attenuation between 200 and 10000
cps mesured along single paths. Technical Report 43, Saclant ASW Research Center, 1965.
Mairal, Julien, Bach, Francis, Ponce, Jean, and Sapiro,
Guillermo. Online dictionary learning for sparse
coding. In ICML ’09.
Nosal, E.-M. and Frazer, L. Track of a sperm whale
from delays between direct and surface-reflected
clicks. Applied Acoustics, 67:1187–1201, 2006.
Paris, Sébastien, Halkias, Xanadu, and Glotin, Hervé.
Efficient bag of scenes analysis for image categorization. In ICPRAM’ 13.
Rabiner, L. and Juang, B.H. Fundamentals of Speech
Recognition. Prentice Hall PTR, 1993.
Tibshirani, Robert. Regression shrinkage and selection via the lasso. Journal of the Royal Statistical
Society, Series B, 58:267–288, 1994.

