Bertrand Rivet

Christian Jutten

WHEN SILENCE IS GOLD

Audiovisual source separation is a fascinating approach to source extraction. Several algorithms have already been proposed for extracting speech sources from audio mixtures by exploiting audiovisual coherence. One of the main property of speech signals is that they are highly non-stationary: there are periods during which speakers do not produce sounds. In this work, the audiovisual coherence is used to estimate such silent periods which are then useful to extract corresponding speech signals.

INTRODUCTION

Have you ever wondered why it could be so difficult to hear the driver of a car from a back seat, particularly in a noisy environment? It is well known that the human brain is able to discriminate sounds: it can enhance the different audio sources so that the interest source(s) can be more easily extracted. However, if the environment noise becomes too loud, this faculty could be deteriorated. As a consequence, how can it be explained that in the same noisy conditions, you can understand what your neighbor says? The answer comes from the cognitive sciences: your brain is able to merge the heard signals and what you see, especially the motion of speech articulators as lips or tongue. Indeed, speech is multimodal [START_REF] Erber | Interaction of audition and vision in the recognition of oral speech stimuli[END_REF][START_REF] Sumby | Visual contribution to speech intelligibility in noise[END_REF][START_REF] Summerfield | Some preliminaries to a comprehensive account of audio-visual speech perception[END_REF]. The bimodal nature of speech is now a basic feature, both for the understanding of speech perception [START_REF] Summerfield | Some preliminaries to a comprehensive account of audio-visual speech perception[END_REF] and for the development of tools for human-machine communication [START_REF] Stork | Speechreadings by Humans and Machines[END_REF][START_REF] Potamianos | Joint Audio-Visual Speech Processing for Recognition and Enhancement[END_REF]. It is thus well known that the vision of a speaker's face influences what is perceived. One of the most famous example is the McGurk's effect [START_REF] Mcgurk | Hearing lips and seeing voices[END_REF] which is an audiovisual illusion: superimposing an audio stimulus [ba] with a video stimulus [ga] leads to perceive [da]. Indeed, most of people are able to lip read (e.g., [START_REF] Benoît | Effects of phonetic context on audio-visual intelligibility of French[END_REF][START_REF] Erber | Interaction of audition and vision in the recognition of oral speech stimuli[END_REF][START_REF] Sumby | Visual contribution to speech intelligibility in noise[END_REF]) even if they are not conscious of doing so: the recognition rate with audiovisual speech increases compared to audio speech only [START_REF] Erber | Interaction of audition and vision in the recognition of oral speech stimuli[END_REF]. In addition, the vision of the speaker's face not only favors understanding but also enhances speech detection in a noisy environment [START_REF] Bernstein | Auditory speech detection in noise enhanced by lipreading[END_REF][START_REF] Grant | The use of visible speech cues for improving auditory detection of spoken sentences[END_REF][START_REF] Kim | Investigating the audio-visual speech detection advantage[END_REF]. All of these studies finally exploit the redundancies and complementarities between visual and audio modalities of speech. It seems actually intuitive that there are strong links between the motion of the speaker's lips and the speech signal.

On the other hand, blind source separation (BSS) is an exciting field of signal processing. It aims at recovering unknown sources from mixtures of them [START_REF] Comon | Handbook of Blind Source Separation Independent Component Analysis and Applications[END_REF]. The lack of prior knowledge about the sources and the mixing process is generally overcame by strong assumptions. For instance, the mutual independence be-tween the sources like in independent component analysis (ICA) e.g., [START_REF] Comon | Independent component analysis, a new concept?[END_REF][START_REF] Cardoso | Blind signal separation: statistical principles[END_REF] or the sparsity of the signals (sparse component analysis, e.g., [START_REF] Gribonval | A survey of Sparse Component Analysis for Blind Source Separation: principles, perspectives, and new challenges[END_REF]). Another additional prior information, close to the sparsity, is the inactivity of some sources in some time windows [START_REF] Rivet | Blind extraction of intermittent sources[END_REF][START_REF] Araki | Blind speech separation in a meeting situation with maximum snr beamformers[END_REF][START_REF] Deville | Temporal and time-frequency correlation-based blind source separation methods. part I: Determined and underdetermined linear instantaneous mixtures[END_REF]. These inactivity periods can be silences in the case of speech signals.

In the last decade, these two fields of signal processing (BSS and audiovisual speech processing) have merged into the audiovisual source separation problem which aims at extracting speech signal(s) from audio mixtures when visual information of the speaker(s) are available. Several approaches (e.g., [START_REF] Sodoyer | Separation of audio-visual speech sources: a new approach exploiting the audiovisual coherence of speech stimuli[END_REF][START_REF] Sodoyer | Developing an audio-visual speech source separation algorithm[END_REF][START_REF] Dansereau | Co-channel audiovisual speech separation using spectral matching constraints[END_REF][START_REF] Rajaram | Bayesian separation of audio-visual speech sources[END_REF][START_REF] Wang | Video assisted speech source separation[END_REF][START_REF] Sigg | Nonnegative CCA for audiovisual source separation[END_REF][START_REF] Rivet | Mixing audiovisual speech processing and blind source separation for the extraction of speech signals from convolutive mixtures[END_REF][START_REF] Rivet | Visual voice activity detection as a help for speech source separation from convolutive mixtures[END_REF][START_REF] Liu | Use of bimodal coherence to resolve spectral indeterminacy in convolutive bss[END_REF][START_REF] Llagostera Casanovas | Blind audiovisual source separation based on sparse redundant representations[END_REF]) have been developped to deal with this problem from the very basic filtering process using enhancement filters estimated with the help of lip shape information [START_REF] Girin | Audio-visual enhancement of speech in noise[END_REF] to the definitively more complex case of moving sources [START_REF] Naqvi | A multimodal approach to blind source separation of moving sources[END_REF]. By exploiting visual information to this problem, one want to provide an easier extraction of speaker(s) than purely audio methods. This paper deals with the extraction of speaker(s) in a multimodal way by exploiting the silence periods in the speech signal. Compared to previous work [START_REF] Rivet | Visual voice activity detection as a help for speech source separation from convolutive mixtures[END_REF], the proposed study provides a better visual detection of silent periods and a more robust estimation of the extraction vector related to the target speaker.

This paper is organized as follows. Section 2 briefly describes an audiovisual voice activity detection (V-VAD) exploited by the proposed method (Section 3) to extract speech signals. Numerical experiments and results are given in Section 4 before conclusions and perspective in Section 5.

AUDIOVISUAL COHERENCE: VISUAL VOICE ACTIVITY DETECTION

Speech is at least a bimodal signal, both acoustic and visual, leading thus to exploit these two modalities in a complementary fashion to improve methods generally derived from audio information only. The possibility of using visual information to detect speech and silence periods in a given audio channel is considered in this section.

Visual voice activity detection (V-VAD) has an advantage over audio only based voice activity detection in that it is not susceptible to the problems associated with the acoustic environment (e.g., noise, simultaneous speakers, reverberations, etc.). The V-VAD used in this study is based on a retinal filter [START_REF] Hérault | Motion processing in the retina: About a velocity matched filter[END_REF][START_REF] Aubrey | Two novel visual voice activity detectors based on appearance models and retinal filltering[END_REF] in order to detect the movement of the speaker's lips. This method requires no prior information: it is based on the assumption that silence frames can be mostly characterized by the lip-shape movements. Indeed, in silence sections, the lip-shape variations are often small. On the contrary, during speech sections these variations are generally quite stronger. Given these observations, the silence sections can be identified with one or several visual parameters of dynamical nature. The used V-VAD applies a retinal filter to the images flow and calculates the change in energy to classify voice activity. The used spatio-temporal filter (Fig. 1) mimics the output of the parvo cells in a human's retina which enhances the moving contours in the image flow [START_REF] Hérault | Motion processing in the retina: About a velocity matched filter[END_REF] as illustrated in Fig. 2 in which white pixels correspond to moving objects in the image. Let denote by R i,j (t) the output of the retinal filter at time t in the (i, j)-th pixel. The change in energy is then simply computed by

η(t) = 1 IJ I i=1 J j=1 R i,j (t), (1)
where I and J are the number of rows and columns pixels respectively. The probability that at time t the speaker is silent given η(t) is finally provided by

P V (t) = .98 exp - η(t) 45 . (2)
The parameters of this expression are estimated from the analysis of the database and they are speaker dependent (see Section 4.1).

AUDIOVISUAL SOURCE EXTRACTION

In this section, the proposed method to extract audio signal related to one (or several) speaker(s) who is (are) filmed is presented (Subsection 3.1) before a brief recall of the used framework (Subsection 3.2).

Extraction of intermittent sources

Let consider the linear instantaneous mixing model

x(t) = As(t), (3)
with as many sources s

(t) = [s 1 (t), • • • , s N (t)] T as mix- tures x(t) = [x 1 (t), • • • , x N (t)] T .
Let us denote the covariance matrix of mixtures x(t) at sample t estimated from a T -length window by

R T,x (t) = E T [x(t)x(t) T ] = N i=1 σ 2 T,i (t)a i a T i ( 4
)
where σ 2 T,i (t) is the average power of the i-th windowed source at sample t and a i is the i-th column of the mixing matrix

A = [a 1 , • • • , a N ]
. Let now suppose that • N 1 sources1 are simultaneously inactive at time τ while they are active at time t, i.e. 1 ≤ i ≤ N 1 , σ 2 T,i (τ ) = 0 and σ 2 T,i (t) = 0, • the N -N 1 other sources have the same average power at times τ and t, i.e.

N 1 + 1 ≤ i ≤ N, σ 2 T,i (τ ) = σ 2 T,i (t) = 0.
It is easy to show that the generalized eigenvalue decomposition [START_REF] Golub | Matrix Computation[END_REF] of pair (R T,x (τ ), R T,x (t)) admits only two distinct generalized eigenvalues: 0 degenerated N 1 times whose eigensubspace E 0 (τ ) is orthogonal to the space spanned by {a N1+1 , • • • , a N } and 1 degenerated N -N 1 times. As a consequence, the projection of observations x(t) onto E 0 (τ ) cancels the contribution of the sources s i (t), N 1 + 1 ≤ i ≤ N . This means that the separation vectors b i , 1 ≤ i ≤ N 1 , lie in E 0 (τ ), where b i is the i-th column of the separation matrix B.

This method allows us to detect how many sources are vanishing by testing the generalized eigenvalues and then to extract the space spanned by the corresponding sources by projecting the observations onto the generalized eigenvectors associated with the eigenvalues equal to zero. Nevertheless, it is worth noting that this method does not allow to identify which sources are inactive. However, even if this difficulty can be overcome by a clustering stage [START_REF] Rivet | Blind extraction of intermittent sources[END_REF], the bimodality of speech is used in this study to extract particular speakers.

Application to speech source extraction

In this section, the proposed method mixing audio and video processes to extract specific speaker(s) from mixtures is presented. It is a two stages method which first estimates the periods of the speaker's silence and then estimates the corresponding separation vector.

To estimate the periods of silence related to the speaker of interest, the audio and video informations are merged into the probability that this specific speaker is silent at time τ by

P AV (τ ) = P A (τ )P V (τ ), (5)
where P V (τ ) is the visual probability that the speaker of interest is silent (2) and P A (τ ) is the audio probability that at least one source is inactive. This latter probability is obtained by mapping the smallest eigenvalue λ 1 (τ ) of pair (R T,x (τ ), R T,x (t)) into [0, 1] thanks to the following sigmoid function

P A (τ ) = 1 - 1 1 + exp -α(log(λ 1 (τ )) -β) (6)
where α and β are parameters that allow to customize the detection of inactive sources (Fig. 3): if log(λ 1 (τ ))β is large (resp. small) compared to 1/α then P A (τ ) is about 0 (resp. 1). The periods of silence are then estimated as the time indices τ such that P AV (τ ) > .5.

In the second step, the separating vector b 1 of the target speaker 2 is estimated from the inactivity periods estimated in the first step of the method. While only one generalized eigenvector is theoretically needed to estimate the separating vector, to obtain a better estimation of this latter vector several eigenvectors are used. If only the target speaker is silent (s 1 (t) = 0), then the corresponding eigenvector is align with the separating vector b 1 . However, if several sources are simultaneously inactive with the target speaker, then the corresponding eigenvectors lying in E 0 (τ ) are not necessary align with the related separating vectors b i . It is worth noting that if dim(E 0 (τ ) ∩ E O (τ )) = 1, the support vector u τ,τ lying in this intersection is then necessary aligned with the separating vector b 1 related to the target speaker. As a consequence, one can estimate the separating vector b 1 as a weighted average of the eigenvectors when only one source is silent and of the support vectors u τ,τ . But as already mentioned, to have a better estimation of b 1 a weighted kernel-PCA is proposed [START_REF] Müller | An introduction to kernel-based learning algorithms[END_REF][START_REF] Rivet | Blind extraction of intermittent sources[END_REF] with kernel

ψ(u τ1,τ2 , u τ3,τ4 ) = |u T τ1,τ2 u τ3,τ4 | -cos(θ 0 ) 1 -cos(θ 0 ) , if |u T τ1,τ2 u τ3,τ4 | ≥ cos(θ 0
) and ψ(u τ1,τ2 , u τ3,τ4 ) = 0 else, with the convention that u τ,τ is a generalized eigenvector when only one source is inactive.

θ 0 is an a priori chosen angle which defines the maximum angle of the beam formed by several estimations of b 1 . It corresponds to the angular resolution. The probability of the inactivity of the target speaker is included in this model by weighting the kernel ψ w (u τ1,τ2 , u τ3,τ4 ) = w τ1,τ2 w τ3,τ4 ψ(u τ1,τ2 , u τ3,τ4 ) with w τi,τj = P AV (τ i )P AV (τ j ). Finally, the separating vector b 1 is estimated by

b 1 = U Ψ w φ 1 , (7)
2 Without loose of generality the target speaker is arbitrary chosen to be s 1 (t). where Ψ w is the symmetric kernel matrix whose entries are ψ w (u τ1,τ2 , u τ3,τ4 ) and φ 1 is the eigenvector of Ψ w associated with the largest eigenvalue. The target speaker is then estimated by

ŝ1 (t) = b T 1 x (8)
for all samples t including those when the target speaker is speaking.

RESULTS

In this study, two distinct databases are used to the target speaker (database D1) and to the concurrent audio sources (database D2). The database D1 [START_REF] Sodoyer | A study of lip movements during spontaneous dialog and its application to voice activity detection[END_REF] has been recorded by two French speakers in a spontaneous dialog condition. Both of them have been simultaneously recorded by a microphone (sampling rate of 16kHz) and a camera (sampling rate of 25 frames per second) focussed on the lips region. The second database D2 is composed of speech signals and of music songs. All these signals have a sampling rate of 16kHz.

Visual voice activity detection

In this section, the analysis of the change in energy (1) is presented (Fig. 4). As one can see, the values of η for the speech and silence classes are not perfectly separated. Besides, even during silence, the speaker's mouth may sometimes slightly move: however the less the lips are moving, the higher the probability of silence is (bottom plot). Indeed, a direct VAD from raw lip shape cannot lead to satisfactory performances because of the intricate relationship between visual and acoustic speech information. Furthermore, the proposed model of the visual probability of silence [START_REF] Aubrey | Two novel visual voice activity detectors based on appearance models and retinal filltering[END_REF] given the change in energy (1) provides a good fitting of the experimental values. As shown on Fig. 5 which compares the timemoving power of audio signal (when there is no additive noise) and the change in energy η(t), these two values are highly correlated as well as the visual probability of silence (2) with the ground truth indexation of silence. This model finally provides a smooth visual decision between silence and non-silence frames. It is worth noting that in real conditions (i.e. with several competitive sources) η(t) remains unchanged while the profile of the audio power is altered by the competitive sources. In the bottom plot, the black line is the ground truth of silence (1 means silence while 0 means speech) and the grey line is the visual probability of silence p V (t) (2), with α = 2 and β = -1.

Audiovisual source extraction

In this experiment, one target speaker from D1 is linearly mixed with two other sources from D2 (Fig. 6).

In this example, the two concurrent sources are speech (s 2 ) and music (s 3 ), respectively. The audio signals have been windowed synchronously with the video signal resulting a 40ms frames without overlapping. To compute the two covariance matrices R T1,x (t) and R T,x (τ ), T 1 is equal to 5 seconds to have a reference of the average powers of the non-stationary sources and T is equal to 40ms to follow the variations of power of the sources. Fig. 6(c) presents the probabilities of silence obtained from audio only P A (t), video only P V (t) and audiovisual P AV (t), respectively. As one can see, the highest values of the audio probability of silence (6) correspond to the silent frames of all sources. In this example the two speech sources (s 1 and s 2 ) present inactive periods. However, as already mentioned the visual probability of silence is highly correlated with the actual periods of silence of the target speaker (s 1 ). As a consequence, the proposed audiovisual probability of silence P AV (t) has two main advantages. First, it allows to select among all the periods of silence of at least one source those which mainly correspond to the target speaker. Secondly, it allows to weight the related eigenvectors which are used by the kernel-PCA to estimate the extraction vector. Finally, the target speaker is correctly extracted from the mixtures by the proposed method.

CONCLUSIONS AND PERSPECTIVES

In this study, nothing more besides the periods of silence of a specific speaker are exploited in a peculiar way to extract this source when it is active. The proposed visual voice activity detection based on a smooth criterion as well as the weighted kernel-PCA provide a soft algorithm to take into account the confidence level of the related silent periods. Moreover, an advantage of this algorithm compared to other methods is that it allows to extract and to identify only the target sources thanks to the visual information. It is worth noting that even if the proposed algorithm is applied on a linear instantaneous mixtures, it can easily be extended to the convolutive case by applying it in the frequency domain. (a) Sources While almost all the classical audiovisual methods to extract speech sources are generally focussed on the periods of activity, conversely the main original assumption of this method is to highlight the inactivity periods of the sources. The latter periods, which are too often ignored, are in fact of a great help to estimate the extraction vector as shown in this study. Future works will focus to propose audiovisual methods to extract moving sources, i.e. with a non-stationary mixing process.

Figure 1 :

 1 Figure 1: Mono-dimensional electrical scheme to model human retina from [15].

Figure 2 :

 2 Figure 2: Illustration of the output (right) of the retinal filter [15] applied on the image flow (left).

α = 2 , 5 Figure 3 :

 253 Figure 3: Probability P A (τ ) of at least one inactive source (6) for several values of parameters α and β.

Figure 4 :

 4 Figure 4: Analysis of change in energy (1) (top plot) and visual probability of silence (2) (bottom plot).

Figure 5 :

 5 Figure 5: Illustration of visual voice activity detector.In the bottom plot, the black line is the ground truth of silence (1 means silence while 0 means speech) and the grey line is the visual probability of silence p V (t) (2), with α = 2 and β = -1.

Figure 6 :

 6 Figure 6: Illustration of the audiovisual speaker extraction.

Without loss of generality, the inactive sources are arbitrary the N 1 first sources.


Gaogang Xie

Peng He

Hongtao Guan

Zhenyu Li

Yingke Xie

Layong Luo

Jianhua Zhang

Yonggong Wang

Kave Salamatian

PEARL: A Programmable Virtual Router Platform

Programmable routers supporting virtualization are a key building block for bridging the gap between new Internet protocols and their deployment in real operational networks. This article presents the design and implementation of PEARL, a ProgrammablE virtuAl Router pLatform with relatively high-performance. It offers high flexibility by allowing users to control the configuration of both hardware and software data paths. The platform makes use of fast lookup in hardware and software exceptions in commodity multi-core CPUs to achieve high speed packet processing. Multiple isolated packet streams and virtualization techniques ensure isolation among virtual router instances.

INTRODUCTION

Deploying, experimenting and testing new protocols and systems over Internet have always been a major issue.

While, one could use simulation tools for the evaluation of a new systems aimed toward large-scale deployment, real experimentation in experimental environment with realistic enough settings, such as real traffic workload and application mix are mandatory. Therefore easy programmable platforms that can support high-performance packet forwarding and enable parallels deployment and experiment of different research ideas are highly demanded. Moreover, we are moving toward a future Internet architecture that seems to be polymorphic rather than monolithic, i.e., the architecture will have to accommodate simultaneous coexistence of several architectures (like the Named Data Network (NDN) [START_REF] Jacobson | Networking Named Content[END_REF], etc.) including the current Internet. Therefore future Internet should be based on platforms running different architectures in virtual slices enabling independent programming and configuration of functions of each individual slice.

Current commercial routers, the most important building blocks for Internet, while attaining very high performance, offer only a very limited access to the researchers and developers to their internal component to implement and deploy innovative networking architecture. In contrast open software based routers naturally facilitate the access and adaptation of almost all their component however with a frequently low packet processing performance.

As an example, recent OpenFlow switches provides flexibility by allowing programmers to configure the 10-tuple of flow table entries, enabling the change of packet processing of a flow. OpenFlow switches are not ready for non-IP based packet flows, such as NDN. Moreover, while the switches allow a number of slices for different routing protocols through the FlowVisor, the slices are not isolated in terms of processing capacity, memory and bandwidth.

Motivated by these facts, we have designed and built a ProgrammablE virtuAl Router pLatform, PEARL that can guarantee high performance. The challenges are two-fold: first to manage the flexibility vs. performance trade-off that translates into pushing functionality to hardware for performance vs. programming them in software for flexibility, second to ensure isolation between virtual router instances both in hardware and software with low performance overhead.

This article describes the PEARL router's architecture, its key components and the performance results. It shows that PEARL meets the design goals of flexibility, high performance and isolation. In the next section, we describe the design goals of the PEARL platform. These goals can be taken as the main features of our designed platform. We then detail the design and implementation of the platform, including both hardware and software platform. In the next section, we evaluate the performance of a PEARL based router using the Spirent TestCenter by injecting into it both IPv4 and IPv6 traffic. Finally, we briefly summarize the related works, and conclude the paper.

DESIGN GOALS

In the past few years, several future Internet architectures and protocols at different layers have been proposed to cope with the challenges that the current Internet faces [START_REF] Jacobson | Networking Named Content[END_REF]. Evaluating the reliability and the performance of these proposed innovative mechanisms is mandatory before envisioning a real scale deployment. Besides theoretically evaluating the performance of a system and simulating these implementations, one needs to deploy them in a production network with real user behavior, traffic workload, resource distribution, and applications mixture.

However, a major principle in experimental deployment is the "No harm" principle that states that normal services on a production network should not be impacted by the deployment of a new service. Moreover, no unique architecture or protocol stack will be able to support all actual and future Internet services and we might need specific packet processing for given services. Obviously, a flexible router platform with high-speed packet processing ability and support of multiple parallel and virtualized independent architectures is extremely attractive for both Internet research and operation. Based on this observation one can define isolation, flexibility, and high performance as the needed characteristics and the design goals of a router platform future Internet.

In particular, the platform should be able to cope with various types of packets including IPv4, IPv6, even non-IP and be able to apply packet routing as well as circuit switching. Various software solutions like Quagga or XORP [START_REF] Handley | XORP: an open platform for network research[END_REF] have provided such flexible platform that is able to adapt their packet-processing components as well as to customize the functionalities of their data, control and management planes. However, these approaches fail to be fast enough to be used in operational context where a wire-speed is needed. Nevertheless, by adding and configuring convenient hardware packet processing resources such as FPGA, CPU cores and memory storage one can hope to meet the performance requirements. Indeed, flexibility and high performance are in conflict in most situations. Flexibility requires more functionalities to be implemented in software to maximize the programmability. On other hand, high performance cannot be reached in software and needs custom hardware. A major challenge for PEARL is to allocate enough hardware and multi-cores in order to achieve both flexibility and high performance.

Another design goal is relative to isolation. By isolation we mean a mechanism that enables different architectures or protocols running in parallel on separate virtual router instances without impacting each other performances. In order to achieve isolation, we should provide a mechanism that will ensure that one instance can only use its allocated hardware (CPU cores and cycles, memory, resources, etc.) and software resources (lookup routing tables, packet queue, etc.) and is forbidden to access resources of other instances even when they are idle. We need also a dispatching component that will ensure that IP or non-IP packet are delivered to specified instances following custom rules defined over MAC layer parameters, protocols, flow label or packet header fields.

The PEARL offers high flexibility through the custom configurations of both hardware data path and software data path. Multiple isolated packet streams and virtualization techniques enable the isolation among virtual router instances, while the fast lookup hardware provides the capacity to achieve high performance

PLATFORM DESIGN AND IMPLEMENTATION

A. System Overview

PEARL uses commodity multi-core CPUs hardware platforms that run generic softwares as well as specialized packet-processing cards for high performance packet processing as shown in Figure 1. The virtualization environment is build using the Linux-based LXC solution. This enables multiple virtual router instances to run in parallel over a CPU core or one router instance over multiple CPU cores. Each virtual machine can be logically viewed as a separate host. The hardware platform contains a FPGA-based packet processing card with embedded TCAM and SRAM. This card enables fast packet processing and strong isolation. LXC takes advantage of a group of the kernel feature (namespace, Cgroup) to ensure isolation in software between virtual router instances. A multi-stream high-performance DMA engine is also used in PEARL, which receives and transmits packets via high-speed PCI Express bus between hardware and software platform. Each IO stream can be either assigned to a dedicated virtual router or shared by several virtual routers using a net IO proxy.

Flexibility. We use TAP/TUN device as the network interface in each virtual machine. Each virtual machine could be considered as a standard Linux host containing multiple network ports. Thus, the IPv4, IPv6, OpenFlow, even Non-IP protocol stack can be easily loaded. Adding new functions to router is also convenient though programming Linux applications. For example, to load IPv4 or IPv6, Quagga routing software suite can be used as the control plane inside each Linux container.

High performance. The operations of routing table lookup and packets dispatch to different virtual machines are always the performance bottleneck. PEARL offloads these two operations into hardware to achieve high speed packet forwarding. In addition, since LXC is a light-weight virtualization technique with low-overhead, the performance is further improved.

B. Hardware Platform

To provide both high performance and strong isolation in PEARL, we design a specialized packet processing card. Figure 2 shows the architecture of hardware data plane. It's a pipeline-based architecture, which consists of two main data paths: the transmitting and the receiving data path. The receiving data path is responsible for processing the ingress packets and the transmitting data path processes the egress packets. In what follows, the processing stages of the pipeline are detailed.

Header Extractor. For each incoming packet, one or many fields are extracted from the packet header. These fields are used for virtual router ID (VID) lookup in the next processing stage. For IP-based protocols, a 10-tuple, defined

following OpenFlow [START_REF] Mckeown | OpenFlow: Enabling innovation in campus networks[END_REF], is extracted, while for non-IP protocols, the MAC address is extracted. Given limited hardware resources, we implement four routing tables in current design. The tables are stored in TCAM as well. We take the VID combined with destination IP address as search key. The VID part of the key is performed exact matching and the IP part is performed the longest prefix matching in TCAM. Once a packet matches in the hardware, it needn't be sent to the kernel for further processing, greatly improving the packet forwarding performance.

For non-IP protocols or the IP-based protocols that are not accommodated in the hardware, we integrate a CPU transceiver module. The module is responsible for transmitting the packet to the CPU directly without looking up routing tables. Whether a packet should be transmitted by the CPU transceiver module is completely determined by software. With the flexibility offered by the CPU transceiver module, it's easy to integrate more flexible software virtual data planes into PEARL.

Action Logic. The result of routing table lookup is the next hop information, including output card number, output port number, the destination MAC address and so on, which is stored in a SRAM-based table. It defines how to process the packet, so it can be considered as an action associated with each entry of the routing tables. Based on the next hop information, this stage performs some decisions such as forwarding, dropping, broadcasting, decrementing TTL or updating MAC address.

Multi-stream DMA Engine. To accelerate the IO performance and to greatly exploit the parallel processing power of multi-core processor, we design a multi-stream high-performance DMA engine in PEARL. It can receive packets of different virtual routers from the network card to different memory regions in the host, and transmit packets in the opposite direction via high-speed PCI Express bus. From a software programmer's perspective, there are multiple independent DMA engines, and the packets of different virtual routers are directed into different memory regions, which is convenient and lockless for programming. Meanwhile, we make a tradeoff between flexibility and high performance of DMA transfer mechanism, and carefully redesign the DMA engine in FPGA. The DMA engine can transfer packets to the pre-allocated huge static buffer at contiguous memory locations. It greatly decreases the number of memory accesses required to transfer a packet to CPU. Each packet transported between the CPU and the network card is equipped with a custom header, which is used for carrying processing information to the destination, such as the matching results.

Output Scheduler. The egress packets sent back by the CPU are scheduled based on their output port number, which is a specific field in the custom header of the packet. Each physical output port is equipped with an output queue. The scheduler puts each packet in the appropriate output queue for transmitting.

C. Software Platform

Our software platform of PEARL consists of several components, as shown in Figure 3. These include vmmd, to provide the basic management functions for the virtual routers and the packet processing cards; nacd, to offer a uniform interface to the underlying processing cards outside the virtual environment; routed, to translate the forwarding rules generated by the kernel or user applications into a uniform format in each virtual router, and install these rules into the TCAM of the processing cards; netio_proxy, to transmit the packets between the physical interfaces and virtual interfaces, and low_proxy, to dispatch packets into low priority virtual routers which shares one

Figure 3 Packet path in the software

There are two types of virtual routers in our PEARL: the high and low priority virtual routers. Each high-priority virtual router is assigned with one pair of DMA Rx/Tx buffers, and an independent TCAM space for lookup. With the high-speed lookup based on TCAM, and efficient IO channels provided by the hardware, the virtual router can achieve the maximum throughput in PEARL platform. For low-priority virtual routers, all the virtual routers share only one pair of DMA Rx/Tx buffers and they can't utilize TCAM for lookup. The macvlan mechanism is adopted to dispatch packets between multiple virtual routers. The applications inside the low priority virtual routers can use the socket APIs to capture packets from the virtual interfaces. Note that each packets needs to go through at least 2 times context switch (system calls) during the transmission to the user application, resulting in a relatively low IO performance.

We take IPv4 and IPv6 as two Internet protocols to show how the virtual routers can be easily implemented on the PEARL.

High-priority IPv4 virtual router.

To create a high-priority IPv4 virtual router in our platform, the vmmd process first starts a new Linux container with several virtual interfaces, and collects the MAC address of each virtual interface, installs these addresses in the underlying cards via the nacd process, so that the hardware can identify the packets heading to this virtual router, and copy the packet into the certain DMA Tx buffer which assigned to this virtual router.

Then, the routed process is launched in the new container. It extracts the routes through the NETLINK socket inside the virtual router and installs routes and the forwarding action in hardware, so the hardware can fill a little structure in the memory to notify the netio_proxy process when a packet match a route in TCAM. The netio_proxy process delivers the packets either to the virtual interface or directly to the physical interface according to the forwarding action in memory. For example, most of time, normal packets will match a route in the hardware. When the netio_proxy receive these packets, it will directly send them through a DMA Tx buffer. An ARP request packet will not match any rules in the TCAM, the netio_proxy process will deliver this packet to the virtual interface, receive the corresponding ARP reply packet from the virtual interface, and then send it to the physical interface.

Low-priority virtual router.

To create low priority virtual routers, a tap device is set up by the vmmd process (tap0).

Low priority virtual routers are configured to share the network traffic of this tap device through the macvlan mechanism. The low_proxy process acts like a bridge, transmitting packets between DMA buffers and the tap device in both directions. Notes that the MAC addresses of the virtual interfaces are generated randomly, we can encode the MAC addresses to identify virtual interface where the packet comes from. For example, we can use the last byte of the MAC address to identify the virtual interfaces, if the low_proxy process receives a packet with source MAC address 02:00:00:00:00:00, it knows that the packet is from the first virtual interface in one of the low priority virtual routers, and transmits the packet to the first physical interface immediately. We adopted this method in the low_proxy process and vmmd process, and use the second byte to identify the different low virtual routers. It not only saves the time consumed by the inefficient MAC hash lookup to determine where the packet comes from, but also saves the space in TCAM, because all the low priority virtual routes only need one rule in TCAM (02:*:00:00:00:*).

IPv4/IPv6 virtual router with kernel forwarding. In this configuration, the routed process does not extract the route from the kernel routing table; instead, it enables the ip_forward options of the kernel. As a result, all packets will match the default route in TCAM without the forwarding action. The netio_proxy process transmits all these packets into the virtual interfaces, so that the kernel will forward the packet instead of the underlying hardware. The tap/tun device is used as the virtual interface. Since the netio_proxy is a user space process, each packet needs two system calls to complete the whole forwarding.

User-defined virtual router. User-define packet process procedure can be implemented as a plug-in loaded by the netio_proxy process, which makes the PEARL extensible. We opened the basic packet APIs to the users, such as read_packet(), send_packet(). Users can write their own process module in C Language, and runs it in the independent virtual routers. For the light-weight applications which do not need to deal with huge amount of network traffic, users can also write a socket program in either high or low priority virtual routers.

EVALUATION AND DISCUSSION

We implemented PEARL prototype using a common server with our specialized network card. The common server is equipped with a Xeon 2.5GHz 64-bit CPU and 16G DDR2 RAM. The OS-level virtualization techniques Linux Containers (LXC) is used to isolate the different virtual routers (VR).

In order to demonstrate the performance and flexibility of PEARL, our implementation is evaluated in three different configurations: high performance IPv4 virtual router, kernel forwarding IPv4/IPv6 virtual router, and IPv4 forwarding in low priority virtual router.

We conducted 4 independent sub-networks with Spirent TestCenter to measure the performance of the three configurations in 1-4 VRs. Three different length of packet (64, 512 and 1518) is used (for IPv6 packet, the packet length is 78, 512 and 1518. 78 bytes is the minimal IPv6 packet length supported by TestCenter).

Figure 4 shows the throughputs of an increasing numbers of VRs using configuration 1. Each virtual data plane has been assigned with a pair of DMA RX/TX buffers, and the independent routing table space in the TCAM, resulting in an efficient IO performance and a high speed IP lookup. The result shows, when the number of the VR reaches to 2, the throughput of minimal IPv4 packet of PEARL is up to 4Gbps, which is the maximum theatrical speed of our implementation. Figure 5 illustrates the throughputs of an increasing numbers of VRs using configuration 2. Each virtual data plane has the same configuration with Configuration 1, except that there is no virtual data plane has its own routing space in TCAM. In Figure 5, the throughput of minimal IPv4/IPv6 packet forwarding is only 100Mbps when there is only one VR. It's because we used the original kernel for the packet forwarding, each packet needs to go through 2 times context switch and 3 times memory copy in our design. We can optimize this by re-implement the forwarding functions as a plug-in in the netio_proxy process in VR. It can be used to verify the applications which handle little traffic (new routing protocols etc). We can see from the results, the total throughputs of the minimal IPv4 packet remains 60Mbps as the number of VR increases. The macvlan mechanism used for the sharing of the network traffic between multiple VRs results in a long kernel path for the packet processing, so the total performance is even lower than the IPv4 kernel forwarding in the Configuration 2. We can improve the performance by developing a new mechanism that suit for our case.

RELATED WORK

Recent researches, such as vRouter Project [START_REF] Egi | Towards high performance virtual routers on commodity hardware[END_REF], RouteBricks [START_REF] Dobrescu | RouteBricks: exploiting parallelism to scale software routers[END_REF] and PacketShader [START_REF] Han | PacketShader: a GPU-accelerated software router[END_REF], have exploited the tremendous power of modern multi-core CPUs and multi-queue network interface cards (NICs) to build high-performance software routers on commodity hardware. Memory latency or IO performance becomes the bottleneck for small packets in such platforms. Due to the complexity of DMA transfer mechanism in commodity NICs, the performance of high speed PCI Express bus is not fully exploited. Meanwhile, the traditional DMA transfer and routing table lookup result in multiple memory accesses. In PEARL, we simplify the DMA transfer mechanism and redesign the DMA engine to accelerate IO performance, and offload the routing table lookup operation to the hardware platform.

OpenFlow [START_REF] Mckeown | OpenFlow: Enabling innovation in campus networks[END_REF] enables rapid innovation of various new protocols, and divides the function into control plane and data plane. OpenFlow provides virtualization through flowvisor, without isolation in hardware. In contrast, PEARL hosts multiple virtual data plane in hardware itself, which could offer both strong isolation and high performance.

Figure 1

 1 Figure 1 Overview of the PEARL architecture Isolation. PEARL implements multiple simultaneous fast virtual data planes by allocating separate hardware resources to each virtual data plane. This facilitates strong isolation among the hardware virtual data planes. Moreover,

Figure 2 :

 2 Figure 2: The architecture of PEARL's hardware data plane VID Lookup. Each virtual router in the platform is marked by a unique VID. This stage classifies the packet based on the fields extracted in the previous stage. We use a TCAM for the storage and lookup of the fields, which can be configured by the users. Due to the special features of TCAM, each field of the rules in VID lookup table can be a

Figure 4

 4 Figure 4 Throughput of high performance IPv4 virtual routers

Figure 5 Figure 6

 56 Figure 5 Throughput of kernel forwarding IPv4/IPv6 virtual routers Figure 6 show the throughputs of an increasing numbers of low priority VRs using configuration 3. Low priority VRs shares only one pair of DMA TX/RX IO channel and can't take advantage of TCAM to speed up the IP Lookup.

Figure 6

 6 Figure 6 Throughput of low priority IPv4 virtual routers

ACKNOWLEDGEMENT

This work was supported by National Basic Research Program of China under grant No. 2007CB310702, by National Natural Science Foundation of China (NSFC) under grant No. 60903208, by NSFC-ANR (Agence Nationale de Recherche, France) under grant No. 61061130562, by the Instrument Developing Project of the Chinese Academy of Sciences under grant No. YZ200926.

SwitchBlade [START_REF] Anwer | SwitchBlade: a platform for rapid deployment of network protocols on programmable hardware[END_REF] builds virtualized data plane in FPGA-based hardware platform. It classified packets based on MAC address. PEARL, on the other hand, dispatches packets into different virtual routers based on the 10-tuple. With the flexibility of the 10-tuple and wildcard, PEARL has the capability to classify packets based on the protocol of any layer. Moreover, in PEARL, all packets are transmitted to the CPU for switch.

The Supercharged Planetlab Platform (SPP) [START_REF] Turner | Supercharging PlanetLab -A High Performance, Multi-Application, Overlay Network Platform[END_REF] is a high performance, multi-application, overlay network platform. SPP utilizes the network processor (NP) to enhance the performance of traditional Planetlab [START_REF] Chun | Planetlab: an overlay testbed for broad-coverage services[END_REF] nodes. SPP divides the whole data plane functions into different pipeline stages on NP. Some of these pipeline stages can be replaced by custom code, resulting in an extensible packet processing. However, the flexibility and isolation of SPP is limited due to the inherent vendor-specific architecture of NP. PEARL takes advantage of the flexibility and full parallelism of FPGA to run multiple high-performance virtual data planes in parallel, while keeping good isolation in fast data path.

CoreLab [START_REF] Corelab | [END_REF] is a new network testbed architecture that supporting full-featured development environment. It enables excellent flexibility by offering researchers a full-virtualized environment on each node for whatever arrangement. PEARL can support the similar flexibility in software data path. Meanwhile, PEARL offers a high performance hardware data path for the performance-critical virtual network applications.

CONCLUSIONS

We aim at flexible routers to bridge the gap between new Internet protocols and the practical test, deployment.

To this end, this work presents a programmable virtual router platform, PEARL. The platform allows users to easily implement new Internet protocols and run multiple isolated virtual data planes concurrently. A PEARL router consists of a hardware data plane and a software data plane with DMA engines for packet transmission. The hardware data plane is built on top of a FPGA based packet processing card with TCAM embedded. The card facilitates fast packet processing and IO virtualization. The software plane is built by a number of modular components and provides easy program interfaces. We have implemented and evaluated the virtual routers on the PEARL.


Pauline Larue

Pierre Jallon
email: pierre.jallon@cea.fr

Bertrand Rivet
email: bertrand.rivet@gipsa-lab.grenoble-inp.fr

MODIFIED K-MEAN CLUSTERING METHOD OF HMM STATES FOR INITIALIZATION OF BAUM-WELCH TRAINING ALGORITHM

published or not. The documents may come

INTRODUCTION

Hidden Markov models (HMM) are widely used for statistical modeling of signals having a temporal structure. They are based on two statistics processes: the state variable process and the observation process. The first one X is the so-called state variable and is a Markov chain, assumed discrete in this paper. It is not observed but it is fully characterized by transitions and initialization probabilities. In many recognition algorithms (e.g., speech [START_REF] Rabiner | A tutorial on hidden markov models and selected applications in speech recognition[END_REF], writing recognition [START_REF] Hu | HMM based online handwriting recognition[END_REF] , gesture recognition [START_REF] Elmezain | A hidden markov model-based continuous gesture recognition system for hand motion trajectory[END_REF], ...), this variable is used to describe the temporal structure of the signals to model. In particular if these signals can be described as sequences of (shorter) stationary signals, a very particular set of models can be used: the left-right models; models whose transition probabilities satisfy the following constraints:

∀(i 1 , i 2 ), i 2 < i 1 , p (X t = i 2 |X t-1 = i 1 ) = 0,
which means that the state variable can only be a, increased sequence. The second process Y is the so-called observation process and is assumed to be an independent process conditionally to X [START_REF] Rabiner | A tutorial on hidden markov models and selected applications in speech recognition[END_REF]. For each state n ∈ {0, • • • , N -1} a probability density function (p.d.f.) p (Y |X = n) is defined, where N is the number of different values taken by X or the number of hidden states. Y is assumed to be a continuous variable and its p.d.f is modeled for each state as a Gaussian distribution.

The following notations will be used in the rest of this paper:

• the initial probability of state is denoted π i = p (X 0 = i), ∀i ∈ {0, • • • , N -1}, • the transition probability a i1,i2 between states i 1 and i 2 is defined by

a i1,i2 = p (X t = i 2 |X t-1 = i 1 ), ∀(i 1 , i 2 ) ∈ {0, • • • , N -1} × {0, • • • , N -1}.
The observation probabilities are described, for state #i with two variables µ i (mean vector) and Σ i (covariance matrix). Finally the whole set of parameters of a HMM with N states is denoted as λ N :

λ N = {π i } i , {a i1,i2 } i1,i2 , {µ i , Σ i } i .
In general λ N is estimated using training sequences.

Given K observation sequences of length T k , Y k,0:T k -1 = {Y k,0 • • • , Y k,T k -1 }, the optimal set of parameters λ N is defined as: λN = argmax λ N K-1 k=0 p Y k,0:T k -1 |λ N .
Without additional assumptions, this problem can not be solved analytically. The mainly used technique to estimate λ N is the expectation maximization (EM) algorithm [START_REF] Dempster | Maximum-likelihood from incomplete data via the EM algorithm[END_REF] trough the forward-backward method so-called Baum-Welch algorithm [START_REF]Inference in Hidden Markov Model[END_REF]. Starting from a preliminary set of parameters λ

N , the algorithm estimates the set of parameters in an iterative manner -denoted as λ (s) N at step #s -such as

K-1 k=0 p Y k,0:T k -1 |λ (s) N increases with respect to s. λN is then estimated as: λN = lim s→∞ λ (s) N .
However convergence to a global maxima is not ensured and this algorithm is well known to be very sensitive to its initialization point λ (0) N [START_REF] Rabiner | Some properties of continuous hidden markov model representations[END_REF]. Several initialization methods have hence been proposed in the literature. A k-mean algorithm on the data can be used to cluster observations Y k,0:T k -1 [START_REF] Rabiner | A tutorial on hidden markov models and selected applications in speech recognition[END_REF][START_REF] Nathan | Initialization of hidden markov models for unconstrained on-line handwriting recognition[END_REF] or several sets of parameters to perform the training can be used [START_REF] Md | A variable initialization approach to the EM algorithm for better estimation of the parameters of hidden markov model based acoustic modeling of speech signals[END_REF]. For left-right models, constrained clustering techniques can also be used with a k-mean algorithm on the observation data [START_REF] Huda | A constraint-based evolutionary learning approach to the expectation maximization for optimal estimation of the hidden markov model for speech signal modeling[END_REF].

In this paper, an alternative approach is proposed to estimate a first set λ (0) N . A set of parameters is first estimated with a very large number of states (M N ) to describe all training sequences in a very accurate way. However, this set of parameters can not be used in practice due to over-learning problems and computational time issues. As a consequence, the number of states has to be reduced to a much smaller number to overcome these difficulties. An unsupervised clustering algorithm (based on the k-mean algorithm) on the observation p.d.f is hence proposed to do this operation while keeping the training signals description accuracy. Transition and initialization probabilities remain to be estimated with the Baum-Welch algorithm. The proposed method is derived in two cases: unconstrained and left-right models.

The paper is structured as follows. Section 2 describes the mainlines of the algorithm. The k-mean algorithm operations are detailed in section 3. The proposed method performance is compared to other methods's based on simulations in section 4. Finally, Section 5 concludes this paper.

INITIALIZATION ALGORITHM DESCRIPTION

The aim of the algorithm is to provide an initialization set of parameters λ

N close enough to λN to ensure a good convergence of the learning algorithm.

It only focuses on the observation probabilities parameters, which aims at finding N sets of observation probabilities parameters {µ

(c) i , Σ (c) i } i∈{1,••• ,N } , i.e. N
Gaussian distributions, to describe in an accurate way the training sequences. It is worth noting that transition probabilities a i1,i2 and initialization probabilities π i are not estimated (although constrained for the left-right models), this operation being done by the Baum-Welch algorithm. The following values are hence used:

• for unconstrained models, transition probabilities are set to a i,i = 0.8, ∀i ∈ {0, • • • , N -1} and ∀i 1 = i 2 , a i1,i2 = 0.2/(N -1), and initialization probabilities

π i = 1/N , ∀i ∈ {0, • • • , N -1}; • for left-right models, transition probabilities are set to a i,i = 0.8, ∀i ∈ {0, • • • , N -2} and a N -1,N -1 = 1, a i,i+1 = 0.
2 and other values are set to 0, finally initialization probabilities

π 0 = 1 and ∀i ∈ {1, • • • , N - 1}, π i = 0.
Concerning the N p.d.f. of hidden states, the estimation is performed in two steps: first M p.d.f. are estimated which accurately describe the signal (Subsection 2.1). These M distributions are then reduced to N (Subsection 2.2).

Initialization step

Given a set of K training sequences, Y k,0:T k -1 , k ∈ {1, • • • , K}, M
Gaussian distributions which accurately describe data can be estimated as follows. Each training sequence Y k is split into several segments of length P , overlapping or not and covering its time support. Each segment is then modeled as a Gaussian signal which parameters µ k,j , Σ k,j are estimated using classical methods: µ k,j being the mean of the related signal segment and Σ k,j its covariance matrix. Without restriction, it is assumed that the distributions are sorted with respect to the Y k segment delay. For instance, with non-overlapping segments, each couple of values (0 ≤ k ≤ K -1 and 0 ≤ j < T k /P ) can be estimated as:

µ k,j = 1 P P -1 p=0 Y k,jP +p , Σ k,j = 1 P P -1 p=0 Y k,jP +p -µ k,j Y k,jP +p -µ k,j T ,
where . T is the transpose operator. P should be chosen small enough so that each signal segment can be considered as stationary, and large enough to ensure a good estimation of mean vector and covariance matrix. It is possible after this step to collect all these values to build the initial set of parameters λ (0)

M . It is worth noting that this initial set λ (0)
M describes in a very accurate way the training sequences: indeed, each training signal segment being described by one of the M states. However, and as already mentioned in the introduction section, this model can not be used in practice because of over-learning problems and of computational time issues.

Clustering step

The second step of the algorithm is hence to approximate these M Gaussian distributions with N ones (N M ). This operation is performed using a k-mean algorithm on these Gaussian distributions.

Let N µ k,j , Σ k,j denotes the Gaussian distribution p.d.f. with mean vector µ k,j and covariance matrix Σ k,j . Moreover, let us refer as 'center' the distribution which characterizes a cluster, of which parameters are labelled by (c): N µ

(c) i , Σ (c) i .
The k-mean algorithm works iteratively as follows.

Given an initial set of N centers N µ

(c) i,0 , Σ (c) i,0 1≤i≤N
, the following two steps are performed at each iteration #s: 1. Association: each distribution N µ k,j , Σ k,j is associated to a center in the set N µ

(c) i,s , Σ (c) i,
(c) i,s+1 , Σ (c) i,s+1 1≤i≤N
.

The N centers are finally estimated as:

∀i , N µ (c) i , Σ (c) i = lim s→∞ N µ (c) i,s , Σ (c)
i,s .

The initial set of N centers depends on the type of HMM (unconstrained or left-right models). For unconstrained models, the M distributions are sorted according to the first component of their mean vector. The set of M distributions is split into N consecutive sets and a center is estimated for each set according to the k-mean algorithm second step as described above. For left-right models, each training sequence distributions are split into N disjointed segments. First segment is then associated to center #0, second one to center #1, and so on. Initialization centers are also estimated using the k-mean algorithm second step.

It is worth noting that the first step requires to adapt the Euclidian distance to associate each distribution to a cluster. This point is discussed in Section 3.1 for unconstrained models and in Section 3.2 for left-right models. Furthermore, the second step is described in Section 3.3.

MODIFIED K-MEAN ALGORITHM

This section details the practical proposed considerations of the proposed modified k-mean clustering algorithm.

State to center distance computation for unconstrained models

The unconstrained model case is first considered: each Gaussian distribution is associated to a center in an independent way. As both distribution and center are Gaussian distributions, the Kullback-Leibler divergence is proposed to estimate the distance between both p.d.f.. This divergence is defined as:

D KL p 1 (u) p 2 (u) = p 1 (u) log p 1 (u) p 2 (u) du.
For Gaussian distributions, this latter expression is proportional to

D KL N µ k,j , Σ k,j N µ (c) i , Σ (c) i ∝ ln det Σ (c) i det Σ k,j + Tr Σ (c) i -1 Σ k,j + µ (c) i -µ k,j T Σ (c) i -1 µ (c) i -µ k,j , (1)
where Tr is the trace operator and det is classically the determinant. Finally, ∀(k, j), the related distribution N µ k,j , Σ k,j is associated to the closest center i * which satisfies:

i * = argmin i D KL N µ k,j , Σ k,j N µ (c) i , Σ (c) i
.

State to center distance computation for left-right models

Compared to previous unconstrained model, left-right assumption sets an additional constraint. For all k, if distribution N µ k,j0 , Σ k,j0 is associated with center #i 0 , then for all j 1 > j 0 , N µ k,j1 , Σ k,j1 cannot be associated to any center #i, i < i 0 . In other words, the distribution -center association has to be jointly done for each training sequence (rather than in an independent way). For sequence #k, let consider therefore the set of N + 1 break indices I k (0), • • • , I k (N ), defining the N states. This set must verify due to the left-right constrain that I k (0) = 0, I k (N ) = M -1 and for two states i 1 and i

2 , ∀ i 1 < i 2 ∈ {0, • • • , N } × {0, • • • , N }, I k (i 1 ) < I k (i 2 ).
Finally, ∀ (k, j), the related distribution N µ k,j , Σ k,j is associated to the closest center i * which satisfies:

i * = i I k (i) ≤ j < I k (i + 1).
The set of break indices is previously estimated as

Îk (1), • • • , Îk (N -1) = argmin I k (1),••• ,I k (N -1) D I k (1), • • • , I k (N -1) ,
where [START_REF] Rabiner | A tutorial on hidden markov models and selected applications in speech recognition[END_REF]. Note that the two extrema I k (0) and I k (N ) do not require to be optimized since I k (0) = 0 and I k (N ) = M -1.

D I k (1), • • • , I k (N -1) = N -1 i=0 j∈I k (i) D KL N µ k,j , Σ k,j N µ (c) i , Σ (c) i , with I k (i) = {I k (i), • • • , I k (i + 1) -1} and D KL (• •) defined by

Center estimation

The third and last issue to solve is the computation of the updated center parameters once each distributions have been associated to one cluster. Consider therefore the center #i and I (c) i the associated set of distributions. A random variable Z in this set can be written as:

Z = (k,j)∈I (c) i δ(H -(k, j))N k,j
where N k,j is a Gaussian random variable with mean and variance {µ k,j , Σ k,j }, δ(u) the Dirac delta function and H is a hidden variable equal to (k, j) if Z shares N k,j distribution.

The i-th center parameters µ

(c) i , Σ (c) i
are estimated as mean and covariance of Z. It is straightforward to check that:

µ (c) i = 1 card I (c) i (k,j)∈I (c) i µ k,j and
Σ (c) i = 1 card I (c) i (k,j)∈I (c) i Σ k,j + µ k,j µ T k,j -µ (c) i µ (c) i T where card I (c) i
is the number of elements in set I (c) i .

NUMERICAL RESULTS

This section presents the results achieved by the proposed method with several kind of configurations: toy example and simulated signals.

Scenario µ 0 , Σ 0 µ 1 , Σ 1 µ 2 , Σ 2 µ 3 , Σ 3 #1 [ 2 2 ], [ 1 0 0 2 ] [ 6 6 ], [ 2 0 0 1 ] [ 1 1 ], [ 1 0 0 2 ] [ 5 5 ], [ 2 0 0 1 ] #2 [ 1 1 ], [ 5 0 0 5 ] [ 3 3 ], [ 1 0 0 1 ] [ 6 6 ], [ 5 0 0 5 ] [ 5 5 ], [ 1 0 0 1 ] #3 [ 1 1 ], [ 5 0 0 5 ] [ 2 2 ], [ 1 0 0 1 ] [ 3 3 ], [ 5 0 0 5 ] [ 4 4 ], [ 1 0 0 1 ]
Table 1: Simulation scenarios.

Evaluation method

The proposed algorithm performance (denoted "k-mean on Gaussian laws") have been estimated using numerical simulations and are compared to the performance of the classical k-mean algorithm applied directly on data (denoted "k-mean on data"). Moreover, the influence of the left-right constrain is investigated (referred as "no constrained" or "left-right"). Nevertheless, constrained k-mean on data for left-right models has not been tested because the number of sets of break indices to test is prohibitive leading thus to an excessive computational time.

As mentioned in the introduction, the aim of these algorithms is to estimate a set of parameters λ(0) N so that the convergence of the Baum-Welch algorithm is improved. This latter algorithm is an iterative algorithm computing in an iterative manner the set of parameters λ (s) N at step #s such as

K-1 k=0 p Y k,0:T k -1 |λ (s) N increases with respect to s.
For a given set of training sequences Y k,0:T k -1 k , the sets of parameters λ(0) N estimated with the different initialization algorithms are hence compared with the value reached and the convergence rate of the following likelihood function

(s) = log K-1 k=0 p Y k,0:T k -1 |λ (s) N .

Toy example

Signals have been generated as the concatenation of 4 segments signals so that Y 0:T = {Y (0)

0:T0-1 , Y (1) 0
:T1-1 , Y (2) 0:T2-1 , Y (3)
0:T3-1 }. Each of these four signals is generated using a Gaussian distribution: for the p-th signal, with mean vector µ p and covariance matrix Σ p . Their time lengths T p , ∀p are randomly chosen using a truncated Gaussian distribution of mean 500 and standard deviation 100 (keeping only positive value).

Three sets of parameters {µ p , Σ p } p have been tested corresponding to 3 scenarios resumed in Table 1. As one can see, these three scenario vary from these difficulty. Indeed, scenario #1 is the easiest since the 4 states are quite separated. Scenario #2 is a little bit more complex than scenario #1 since states 1 and 2, and 3 and 4 overlap. Finally, scenario #3 is the more confusing model since all the state largely overlap. For each scenario, 100 realizations have been generated, and for each realization 5 sequences have been used for the initialization and the training algorithms. P = 20 has been used for estimation of M Gaussian distributions prior to k-mean algorithm.

Fig. 1 presents the achieved results as the averaged values of (s) over the 100 realizations. For the first scenario (Fig. 1 The proposed method hence improves the estimation of the set of parameters used to initiate the Baum Welch algorithm: a better and a faster convergence is shown by simulations.

Recognition problem

In a second set of numerical experiments (Tab. 2) classification accuracy (CA) of the proposed method is compared to CA achieved by the reference 'dynamic time warping' method (DTW) [START_REF] Sakoe | Dynamic programming algorithm optimization for spoken word recognition[END_REF]. Two configurations are compared, and for each of them 50 bi-dimensional observations (Y k,i ∈ R 2 ) are generated from 25 different models. Each model is composed from the concatenation of 2 to 4 states (Configuration #1) and of 4 to 8 states (Configuration #2). In all this experiment (configurations #1 and #2), a 5-state left-right HMM with the proposed initialization procedure is considered (referred as HMM). As one can see, the proposed method outperforms the classical DTW by about 10% of classification accuracy computed in a 10-fold cross validation procedure: i.e. the 50 observations sequences are partitioned into 10 sets and each of them is sequentially used as the test database while the other 9 sets are used to train the HMM parameters. It is worth noting that the good behavior of proposed method is achieved without optimization of the number of states. As a consequence it is not surprising that the best CA is obtained on the simplest configuration (Conf. #1), however the gap between the two detection method remains of about 10% even with the more complex configuration (Conf. #2).

CONCLUSION

In this study, an initialization of the Baum-Welch training algorithm based on a modified k-mean clustering method of HMM states is presented. The proposed procedure differs from classical implementations by a clustering of the states rather than on the training data. The simulations results have shown that the proposed method improves the initialization of the Baum-Welch algorithm since the value of the log-likelihood achieved by our method is higher than value achieved by the classical initialization. Moreover, in a recognition prob-lem, the proposed method outperforms the reference dynamic time warping method showing its good behavior.

Future works will include a deeper analysis of this method as well as an automatic procedure to adjust jointly the number of hidden states.

  Averaged log-likelyhood values no constrained k-mean on gaussian laws Left-right k-mean on gaussian laws no constrained k-mean on data Averaged log-likelyhood values no constrained k-mean on gaussian laws Left-right k-mean on gaussian laws no constrained k-mean on data Averaged log-likelyhood values no constrained k-mean on gaussian laws Left-right k-mean on gaussian laws no constrained k-mean on data (c) Scenario #3

Figure 1 :

 1 Figure 1: Performance evaluation achieved on the toy example.

2 DTWTable 2 :

 22 (a)), the 3 tested algorithms ('left-right kmean on Gaussian distributions', 'unconstrained k-mean on Gaussian distributions' and 'on data') converge to the same value, but k-mean on Gaussian distributions Conf. #1 Conf. #Classification accuracy in percentage [%].algorithms converge much faster than classical k-mean initialization. For the second scenario (Fig.1(b)), the two proposed k-mean on Gaussian distributions algorithms (left-right and unconstrained) converge to the same point with a few iterations whereas the classical k-mean on data algorithm fails into a local maxima. Finally for the most complex scenario (Fig.1(c)), the constrained left-right proposed k-mean algorithm on Gaussian distribution has better performance than both other algorithms.