On considère la diagonalisation conjointe de tenseurs complexes d'ordre trois symétriques. Nous proposons un algorithme de type Jacobi pour la réaliser et la solution optimale dans le cas 2 × 2 est donnée. Nous illustrons le lien entre cette décomposition conjointe et la séparation de sources non circulaires par l'intermédiaire d'une fonction de contraste. Des simulations illustrent l'approche proposée

DE LUIGI, Christophe

MOREAU, Eric

I-Revues

Diagonalisation conjointe de tenseurs d’ordre 3
Application a` la se´paration de signaux non-circulaires
Christophe De Luigi et Eric Moreau
LSEET, UMR CNRS 6017, Universite´ du Sud Toulon Var
avenue G. Pompidou, BP56, 83162 La Valette du Var Cedex
deluigi@univ-tln.fr, moreau@univ-tln.fr
Re´sume´ – On conside`re la diagonalisation conjointe de tenseurs complexes d’ordre trois syme´triques. Nous proposons un
algorithme de type Jacobi pour la re´aliser et la solution optimale dans le cas 2× 2 est donne´e. Nous illustrons le lien entre cette
de´composition conjointe et la se´paration de sources non circulaires par l’interme´diaire d’une fonction de contraste. Des simulations
illustrent l’approche propose´e.
Abstract – We consider the joint diagonalization problem of complex symmetric third order tensors. A Jacobi like algorithm is
presented where the optimal solution in the 2× 2 case is derived. In case of the source separation problem of non circular signals,
a link between the decomposition and a given contrast function is enlightened. Computer simulations illustrate the proposed
approach.
1 Introduction
Dans le domaine de la se´paration de sources, l’utilisation
de statistiques d’ordre supe´rieur a` deux par l’interme´diaire
de cumulants est devenue relativement classique. La pro-
prie´te´ de multi-line´arite´ des cumulants fait qu’on les quali-
fie souvent de tenseurs. L’ordre du tenseur peut alors eˆtre
directement l’ordre des cumulants conside´re´s. Lorsque que
les signaux sources sont suppose´s statistiquement inde´pen-
dants, le tenseur cumulant associe´ a` ces signaux est diago-
nal. Cette proprie´te´ e´tant perdue par l’action du me´lange,
la restauration de celle ci permettra d’estimer le me´lange.
On parle alors de diagonalisation de tenseur. L’approche
peut-eˆtre la plus connue de ce type a e´te´ propose´e dans
[2] ou`, pour un me´lange instantane´ et apre`s blanchiment
des observations, la matrice de me´lange alors unitaire est
estime´e par un algorithme de type Jacobi.
Dans le meˆme cadre, il a aussi e´te´ montre´ que l’utilisa-
tion de tranches de tenseur cumulant pouvait eˆtre parti-
culie`rement inte´ressante. Ces tranches d’ordre infe´rieur au
tenseur initial peuvent eˆtre utilise´es conjointement pour
estimer la matrice de me´lange. On parle alors de diago-
nalisation conjointe de tenseurs dont l’ordre correspond
a` celui des tranches. Lorsque qu’il s’agit de matrice, cela
a e´te´ propose´e dans [1] pour un tenseur cumulant initial
d’ordre 4 et lorsqu’il s’agit de tenseur d’ordre 3, cela a e´te´
propose´e dans [3] pour un tenseur cumulant initial d’ordre
3 ou 4. La ge´ne´ralisation de ces approches pour des ten-
seurs cumulant initiaux d’ordre quelconque (supe´rieur ou
e´gal a` 3) a e´te´ propose´ dans [4] ou` il a aussi e´te´ montre´
que l’on pouvait tre`s simplement conside´rer des tranches
tensorielles d’ordre fixe de cumulants d’ordre diffe´rent.
Dans le cas de signaux complexes, les cumulants ne sont
pas de´finis de manie`re unique. Pour un ordre donne´, ils de´-
pendent du nombre relatif de conjugaison complexe consi-
de´re´. Le choix du nombre pre´ce´dent implique le type de
de´composition des tenseurs. Pour les matrices, il y a en
fait deux types de de´composition. La de´composition ma-
tricielle classique a e´te´ conside´re´e dans [1], tandis que la
seconde a e´te´ conside´re´e dans [5]. Pour les tenseurs d’ordre
3, il y a encore deux types de de´composition. Une premie`re
de´composition a e´te´ conside´re´e dans [3] tandis que la se-
conde ne semble pas avoir e´te´ conside´re´e jusqu’a` mainte-
nant.
L’objectif de cette communication est double. Le pre-
mier concerne l’e´tude de la deuxie`me de´composition pour
des tenseurs d’ordre 3 dans le cadre des algorithmes de
type Jacobi. Le second objectif concerne l’e´tude de la prise
en compte de l’ensemble des statistiques disponibles pour
la se´paration de sources non-circulaires de type te´le´com-
munications nume´riques, [6] [7] [8].
2 Formulation du proble`me
Dans la suite on ne va conside´rer que des tenseurs d’ordre
3 carre´s de dimension (N,N,N) avec N ∈ N∗ \ {1}. On
conside`re un ensemble T de Nt ∈ N
∗ tenseurs complexes
d’ordre 3 T(`), ` ∈ {1, . . . , Nt} ≡ INt , admettant tous
la de´composition suivante : il existe une matrice (N,N)
unitaire V et Nt tenseurs (N,N,N) d’ordre 3 diagonaux
Λ(`), ` ∈ INt , tels que
Ti,j,k(`) =
N∑
a,b,c=1
V ∗a,iV
∗
b,jV
∗
c,kΛa,b,c(`) , ∀` ∈ INt (1)
ou` (·)∗ est l’ope´rateur de conjugaison complexe. Pour des
raisons de simplicite´, on re´e´crit la relation (1) sous la forme
compacte
T(`) = Λ(`)×3 V
H , ∀` ∈ INt (2)
Colloque GRETSI, 11-14 septembre 2007, Troyes 1009
ou` (·)H est l’ope´rateur de transposition conjugaison com-
plexe des matrices.
E´tant donne´ les tenseurs T(`), ` ∈ INt , (typiquement
a` partir de l’estimation de grandeurs statistiques comme
les cumulants) le proble`me principal consiste en l’estima-
tion de la matrice unitaire V. Pour cela il est classique
d’introduire le crite`re quadratique d’ajustement de mo-
de`le suivant
CT (U, {D}) =
Nt∑
`=1
∥∥T(`)−D(`)×3 UH∥∥2 (3)
ou` ‖ · ‖2 correspond a` la norme Euclidienne du tenseur
conside´re´, i.e. ‖T‖2 =
∑N
l1,l2,l3=1
|Tl1,l2,l3 |
2. Pour des rai-
sons de simplicite´, dans la de´finition pre´ce´dente (3), {D}
correspond a` l’ensemble des tenseurs diagonaux de la de´-
composition. Il s’agit donc maintenant de re´soudre le pro-
ble`me de minimisation suivant
min
U,{D}
CT (U, {D}) . (4)
A cette fin, CT (U, {D}) peut eˆtre d’abord minimise´ par
rapport aux tenseurs diagonaux D(`), ` ∈ INt , tout en
fixant la matrice unitaire U. Pour cela, on prendra en
compte les deux re´sultats suivants dont la preuve, relati-
vement simple, est omise. Etant donne´ un tenseur T, si U
est une matrice unitaire alors on a
‖T×3 U‖
2
= ‖T‖2 (5)(
T×3 U
H
)
×3 U = T . (6)
En utilisant la proprie´te´ (5), CT (U, {D}) peut eˆtre re´e´crit
comme
CT (U, {D}) =
Nt∑
`=1
∥∥T(`)×3 U− (D(`)×3 UH)×3 U∥∥2
(7)
et en utilisant maintenant la proprie´te´ (6), on obtient
CT (U, {D}) =
Nt∑
`=1
‖T(`)×3 U−D(`)‖
2
. (8)
Ainsi la minimisation de CT (U, {D}) par rapport aux ten-
seurs diagonauxD(`), ` ∈ INt , est e´quivalente a` la minimi-
sation de chacun des termes de la somme. Donc, pour une
valeur donne´e de ` et pour une matrice unitaire fixe´eU, la
valeur optimale deD(`) note´e D̂(`) minimisant CT (U, {D})
vaut directement
D̂(`) = Tdiag (T(`)×3 U) (9)
ou` Tdiag(T) est le tenseur diagonal de diagonal e´gale a`
celle de T, i.e. on a
(Tdiag(T))l1,l2,l3 = Tl1,l2,l3δl1,l2,l3 .
Maintenant en remplac¸ant D(`) par D̂(`) dans (8), on a
CT (U, {D̂}) =
Nt∑
`=1
‖T(`)×3 U− Tdiag (T(`)×3 U)‖
2
(10)
et ainsi
CT (U, {D̂})
=
∑Nt
`=1 ‖T(`)×3 U‖
2 −
∑Nt
`=1 ‖Tdiag (T(`)×3 U)‖
2
=
∑Nt
`=1 ‖T(`)‖
2 −
∑Nt
`=1 ‖Tdiag (T(`)×3 U)‖
2
.
(11)
Donc, comme
∑Nt
`=1 ‖T(`)‖
2 est un terme constant, la mi-
nimisation de CT (U, {D̂}) par rapport a`U est e´quivalente
a` la maximisation de
CT (U) =
Nt∑
`=1
‖Tdiag (T(`)×3 U)‖
2
. (12)
Finalement, le proble`me initial consiste maintenant en la
recherche d’une matrice Û satisfaisant
Û = argmax
U
CT (U) . (13)
3 Algorithmes de type Jacobi
Nous conside´rons maintenant l’optimisation du crite`re
pre´ce´dent par un algorithme de type Jacobi. On de´com-
pose ainsi la matrice unitaire recherche´e en un produit de
rotations planes. Le but revient alors a` la de´termination de
chacune de ces rotations. Dans ce cadre, on montre alors
que l’on peut simplement conside´rer le cas N = 2 [2]. Dans
ce cas, la matrice unitaire U peut eˆtre parame´tre´e par
U =
(
cos(α) − sin(α) exp(jφ)
sin(α) exp(−jφ) cos(α)
)
(14)
ou` α et φ sont deux angles. Nous avons alors le re´sultat
suivant :
Proposition 1 Le crite`re dans (12) peut s’e´crire sous la
forme
CT (U) = u
TB1u (15)
avec
u = (cos(2α) sin(2α) sin(φ) sin(2α) cos(φ)) T (16)
et B1 est un matrice dont l’expression est donne´e en fin
de de´monstration.
De´monstration. Nous pouvons e´crire, avec la proprie´te´ de
syme´trie du tenseur T(`) soit Tppk(`)=Tpkp(`) =Tkpp(`),
les 2 e´le´ments de la diagonale d’un tenseur (2×2×2) D(`)
D111(`) = T111(`) cos
3(α)− 3T112(`) sin(α) cos
2(α)ejφ · · ·
+3T122(`) sin
2(α) cos(α)ej2φ −T222(`) sin
3(α)ej3φ
D222(`) = T111(`) sin
3(α)e−j3φ + 3T112(`) cos(α) sin
2(α) · · ·
e−2jφ + 3T122(`) cos
2(α) sin(α)e−jφ +T222(`) cos
3(α).
(17)
1010
Nous obtenons alors la norme des e´le´ments diagonaux de
ce tenseur
|D111(`)|
2 + |(D222(`)|
2
=
(
|T111(`)|
2 + |T222(`)|
2
) (
cos6(α) + sin6(α)
)
+ 2.25
(
|T112(`)|
2 + |T)222(`)|
2
)
sin2(α)
+ 3 Re
{
T222(`)(T122(`)
∗ejφ −T111(`)T112(`)
∗e−jφ
}
sin(2α) cos(2α)
+ 1.5 Re
{
T111(`)T122(`)
∗e−j2φ +T222(`) T112(`)
∗ejφ
}
sin2(2α)
(18)
qui peut s’e´crire
|D111(`)|
2 + |D222(`)|
2 = uT B(`) u (19)
ou` u est un vecteur re´el (3 × 1) tel que uTu = 1 et de´fini
par (16), et B(`) est une matrice re´elle syme´trique (3× 3)
donne´e par
B11(`) = |T111(`)|
2 + |T222(`)|
2
B12(`) = −1.5 Im {T222(`)T122(`)
∗ +T111(`)T112(`)
∗}
B13(`) = 1.5 Re {T222(`)T122(`)
∗ −T111(`)T112(`)
∗}
B22(`) = 0.25
(
|T111(`)|
2 + |T222(`)|
2
)
+2.25
(
|T112(`)|
2 + |T122(`)|
2
)
−1.5 Re {T111(`)T122(`)
∗ +T222(`)T112(`)
∗}
B23(`) = 1.5 Im {T111(`)T122(`)
∗ −T222(`)T112(`)
∗}
B33(`) = 0.25
(
|T111(`)|
2 + |T222(`)|
2
)
+2.25
(
|T112(`)|
2 + |T122(`)|
2
)
+1.5 Re {T111(`)T122(`)
∗ +T222(`)T112(`)
∗}
(20)
La maximization de la somme des normes des e´le´ments
diagonaux de l’ensemble des tenseurs {D(`)} est obtenue
par
Nt∑
l=1
|(D(`))111|
2 + |(D(`))222|
2 = uT B1 u (21)
ou` B1 est une matrice re´elle syme´trique de´finie par
B1 =
Nt∑
l=1
B(`) (22)
ce qui comple`te la de´monstration. ♦
La maximisation du crite`re (15) peut facilement eˆtre
obtenu via une de´composition en e´le´ments propres de la
matrice B1. A partir du vecteur propre associe´ a` la plus
grande valeur propre, les angles sont alors obtenus en uti-
lisant
cos(α) =
√
1 + 12u(1)
sin(α) exp(jφ) = 12 cos(α) (u(3) + j u(2))
(23)
avec α ∈ [−pi/4, pi/4].
Finalement, la matrice unitaire inconnue U est obtenue
par multiplication successive des rotations planes la de´-
composant.
En ge´ne´ral, il peut eˆtre a priori inte´ressant d’utiliser
l’ensemble des informations disponibles pour re´soudre un
proble`me. Dans le cadre de l’algorithme STOTD [3], la de´-
composition suivante de tenseurs d’ordre 3 est conside´re´e
T ′i,j,k(`) =
N∑
a,b,c=1
Va,iV
∗
b,jV
∗
c,kΛ
′
a,b,c(`) , ∀` ∈ INt′ (24)
ou` les Λ′(`) sont des tenseurs diagonaux complexes d’ordre
3. Il a e´te´ montre´ dans [3], que dans le cas de N = 2, le
crite`re (12) adapte´ a` cette de´composition peut s’e´crire sous
la forme
CT (U) = u
TB2u (25)
ou` B2 est une matrice re´elle syme´trique.
Ainsi une fac¸on de combiner ces deux types de de´com-
position de tenseurs revient simplement a` conside´rer la
matrice B1 +B2 dans le proble`me d’optimisation.
L’algorithme de type Jacobi prenant en compte la seule
matriceB1 est appele´ NC-STOTD tandis que l’algorithme
prenant en compte B1 +B2 est appele´ G-STOTD.
4 Lien avec la Se´paration de Sources
La se´paration de sources repose sur le mode`le de signaux
suivant
x[n] = As[n] (26)
ou` n ∈ Z est le temps discret, s[n] le (N, 1) vecteur des
N 6= 2 signaux sources complexes inobservables si[n], i ∈
{1, . . . , N}, x[n] le (N, 1) vecteur des observations xi[n],
i ∈ {1, . . . , N} et A la matrice de me´lange carre´e (N,N)
suppose´e inversible.
De plus les signaux sources sont suppose´s eˆtre des si-
gnaux ale´atoires de moyenne nulle, de puissance unite´, de
cumulant d’ordre 4 non nul et statistiquement mutuelle-
ment inde´pendant. Par ailleurs on supposera qu’un blan-
chiment des observations a e´te´ re´alise´ de telle sorte que la
matrice A est maintenant conside´re´e unitaire.
La se´paration aveugle de sources consiste en l’estimation
d’une matrice unitaire H de telle sorte que le vecteur
y[n] = Hx[n] (27)
restore un signal source diffe´rent sur chacune de ses com-
posantes.
Peut eˆtre l’une des fac¸ons la plus pratique de re´soudre
le proble`me de la se´paration de sources passe par l’utilisa-
tion de fonctions de contraste [2]. Elles correspondent a` des
fonctions de couˆt de´pendant de statistiques des sorties du
syste`me de se´paration dont la maximisation globale suf-
fit pour obtenir une solution. Nous proposons le re´sultat
suivant :
Proposition 2 La function
J (y) =
N∑
i,j=1
|Cum{yi, yi, yi, yj}|
2
(28)
est une fonction de contraste.
Par manque de place, la preuve est reporte´e dans une com-
munication future.
Plus important est le lien existant entre ce contraste et
un crite`re de diagonalisation conjointe. Il est donne´e par
le re´sultat suivant :
Colloque GRETSI, 11-14 septembre 2007, Troyes 1011
Proposition 3 Conside´rons l’ensemble T c des N tenseurs
d’ordre 3
Tc(`) =
(
T ci,j,k(`)
)
, ∀` ∈ IN
de´finis par
T ci,j,k(`) = Cum{xi, xj , xk, x`} . (29)
Si H est une matrice unitaire, on a alors
J (Hx) = CT c(H) . (30)
5 Simulations
Dans le cadre de la se´paration de sources nous illus-
trons les performances obtenues a` partir des algorithmes
JADE [1], STOTD [3] et les deux approches pre´sente´es
NC-STOTD et G-STOTD. Les signaux sources utilise´s
prennent leurs valeurs dans {−1;−j; 0;β; jβ} avec les pro-
babilite´s respectives
{
1
2(1+β) ;
1
2(1+β) ;
β−1
β
; 12β(1+β) ;
1
2β(1+β)
}
.
On conside`re β = 1, 4. Le nombre d’e´chantillons sera 10000.
Nous allons conside´rer successivement deux cas de me´-
lange. Dans le premier cas, il s’agit d’un me´lange 2 × 2
unitaire fixe´. Le blanchiment n’est pas re´alise´. Dans le
deuxie`me cas, il s’agit d’un me´lange 4 × 2 dont les com-
posantes sont tire´es de manie`re ale´atoire suivant une loi
Gaussiennes, centre´es et de puissance 1. Le blanchiment
est re´alise´ avec une projection sur l’espace signal de dimen-
sion 2. Les me´langes sont bruite´s avec un rapport signal
sur bruit (RSB) allant de -5 dB a` 30 dB.
L’indice de performance utilise´ est celui donne´ dans [4].
La moyenne sur 1000 tirages (pour les sources, le bruit et,
dans le deuxie`me cas, le me´lange) est conside´re´e.
Pour les deux cas de me´lange, les re´sultats sont donne´s
par les figures 1 et 2 respectivement. Nous constatons que
l’on obtient les meilleures performances avec l’algorithme
G-STOTD, quelque soit le RSB. De plus les performances
obtenues par l’algorithme NC-STOTD sont meilleures que
celle obtenues par STOTD et JADE (qui ont les meˆmes
performances, ce qui va dans le sens de [3]).
−5 0 5 10 15 20 25 30
10−6
10−5
10−4
10−3
10−2
10−1
RSB
Pe
rf
or
m
an
ce
JADE
STOTD
NC−STOTD
G−STOTD
Fig. 1 – Performance des diffe´rents algorithmes en fonc-
tion du RSB dans le cas d’un me´lange unitaire.
−5 0 5 10 15 20
10−4
10−3
10−2
RSB
Pe
rf
or
m
an
ce
JADE
STOTD
NC−STOTD
G−STOTD
Fig. 2 – Performance des diffe´rents algorithmes en fonc-
tion du RSB dans le deuxie`me cas de me´lange.
Ainsi l’association de plusieurs grandeurs statistiques
associe´es ici a` des signaux sources non circulaires permet
d’obtenir, dans le cadre pre´sente´, de meilleures perfor-
mances que l’approche “classique”.
Re´fe´rences
[1] J.-F. Cardoso and A. Souloumiac, “Blind Beamfor-
ming for Non Gaussian Signals”, IEE Proceedings-F,
Vol. 140, pp. 362-370, 1993.
[2] P. Comon,“Independent Component Analysis, a New
Concept ?”, Signal Processing, Vol. 36, pp. 287-314,
1994.
[3] L. De Lathauwer, B. De Moor and J. Vanderwalle,
“Independent Component Analysis and (Simulta-
neous) Third-Order Tensor Diagonalization”, IEEE
Transactions on Signal Processing, Vol. 49, pp. 2262-
2271, 2001.
[4] E. Moreau, “A Generalization of Joint-
Diagonalization Criteria for Source Separation”,
IEEE Transactions on Signal Processing, Vol. 49, pp
530-541, 2001.
[5] L. De Lathauwer, B. De Moor and J. Vanderwalle,
“ICA Techniques for More Sources than Sensors”,
Proceeding of the IEEE Signal Processing Workshop
on Higher-Order Statistics (HOS’99), pp. 121-124,
June 1999, Caesarea, Israel.
[6] L. De Lathauwer and B. De Moor, “On the Blind
Separation of Non-Circular Sources”, Proceeding of
EUSIPCO-02, Vol. II, pp. 99-102, Sept. 2002, Tou-
louse, France.
[7] P. Chevalier, “Optimal Array Processing for Non Sta-
tionary Signals”, in Proc. ICASSP’96, pp. 2868-2871,
Atlanta, USA, May 1996.
[8] P. Chevalier, “New Insights into Optimal Widely Li-
near Array Receivers for the Demodulation of BPSK,
MSK, and GMSK Signals Corrupted by Noncircular
Interferences - Application to SAIC”, IEEE Transac-
tions on Signal Processing, Vol. 54, pp. 870-883, 2006.
1012


Diagonalisation conjointe de tenseurs d'ordre 3 Application à la séparation de signaux non-circulaires

http://documents.irevues.inist.fr/bitstream/2042/17675/1/GRETSI_2007_1009.pdf

Diagonalisation conjointe de tenseurs d'ordre 3 Application à la séparation de signaux non-circulaires

Abstract

Similar works

Full text

Available Versions

I-Revues