Andrea Cherubini

Fabien Spindler

Franc ¸ois Chaumette

A New Tentacles-based Technique for Avoiding Obstacles during Visual Navigation

Keywords: Visual Navigation, Visual Servoing, Collision Avoidance

In this paper, we design and validate a new tentacle-based approach, for avoiding obstacles during appearance-based navigation with a wheeled mobile robot. In the past, we have developed a framework for safe visual navigation. The robot follows a path represented as a set of key images, and during obstacle circumnavigation, the on-board camera is actuated to maintain scene visibility. In those works, the model used for obstacle avoidance was obtained using a potential vector field. Here, a more sophisticated and efficient method, that exploits the robot kinematic model, and predicts collision at look-ahead distances, is designed and integrated in that framework. Outdoor experiments comparing the two models show that the new approach presents many advantages. Higher speeds and precision can be attained, very cluttered scenarios involving large obstacles can be successfully dealt with, and the control inputs are smoother.

I. INTRODUCTION

One of the main objectives of robotics research is the development of vehicles capable of autonomously navigating in unknown environments [START_REF] Nunes | Introducing perception, planning, and navigation for Intelligent Vehicles[END_REF], [START_REF] Buehler | Special Issue on the 2007 DARPA Urban Challenge, Part I-III[END_REF]. In this field, an important task is obstacle avoidance: if possible, a collision-free trajectory to the goal should be generated; otherwise, the vehicle must brake to prevent collision [START_REF] Wada | A deceleration control method of automobile for collision avoidance based on driver's perceptual risk[END_REF].

The task that we focus on is outdoor visual navigation: a wheeled vehicle, equipped with an actuated pinhole camera and with a forward-looking range scanner, must follow a path represented by key images, without colliding with the ground obstacles. The camera detects the features required for navigating, while the scanner senses the obstacles (in contrast with other works, such as [START_REF] Kato | An obstacle detection method by fusion of radar and motion stereo[END_REF], only one sensor is used to detect the obstacles). In the past, obstacle avoidance has been integrated in visual navigation [START_REF] Ohya | Vision-based navigation by a mobile robot with obstacle avoidance using a single-camera vision and ultrasonic sensing[END_REF], [START_REF] Yan | Mobile robot indoor navigation using laser range finder and monocular vision[END_REF] and path following [START_REF] Lamiraux | Reactive path deformation for nonholonomic mobile robots[END_REF], [START_REF] Lee | Mobile robot navigation with reactive free space estimation[END_REF], by using the path geometry or the environment 3D model (including, for example, walls and doors). However, since our task is defined in the image space, we seek a merely sensor-based solution, which does not need a global model of the environment and trajectory.

In our recent work [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF] we adopted a method based on potential fields [START_REF] Khatib | Real-time obstacle avoidance for manipulators and mobile robots[END_REF] built on an occupancy grid. This method, however, suffered from its simplicity, which led to very conservative collision avoidance, and strongly varying control inputs. Alternative reactive strategies include: the vector field histogram [START_REF] Borenstein | The Vector Field Histogram -Fast obstacle avoidance for mobile robots[END_REF], dynamic window [START_REF] Fox | The Dynamic Window approach to obstacle avoidance[END_REF], obstacle-restriction method [START_REF] Minguez | The Obstacle-Restriction Method (ORM) for robot obstacle avoidance in difficult environments[END_REF], and closest gap [START_REF] Mujahad | Closest Gap based (CG) reactive obstacle avoidance navigation for highly cluttered environments[END_REF]. The first two methods rely on a candidate set of commands; however, trap situations and oscillations may arise. In [START_REF] Minguez | The Obstacle-Restriction Method (ORM) for robot obstacle avoidance in difficult environments[END_REF], these issues are solved, but, since the robot is assumed holonomic and a 3D subgoal A. Cherubini was with INRIA Rennes -Bretagne Atlantique, IRISA, Campus de Beaulieu 35042, Rennes, France, and is now with the Laboratory for Computer Science, Microelectronics and Robotics LIRMM -Université de Montpellier 2 CNRS, 161 Rue Ada, 34392 Montpellier, France. {Andrea.Cherubini}@lirmm.fr. F. Spindler and F. Chaumette are with INRIA Rennes -Bretagne Atlantique, IRISA. {firstname.lastname}@inria.fr is required, the approach is not suitable for our problem.

Similarly, [START_REF] Mujahad | Closest Gap based (CG) reactive obstacle avoidance navigation for highly cluttered environments[END_REF] uses the robot pose, which is noisy in our framework. Instead, we take inspiration from [START_REF] Bonnafous | Motion generation for a rover on rough terrains[END_REF] and [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF], where a set of trajectories (arcs of circles) is evaluated for navigating. However, in [START_REF] Bonnafous | Motion generation for a rover on rough terrains[END_REF], a sophisticated probabilistic elevation map is used, and the selection of the optimal arc of circle is based on its risk and interest, which both require accurate pose estimation. Similarly, in [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF] the candidate trajectories ('tentacles') are also used for free navigation, which relies on GPS way points, hence -once more -on the robot pose. A deeper comparison with this work is carried out in Sect. IV.

The main contribution of our work is the development of a novel, reactive, pose-independent, tentacle-inspired obstacle avoidance technique, which is perfectly suitable for appearance-based tasks, such as visual navigation. The main strengths of this approach are that it is reactive and sensorbased (it does not require any model of the environment), whilst providing useful "look-ahead" information, by exploiting the robot geometric and kinematic characteristics. Our method is experimentally validated in the framework designed in [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF], and compared with the previous potential field method. In [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF], we guaranteed that obstacle avoidance had no effect on the visual task, by discerning two contexts (safe and unsafe), and pursuing the corresponding tasks. Here, the tentacle-based approach intervenes both in assessing the context, and in designing the task in case of danger.

The article is organized as follows. In Sect. II, the control law from [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF] is recalled. Then, two alternative obstacle models are presented: vortex potential fields (Sect. III) and tentacles (Sect. IV). Experimental results are reported in Section V, and summarized in the conclusion.

II. CONTROL SCHEME

A. General Definitions

The reader is referred to Fig. 1. We define the robot frame F R (R, X, Y ) (R is the robot center of rotation) and image frame F I (O, x, y) (O is the image center); C is the camera optical center. The robot control inputs are:

u = (v, ω, φ) .
These are the translational and angular velocities of the vehicle, and the camera pan angular velocity. We use the normalized perspective camera model, and we assume that the sequence of images that defines the path can be tracked with continuous v (t) > 0. This ensures safety, since only obstacles in front of the robot can be detected by our scanner.

The path that the robot must follow is represented as a database of ordered key images, such that successive pairs contain some common static visual features (points). First, the vehicle is manually driven along a taught path, with the camera pointing forward (ϕ = 0), and all the images are saved. Afterwards, a subset (database) of N key images I 1 , . . . , I N representing the path (Fig. 1(b)) is selected. Then, during autonomous navigation, the current image, noted I, is compared with the next key image I d ∈ {I 1 , . . . , I N }, and a relative pose estimation between I and I d is used to check when the robot passes the pose where I d was acquired. For key image selection, and visual point detection and tracking, we use the algorithm in [START_REF] Royer | Monocular vision for mobile robot localization and autonomous navigation[END_REF]. The output of this algorithm, which is used by our controller, is the set of points visible both in I and I d . Then, navigation consists of driving the robot forward, while I is driven to I d . We maximize similarity between I and I d using only the abscissa x of the centroid of the points matched on I and I d . When I d has been passed, the next image in the set becomes the desired one, and so on, until I N is reached.

Along with the visual path following problem, we consider obstacles which are on the path, but not in the database, and sensed by the range scanner in a plane parallel to the ground. For obstacle modeling, we use the occupancy grid in Fig. 2(a): it is linked to F R , with cell sides parallel to X and Y . Its extension is limited (X m ≤ X ≤ X M and Y m ≤ Y ≤ Y M ), to ignore obstacles that are too far to jeopardize the robot. Any grid cell c centered at (X, Y ) is considered occupied if an obstacle has been sensed in c. For cells entirely lying in the scanner area, only the current scanner reading is considered. For all other cells, we use past readings, which are progressively displaced using odometry.

B. Control Design

The desired behaviour of the robot is related to the surrounding obstacles. When the environment is safe, the vehicle should progress forward while remaining near the taught path, with camera pointing forward (ϕ = 0). If avoidable obstacles are present, we apply a robot rotation for circumnavigation with an opposite camera rotation to maintain visibility. Finally, if collision is inevitable, the vehicle should simply stop. To select the behaviour, we assess the danger at time t with a situation risk function

H : IR * + → [0, 1].
Stability of these tasks has been guaranteed in [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF] by:

       v = (1 -H) v s + Hv u ω = (1 -H) λx(x d -x)-jvvs+λϕj φϕ jω + Hω u φ = H λx(x d -x)-jvvu-jωωu j φ -(1 -H) λ ϕ ϕ (1)
In the above equations:

• H ∈ [0, 1]
is the situation risk function introduced above; two alternative definitions, depending on the obstacle model, will be given in Sect. III and IV. • v s is the translational velocity in the safe context (i.e., when H = 0). It must be reduced when the features are moving quickly in the image, making tracking difficult. This is typically the case at sharp robot turns, and when the camera pan angle ϕ is strong. We define v s as:

v s (ω, ϕ) = v m + v M -v m 4 σ (2)
with function σ defined as:

σ : IR × -π 2 , π 2 → [0, 4] (ω, ϕ) → [1+tanh (π-k ω |ω|)] [1+tanh (π-k ϕ |ϕ|)] .
Function (2) has an upper bound v M > 0 (for ϕ = ω = 0), and smoothly decreases to the lower bound v m , as either |ϕ| or |ω| grow. The decreasing trend of v s is determined by empirically tuned positive parameters v M , v m , k ω and k ϕ . This definition of v s yields better results than the one in [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF], which was only characterized by the image x variation.

• v u ∈ [0, v s ] is the translational velocity in the unsafe context (H = 1). It guarantees that the vehicle slows down (and eventually stops) in dangerous situations. • x and x d are abscissas of the feature centroid respectively in the current and next key image. • λ x > 0 and λ ϕ > 0 are empirical gains determining the convergence trend of x to x d and of ϕ to 0. • j v , j ω and j φ are the components of the Jacobian relating ẋ and the control inputs u:

j v = -sin ϕ+x cos ϕ ζ j ω = RC(cos ϕ+x sin ϕ) ζ + 1 + x 2 j φ = 1 + x 2 ,
with RC and ζ (the centroid depth, hand-tuned according to the environment characteristics) depicted in Fig. 1

(a).

• ω u is the angular velocity that makes the robot avoid collisions while advancing. In [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF], we proved that (1) is well defined by setting ζ > RC/2, and that ∀H ∈ [0, 1] obstacle avoidance has no effect on the visual task. Our main contribution here will be in the obstacle model, i.e., in the definition of the variables H, v u and ω u according to the danger. Two alternative models will be presented and confronted experimentally. To carry out the comparison, the vortex potentials [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF] are briefly recalled in Sect. III. Then, our novel tentacle-based approach is designed and presented in Sect. IV.

III. VORTEX POTENTIAL FIELDS

This approach is illustrated in Fig. 2(a). Given an arbitrary integer K, obstacles are modeled using the latest 2K + 1 scans. For each cell c = (X, Y ), we define the 2K + 1 occupancies r at the j-th oldest iteration as:

r j (c) = {0, 1}, j = 0, . . . , 2K + 1.
We set r j = 1 if an obstacle has been sensed in c at the j-th iteration prior to the current one, and 0 otherwise. Then, we associate to each cell a coefficient µ (c), obtained by linear combination of the occupancies, weighted with a normalized Gaussian filter that smoothens the cell effect over time:

µ (c) = 2K+1 j=0 e -(j-K) 2 /K √ Kπ r j (c) .
The maximum weight is at the K-th latest scan, to avoid overshoot at a new obstacle detection. If the robot velocity is negligible with respect to the scanner frequency, and K is small, the effect of motion on the occupancies can be ignored. In [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF], we showed that these assumptions are appropriate in our setup.

The potential associated to each cell c = R, is defined as:

U c = µ (c) c ,
where c is the distance from R to c. We define the vortex field for each cell as the rotor of U c :

f c = f c,X f c,Y =    ± ∂U c ∂Y ∓ ∂U c ∂X    = µ (c)      ∓ Y c 3 ± X c 3     
.

The signs of f c,X and f c,Y depend on X, so that the field always points forward, and the fields f c,i generated by all cells are superimposed to obtain the total field:

f = i f c,i .
We then use the magnitude and phase of f (denoted |f | and ∠f ), and two tuned thresholds ∆ d and ∆ s such that 0 < ∆ d < ∆ s , to design the situation risk function as:

H=              0 if |f | = 0 k f |∠f | if 1 |f | ≥ ∆ s 1 if 1 |f | ≤ ∆ d 1+k f |∠f | 2 + 1-k f |∠f | 2 tanh ∆ d 1-|f |∆ d + ∆s 1-|f |∆s otherwise.
Note that H = 0 if no obstacle is detected, and it is bounded by 1. For small |f |, H is determined by ∠f : the obstacles are far enough to be circumnavigated, and parameter k f ∈ 0, 2 π weighs the danger provoked by the field orientation. Instead, for large |f |, the obstacles are 'too near', thus H = 1. A hyperbolic function is used to interpolate in between.

In the presence of dangerous obstacles (i.e., for large H), the robot should slow down, and eventually stop. Imposing v = 0 when H = 1 in (1), yields:

v u = 0.
Hence, in the general case, the first equation in [START_REF] Nunes | Introducing perception, planning, and navigation for Intelligent Vehicles[END_REF] becomes:

v = (1 -H) v s .
The velocity is reduced from v s to 0, as the risk increases.

To ensure that to circumnavigate the obstacles in dangerous contexts, the robot aligns its heading with f , we use a gain λ f > 0, and set the desired unsafe angular velocity to:

ω u = λ f ∠f .

IV. TENTACLES A. Related work

Our new technique is mainly inspired from [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF], although some differences, listed below, have been introduced, to deal with the specific constraints of appearance-based navigation.

The design in [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF], although strengthened by a dynamic analysis, is justified only empirically. For example, the danger of the cells is associated to their orthogonal projection on the arc of circle, without measuring the actual distance to collision. Instead, in the case of two cells with same orthogonal projection (e.g., c 1 and c 3 in Fig. 2(e)), the internal cell should be considered more dangerous. Besides, in [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF], the requested trajectory is defined by GPS way points. Hence, in the absence of GPS, the approach in [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF] becomes purely reactive, while our approach can still follow a path, defined by the image database.

In contrast with that work, our approach does not require a 3D notion of the goal. It must be purely reactive, since the visual task and the obstacle avoidance task are defined in different state spaces (respectively in the image and in the local planar surroundings). As in [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF], our grid is built locally at each new scanner acquisition, without accumulating past data, and each tentacle is associated to some classification areas overestimating the vehicle encumbrance. However, these areas are defined and used differently. Their definition is associated to the rigid body kinematics of the boxes representing the vehicle encumbrances, to consider distances to collision. Then, with the largest areas we select the safest tentacle and its risk, and with the thinnest one, the eventual deceleration. Finally, since our vehicle navigates at slightly varying speeds, we do not relate the tentacle sets to v (as in [START_REF] Hundelshausen | Driving with tentacles-Integral structures of sensing and motion[END_REF]), so that a reduced number of candidate paths is sufficient. Our approach operates locally and instantaneously, without planning nor deriving the robot pose, to determine the values of H, v u and ω u in (1).

B. Classification areas and metrics

We use, along with the set of all occupied grid cells:

O = {c 1 , . . . , c n } ,
a set of drivable paths (tentacles). Each tentacle j is a semicircle that starts in R, is tangent to X, and is characterized by its curvature (i.e., inverse radius) κ j , which belongs to K, a uniformly sampled set:

κ j ∈ K = {-κ M , . . . , κ M } .
The maximum desired curvature κ M > 0, must be feasible considering the robot kinematics. We consider an odd number of tentacles, so that a straight tentacle (with κ = 0) also exists. To illustrate our method, in Fig. 2(b-e), the straight and the sharpest counterclockwise (κ = κ M ) tentacle are shown in dashed black. When a total of 3 tentacles is used, these correspond respectively to j = 2 and j = 3.

Each tentacle j is characterized by three classification areas (collision, dangerous central, and dangerous external), obtained by rigidly displacing, along the tentacle, three rectangular boxes (red, blue, and dark green in Fig. 2), with increasing width. The boxes, with same height and different widths, are all overestimated with respect to the robot dimensions, to ensure that, in the presence of disturbances, the actual path is included in the classification area. All three areas are delimited by the box and by three arcs of circle (or lines, in the particular case κ = 0) concentric with the tentacle, and starting from three points (denoted A, B and C) on the box perimeter. In all cases, these points are: the two corners on one side of the box (the outer side must be used if κ = 0), and the intersection between the other side (intern if κ = 0) and the rear wheel axis. For the collision areas, A, B and C are shown in Fig. 2. We then associate each area to the set of all cells (pink, cyan and light green in Fig. 2) whose center lies within the area. For tentacle j, the three sets of cells are noted C j , D j and E j . The only exception is with cells belonging to the dangerous central set, which are not considered in the external set as well:

D j E j = ∅.
During navigation, O is used, along with the sets just defined, to calculate a candidate risk function H j ∈ [0, 1] for each tentacle j, and select the best tentacle accordingly. Then, the unsafe translational velocity on the best tentacle, v u , is calculated to adapt the speed to the potential danger, and to finally derive ω u . All these steps are detailed below.

C. Tentacle risk function

The tentacle risk function H j is derived from the risk distance of all occupied cells in the dangerous areas. This distance is denoted ∆ ij ≥ 0 for each c i ∈ O (D j E j ).

For occupied cells in the central set D j , ∆ ij is the distance that the middle boundary box (blue) would cover along tentacle j before touching the cell center.

For the external set, we consider only the subset Ēj ⊆ O E j of cells which reduce the clearance in the tentacle normal direction. For each external occupied cell, we denote Γ ij the ray starting at the tentacle center and passing through c i . Cell c i is added to Ēj if and only if, in D j E j , there is at least an occupied cell crossed by Γ ij on the other side of the tentacle. In the example of Fig. 2(e), O E 3 = {c 1 , c 3 , c 4 }, whereas Ē3 = {c 1 , c 3 }. Then, for cells in Ēj , ∆ ij is the sum of two terms: the distance from the center of c i to its normal projection on the perimeter of the dangerous central area, and the distance that the middle boundary box would cover along tentacle j before reaching the normal projection. The derivation of ∆ ij is illustrated, in Fig. 2, for 4 cells. Note that for a given cell, ∆ ij may have different values (or even be undefined) according to the tentacle.

When all risk distances on the tentacle are calculated, we compute ∆ j as their minimum:

∆ j = inf ci∈(O∩Dj )∪ Ēj ∆ ij . If (O D j ) Ēj ≡ ∅, ∆ j = ∞.
In the example of Fig. 2, ∆ 2 = ∆ 12 and ∆ 3 = ∆ 33 . Obviously, overestimating the bounding box sizes leads to more conservative ∆ j .

We then use ∆ j and two tuned thresholds ∆ d and ∆ s (0 < ∆ d < ∆ s ), to design the tentacle risk function:

H j =        0 if ∆ j ≥ ∆ s 1 2 1 + tanh 1 ∆j -∆ d + 1 ∆j -∆s if ∆ d < ∆ j < ∆ s 1 if ∆ j ≤ ∆ d .
(3) Note that H j smoothly varies from 0, when the dangerous cells on the tentacle (if any) are far, to 1, when they are near. If H j = 0, the tentacle is tagged as clear. All the H j s are compared (with the strategy explained below), to determine H in (1) and select the best tentacle for navigation.

D. Situation risk function and best tentacle

Here we detail our strategy for determining the best tentacle curvature κ b for navigation, and therefore u in [START_REF] Nunes | Introducing perception, planning, and navigation for Intelligent Vehicles[END_REF].

Initially, we calculate the path curvature κ = ω/v ∈ IR that the robot would follow if there were no obstacles. Replacing H = 0 in (1), it is:

κ = λ x x d -x -j v v s + λ ϕ j φϕ /j ω v s ,
which is always well-defined, since j ω = 0 and we have set v s > v m > 0. We obviously constrain κ to the interval of feasible curvatures [-κ M , κ M ]. Then, we derive the two neighbors of κ among all the existing tentacle curvatures:

κ n , κ nn ∈ K such that κ ∈ [κ n , κ nn ) .
Let κ n be the nearest one, i.e., the curvature of the tentacle that best approximates the safe path 1 . We denote it as the visual task tentacle. The situation risk function H v of that tentacle, which measures the risk on the visual path, is obtained by linear interpolation of its neighbours:

H v = (H nn -H n ) κ + H n κ nn -H nn κ n κ nn -κ n . ( 4
)
If H v = 0, the visual task tentacle can be followed: we set κ b = κ n , and we apply [START_REF] Nunes | Introducing perception, planning, and navigation for Intelligent Vehicles[END_REF] with H = 0. Instead, if H v = 0, we seek a clear tentacle (H j = 0). First, to avoid abrupt control changes, we only search among the tentacles between the visual task one and the best one at the previous iteration2 , noted κ pb . If many clear ones are present, the closest to the visual task tentacle is chosen. If none of the tentacles with curvature in [κ n , κ pb ] is clear, we search among the others. Again, the best tentacle will be the clear one that is closest to κ n and, in case of ambiguity, the one closest to κ nn . If a clear tentacle has been found, we select it and set H = 0. Instead, if no tentacle in K is clear, the one with minimum H j calculated using (3) is chosen, and H is set equal to that H j . Ambiguities are again solved first with the distance from κ n , then from κ nn .

In all cases, the unsafe translational velocity v u is derived from the obstacles on the best tentacle, as explained below.

E. Unsafe translational velocity

The unsafe translational velocity is derived from the collision distance on the best tentacle, δ b , which is a conservative approximation of the maximum distance that the robot can travel along the best tentacle without colliding. Since the thinner (red) box contains the robot, if R follows the best tentacle, collisions can only occur in occupied cells in C b .

In fact, the collision with cell c i will occur at the distance, denoted δ ib ≥ 0, that the red box would cover along the best tentacle, before touching the center of c i . The derivation of δ ib is illustrated in Fig. 2 for four occupied cells.

Then, we define δ b as the minimum among the collision distances of all occupied cells in C b : The translational velocity must be designed accordingly. Let δ d and δ s be two tuned thresholds such that 0 < δ d < δ s . If the probable collision is far enough (δ b ≥ δ s ), v can be maintained at the safe value defined in [START_REF] Buehler | Special Issue on the 2007 DARPA Urban Challenge, Part I-III[END_REF]. Instead, if it is near (δ b ≤ δ d ), the robot should stop. To comply with the boundary conditions v u (δ d ) = 0 and v u (δ s ) = v s , in between we apply a constant deceleration:

δ b = inf
a = v 2 s /2(δ d -δ s ) < 0.
Since the distance required for braking at velocity v u (δ b ) is:

δ b -δ d = -v 2
u /2a, the expression of the unsafe translational velocity becomes:

v u (δ b ) =    v s if δ b ≥ δ s v s δ b -δ d /δ s -δ d if δ d < δ b < δ s 0 if δ b ≤ δ d , (5)
in order to decelerate as the collision distance decreases.

F. Unsafe angular velocity

Once the best tentacle (with curvature κ b ) is chosen, and the corresponding v u is obtained with (5), we set:

ω u = v u κ b . (6)
Setting this value of unsafe angular velocity in (1) guarantees that when H = 1, the robot precisely follows the best tentacle, with translational velocity v u .

G. Discussion

The values of H, v u and ω u derived as explained above are inserted in (1) to derive the control inputs u.

When H = 0, the robot tracks at its best the taught path: the image error is regulated by ω, while v is set to v s to improve tracking, and the camera is driven forward (ϕ = 0). This occurs if the 2 neighbour tentacles are clear, while in the vortex approach, even a single occupied cell would generate H > 0. Thus, one advantage of the new approach is that only obstacles on the visual path are taken into account. When H = 1, φ ensures the visual task, and the two other inputs guarantee that the best tentacle is followed: ω/v = κ b . The applied curvature fluctuates less than with the potential fields, where it is driven by strongly varying f .

In general, the robot navigates between the taught path, and the best path considering obstacles. Only the transition, but not the speed, is driven by H. In fact, note that, for all H ∈ [0, 1], when δ b ≥ δ s : v = v s . This is another advantage of tentacles: a high velocity can be applied if the path is clear up to δ s . With vortex fields, instead, the vehicle stops frequently (as soon as H = 1). Finally, one may object that processing can be costly. However, since ∆ ij and δ ij are invariant geometric characteristics related to the cell positions, their values at each cell can be computed and stored offline, for use when the cell becomes occupied. Control inputs using tentacles in scenario A.1 (top) and A.2 (bottom): v (black, in ms -1 ), ω (green, in rads -1 ), and φ (red, in rads -1 ). The iterations with strong H are highlighted in yellow.

V. EXPERIMENTS Here, we report the experiments (also shown in the video attached to this paper) that we performed to validate our approach. All experiments have been carried out on our CyCab robot, equipped with a 70 • field of view, B&W Marlin (F-131B) camera mounted on a TRACLabs Biclops Pan/Tilt head, and with a 2-layer, 110 • scanning angle, laser SICK LD-MRS. The grid is built by projecting the readings from the 2 layers on the ground, and by using: First, we have compared the two obstacle avoidance methods (vortex potential fields and tentacles) in the experiments in Fig. 3. The taught path (denoted A, and white in the figure) is 60 m long. By placing various obstacles, we design the scenarios A.1 (with a long obstacle perpendicular to the path) and A.2 (with two avoidable obstacles and a blocking barrier). Then, we attempt to replay path A in each scenario, using either technique. The replayed paths, estimated from odometry and snapshots, are drawn in black (potentials) and red (tentacles). The control inputs u using tentacles are plotted in Fig. 4, with dangerous iterations (i.e., with strong H) highlighted in yellow. The smooth trend of u at the beginning and end is due to the acceleration saturation carried out by the CyCab low-level control. In Fig. 5, we have plotted the curvatures ω/v that are applied in the 4 experiments. Since our preliminary tests, it was clear that tentacles could be implemented at higher speeds than potentials, since they provide look-ahead information, and a stabler model. Hence, in (2), we set: v M = 1 ms -1 with tentacles, v M = 0.4 with potentials, and v m = 0.3 for both. For v M > 0.4, potentials fail to avoid the obstacles, while with tentacles, the velocity has been limited to reduce the motion of features between successive images; the maximum speed attainable by the CyCab is 1.3 ms -1 anyway.

X M = Y M =
In scenario A.1, the long obstacle cannot be circumnavigated using potentials. It induces oscillations on the orientation ∠f , hence on the applied curvature (solid black curve in Fig. 5), until the robot is too near to the obstacle, and eventually stops. Instead, using tentacles, the obstacle is overtaken on the left, while the camera rotates right to maintain scene visibility (green and red curves in Fig. 4, top). The robot successfully reaches the final key image and completes navigation, although it is driven over 5 meters away from the taught path. In practice, soon after the obstacle is detected, tentacles with first positive (5 -16 s), and then negative (16 -25 s) curvatures are selected. Since these tentacles are clear, v is reduced only for visual tracking, by (2) (black curve in Fig. 4, top). This is an advantage over the too conservative potential fields, which do not consider the real collision risk, and force the stop. After 25 s, the environment returns clear (H = 0), so the visual tentacle can be followed again, and the robot is driven back to the path. Then (38 -52 s) a small bush on the left triggers H and causes a rotation along with a slight decrease in v. Then the context returns safe, and the visual path can be followed until the end. The translational velocity averaged over the experiment is 0.79 ms -1 , which is more than twice the speed reached in [START_REF] Cherubini | Visual navigation with obstacle avoidance[END_REF].

In scenario A.2, with both methods Cycab navigates without colliding and stops at the barrier. However, as aforementioned, the average navigation velocities v are very different, i.e., 0.83 ms -1 with tentacles, and only 0.35 ms -1 with potentials. Moreover, the applied curvature varies less in the first case (see dashed curves in Fig. 5). This leads to smoother and faster navigation than with potentials. When the barrier is reached, the vortex method behaves as in A.1 to force the stop. Instead, the new method seeks a feasible path until all tentacles are occupied. Then, (5) makes v gradually decrease to zero (black curve in Fig. 4, bottom).

After these results, which confirmed the advantages of the new over the old method, we have run some experiments, uniquely with tentacles, on longer and more crowded paths, shown in Fig. 6(a). The Cycab completed all paths (including 650 m long path E), while dealing with various natural and unpredictable obstacles, such as parked and driving cars, pedestrians, and even bicycles. A first result that emerged from the experiments was that, by assessing the collision risk only along the visual path, non-dangerous obstacles (e.g., walls or cars parked on the sides) are not taken into account. This is clear from Fig. 6(b): the cars parked on the right (which were not present during teaching) do not belong to any of the visual task tentacle classification areas. Hence, they are irrelevant, and do not deviate the robot from the path. Another nice behavior is shown in Fig. 6(c): if a stationing car is unavoidable, the robot decelerates and stops with (5), but, as soon as the car departs, it gradually accelerates (again with ( 5)), to resume navigation. An experiment with a crossing pedestrian is presented in Fig. 6(d). The pedestrian is considered irrelevant, until it enters the visual task tentacle. Then, the clockwise tentacles (blue in the figure) are selected to pass between the person and the right sidewalk. When the path is clear again, the robot turns left to recover it.

VI. CONCLUSIONS

We have presented a novel, robust and reactive technique for avoiding obstacles with a wheeled robot. By exploiting the robot kinematics, we can predict collisions at look-ahead distances along candidate circular paths (tentacles). Since our method is sensor-based and pose-independent, it is perfectly suited for visual navigation. Extensive experiments show that it generates smoother control inputs than its predecessor, and that it can be applied in realistic situations. To our knowledge, this is the first time that outdoor visual navigation with obstacle avoidance is carried out at approximately 1 ms -1 on over 500 m, using neither GPS nor maps. Perspective work includes using a more sophisticated obstacle model (e.g., with shape, dimension and velocity), to design an optimal/complete rather than a worst case algorithm.

Fig. 1 .

 1 Fig. 1. General definitions. (a) Top view of the robot (orange), with actuated camera (blue). (b) Current and next key images, and key image database.

Fig. 2 .

 2 Fig. 2. Obstacle models with 4 occupied cells c 1 , . . . , c 4 (purple). (a) Singular cell (blue vectors) and total (red vector) vortex potential fields. (b -e) Tentacles (dashed black), with classification areas (collision in pink, dangerous central in cyan, dangerous external in green), corresponding boxes and delimiting arcs of circle (red, blue and dark green), and cell risk and collision distances (∆ ij , δ ij ); ray of cell c 1 with respect to the third tentacle, Γ 13 (dashed red).

  If all cells in C b are free, δ b = ∞. In the example of Fig. 2, assuming the best tentacle is the straight one (b = 2), δ b = δ 12 . Again, oversizing the box leads to more conservative δ b .

Fig. 3 .

 3 Fig.3. Two obstacle scenarios with taught (white) and replayed (red: using tentacles, black: using potentials) paths.

Fig. 4 .

 4 Fig. 4.Control inputs using tentacles in scenario A.1 (top) and A.2 (bottom): v (black, in ms -1 ), ω (green, in rads -1 ), and φ (red, in rads -1 ). The iterations with strong H are highlighted in yellow.

  10 m, X m = -2 m, Y m = -10 m. The cells have size 20 × 20 cm. We use 21 tentacles, with κ M = 0.35 m -1 (the Cycab maximum applicable curvature).

Fig. 5 .

 5 Fig. 5. Applied curvature ω/v (in m -1 ) in scenario A.1 (solid) and A.2 (dashed) using tentacles (red) and potentials (black).

Fig. 6 .

 6 Fig. 6. Map of the four navigation paths B, C, D, E (a), and validation with: irrelevant obstacles (b), traffic (c) and a moving pedestrian (d). The visual task tentacle is shown in red, and the best tentacle (when different) in blue; only cells that can activate H (i.e., cells at distance ∆ < ∆s) have been drawn, with brightness proportional to ∆.

We consider that intervals are defined even when the first endpoint is greater than the second: [κn, κnn) must be read (κnn, κn] if κn > κnn.

At the first iteration, we set κ pb = κn.

ACKNOWLEDGMENT

This work has been supported by ANR (French National Agency) CityVIP project under grant ANR-07 TSFA-013-01.


Etienne Côme

Latifa Oukhellou
email: oukhellou@u-pec.fr

Thierry Denoeux

Patrice Aknin

Fault diagnosis of a railway device using semi-supervised independent factor analysis with mixing constraints

Keywords: Independent factor analysis, Semi-supervised learning, Mixing constraints, Mixture models, Maximum likelihood, Fault diagnosis

Independent factor analysis (IFA) defines a generative model for observed data that are assumed to be linear mixtures of some unknown non-Gaussian, mutually independent latent variables (also called sources or independent components). The probability density function of each individual latent variable is modelled by a mixture of Gaussians (MOG). Learning in the context of this model is usually performed within an unsupervised framework in which only unlabelled samples are used. Both the mixing matrix and the parameters of latent variable densities are learned from observed data. This paper investigates the possibility of estimating an IFA model in a noiseless setting when two kinds of prior information are incorporated; namely, constraints on the mixing process and partial knowledge on the cluster membership of some training samples. Semi-supervised or partially-supervised learning frameworks can thus be handled. The investigation of these two kinds of prior information was motivated by a real-world application concerning the fault diagnosis of railway track circuits. Results from both this application, simulated data are provided to demonstrate the capacity of our approach to enhance estimation accuracy and remove the indeterminacy commonly encountered in unsupervised IFA, such as source permutations.

Introduction

Generative latent variable models aim to describe observed variables in terms of a set of unobservable (or latent) variables. Depending on the assumptions made on the latent and observed variable distributions, different kinds of models can be distinguished, including principal component analysis (PCA) (21; 32), factor analysis (FA) (6; 17), and independent component analysis (ICA) (22; 7; 9). This paper deals with a particular model of this last family recently proposed by (3; 26) in the context of signal processing known as independent factor analysis (IFA).

The generative model involved in IFA assumes that observed variables are generated by a linear mixture of independent and non-Gaussian latent variables, as in the ICA model. Furthermore, it assumes that each individual latent variable has its own distribution, which is modelled semi-parametrically by a mixture of Gaussians (MOG). The general IFA model includes data that are noisy, or it can be considered in its noiseless setting. In the latter case, the model is similar to standard ICA but it includes a mixture of Gaussians model for the source densities.

Learning in the context of the IFA and ICA models is often considered within an unsupervised or blind framework. The model parameters and latent variables are learned exclusively from the observed data. These models yield reliable results, provided that the independence assumption is satisfied and the postulated mixing model is suited to the physics of the system. Otherwise, they may fail to recover the sources. Several extensions of the basic ICA model have been proposed to improve its performance, they are commonly referred to as semi-blind independent component analysis (see [START_REF] Babaie-Zadeh | Semi-Blind Approaches for Source Separation and Independent component Analysis[END_REF] for a review). The main approaches exploit nonlinear mixtures [START_REF] Honkela | Blind Separation of Nonlinear Mixtures by Variational Bayesian Learning[END_REF], temporal correlation (4), positivity [START_REF] Cichocki | Nonnegative matrix and tensor factorizations: Applications to exploratory multi-way data analysis and blind source separation[END_REF] or sparsity (19; 23; 34). An approach to modelling class-conditional densities based on an IFA model was also recently introduced in [START_REF] Montanari | Independent factor discriminant analysis[END_REF].

The method proposed here falls in this class of semi-blind approaches and may take advantage of two kinds of prior information through two extensions of the basic noiseless IFA model. The first one concerns the possibility of incorporating independence hypotheses with respect to a subset of latent and observed variables. Such hypotheses can be derived from physical knowledge available on the mixing process. This kind of approach has not been applied within the framework of IFA, but it has been widely considered in factor analysis [START_REF] Bartholomew | Latent variable models and factor analysis[END_REF] and, more specifically, in the domain of structural equation modelling [START_REF] Bollen | Structural Equations with Latent Variables[END_REF]. The second extension incorporates additional information on the cluster membership of certain samples to estimate the IFA model. In this way, semi-supervised learning can be addressed.

Considering the graphic model of IFA shown in Figure 1, the mixing process prior consists in omitting some connections between some set of observed variables X and some set of latent variables Z. The second prior means that additional information on the discrete latent variables Y is taken into account. Indeed, it can be seen that the IFA model provides two levels of interpretation corresponding to discrete and continuous latent variables. For each one of the S latent variables, the discrete latent variable Y si encodes the cluster from which each sample i is drawn. Therefore, partial labelling over these discrete latent variables may help to improve the performance. Such an approach was very recently also explored by [START_REF] Wang | Non-Negative Semi-Supervised Learning[END_REF] in the context of non-negative matrix factorisation using quite different tools, but it is up to our knowledge the first attempt to use such information in the context of IFA. This article is organised as follows. We first review the noiseless IFA model and show how it can be estimated using maximum likelihood. In Section 3, the problem of learning in the context of the IFA model with prior knowledge on the mixing process is addressed. Section 4 focuses on the problem of incorporating additional information on the cluster membership of some training samples to learn the model. A generalised maximum likelihood criterion is defined, and the algorithm for its optimisation is also detailed. Experiments with simulated and standard data sets are presented in Section 5. The railway track circuit diagnosis problem is then addressed in Section 6, and Section 7 concludes the paper.

2 Noiseless IFA 2.1 Background on noiseless IFA IFA was introduced by (26; 3). It originates from both ordinary factor analysis (FA) in applied statistics (29; 31) and independent component analysis (ICA) in signal processing (22; 7). IFA aims to recover independent latent variables from their observed linear mixtures. The latent variables are assumed to be mutually independent and non-Gaussian. In the noiseless form that is used throughout this paper, the IFA model can be expressed as:

x = A z, (1)
where A is a square matrix of size S × S, x is the random vector where elements (x 1 , . . . , x S ) are the mixtures and z is the random vector where elements (z 1 , . . . , z S ) are the latent variables. In this noiseless setting, the un-mixing matrix plays an important role, and we note that W = A -1 in this paper. Thanks to the noiseless setting, a deterministic relationship between the distributions of observed and latent variables can be expressed as:

f X (x) = 1 | det(A)| f Z (A -1 x), (2)
Unlike the ICA model in which the probability density functions of the latent variables are fixed using prior knowledge or according to some indicator that allows switching between sub-and super-Gaussian densities [START_REF] Hyvärinen | Independent Component Analysis[END_REF], each latent variable density in IFA is modelled as a mixture of normally-distributed components (i.e., mixture of Gaussians MOG) so that a wide class of densities can be approximated (26; 3):

f Zs (zs) = Ks k=1 π s k ϕ(zs; µ s k , ν s k ), (3)
Note that ϕ(.; µ, ν) denotes a normal density function with vector mean µ and variance ν. The Equation (3) means that each latent variable is described as a mixture of Ks Gaussians with mean µ s k , variance ν s k and mixing proportions π s k .

Noiseless IFA model Learning

The learning problem associated with the IFA model consists in estimating both the mixing matrix A and the MOG parameters from the observed variables alone. Consider an i.i.d random sample of size N . Using the Equation (2) under the latent variable independence hypothesis, the log-likelihood has the form:

L(A; X) = -N log(| det(A)|) + N i=1 S s=1 log f Zs (A -1 x i )s , (4)
By substituting the density distribution by its expression given in (3), the log-likelihood can be rewritten as:

L(ψ; X) = -N log(| det(A)|) + N i=1 S s=1 log Ks k=1 π s k ϕ (A -1 x i )s, µ s k , ν s k . ( 5
)
Note ψ is the IFA parameter vector ψ = (A, π 1 , . . . , π S , µ 1 , . . . , µ S , ν 1 , . . . , ν S ), A is the mixing matrix; π j is the vector of mixing proportions of source j that sum to 1; and µ j and ν j are vectors of size K j containing the means and the variances of each cluster, respectively. All IFA parameters can therefore be estimated by maximising the likelihood function. When the latent variable densities are known, the mixing matrix estimation is based on gradient methods that maximize the likelihood. The stochastic gradient of the log-likelihood defined in (5) can be derived as:

∆A ∝ ∂L(A; X) ∂A ∝ -N.(A -1 ) t + N i=1 (A -1 ) t g A -1 .x i (A -1 x i ) t ∝ -(A -1 ) t + 1 N N i=1 (A -1 ) t g (z i ) z t i ∝ (A -1 ) t 1 N N i=1 g (z i ) z t i -I , (6)
where

z i = A -1 x i and g(z) = -∂ log(f Z 1 (z1)) ∂z1 , . . . , -∂ log(f Z S (zS )) ∂zS t
. The update rule of the mixing matrix is thus given by:

A (q+1) = A (q) + τ ∆A (q) , (7)
Note that τ is the gradient step that can be adjusted by means of line search methods such as the backtracking line search method (27, page 37), which decreases the learning rate τ by a factor ρ ∈ [0, 1] until sufficient increase conditions are met. The convergence of this algorithm can be improved by using the natural gradient, which is based on the geometrical structure of the parameter space. This can be obtained by multiplying the left-hand side by the matrix AA t , which leads to the following mixing matrix update rule [START_REF] Hyvärinen | Independent Component Analysis[END_REF]:

∆ nat A = AA t ∆A = A 1 N N i=1 g (z i ) z t i -I . (8)
with z

(q) i = W (q) x i and I the identity matrix. Maximization of the likelihood with respect to all model parameters can be achieved by an alternating optimisation strategy. The gradient algorithm is indeed well suited to optimise the log-likelihood function with respect to the mixing matrix A when the parameters of the source marginal densities are held constant. In contrast, when A is fixed, an EM algorithm can be used to optimise the likelihood function with respect to the parameters of each source. These observations have led to the development of a generalised EM algorithm (GEM) that simultaneously maximises the likelihood function with respect to all model parameters (16; 24).

Constraints on the mixing process

This section investigates the possibility of incorporating independence hypotheses concerning relationships between some subset of latent and observed variables in the ICA model. Such hypotheses are often deduced from physical knowledge of the mixing process. The hypothesis that we consider in this section has the following form, namely, X h ⊥⊥ Zg, which means that X h is statistically independent from Zg. This kind of hypothesis constrains the form of the mixing matrix, as shown by the following proposition: Proposition 1 In the noiseless ICA model,

X h ⊥⊥ Zg ⇔ A hg = 0. ( 9
)
Proof. Independence can be defined as:

X h ⊥⊥ Zg ⇔ f X h ×Zg (x h , zg) = f X h (x h ) f Zg (zg). (10)
In the case of the noiseless IFA model, the joint probability density function on X h × Zg is given by:

f X h ×Zg (x h , zg) = R S-1 f X h ×Z1ו••×ZS (x h , z 1 , . . . , z S ) S l=1,l =g dz l (11) = R S-1 S s=1 f Zs (zs) × δ(x h -A h. z) S l=1,l =g dz l = f Zg (zg) ×   R S-1 S l=1,l =g f Z l (z l ) × δ(x h -A h. z)dz l   , (12)
Note that A h. is the h th row of the mixing matrix A. Using the independence assumption, we identify:

X h ⊥⊥ Zg ⇔ f X h (x h ) = R S-1 S l=1,l =g f Z l (z l ) × δ(x h -A h. z)dz l , (13)
Note that δ is the Dirac function. The integral must not depend on zg (the g th row of z), which is possible only if A satisfies A hg = 0.

The log-likelihood must be maximised under the constraint that some of the mixing coefficients are null, and the gradient ascent is only performed with respect to the non-null coefficients. In this case, the initialisation and the update rules of the mixing matrix are given by:

A (0) = C • A (0) A (q+1) = A (q) + τ C • ∆A (q) , (14)
Note that • denotes the Hadamard product between two matrices, and C is a binary matrix defined by C hk = 0 if Z k ⊥⊥ X h , and C hk = 1, otherwise.

4 Semi-supervised learning in noiseless IFA

The expectation maximization (EM) algorithm provides a general solution to problems involving missing data [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF]. Here, we are interested in one of its most classical applications, which concerns mixture estimation problems.

In a semi-supervised learning context, the IFA model is built from a combination of M labelled and N -M unlabelled samples (12; 13). The maximum likelihood criterion can be decomposed into two parts respectively corresponding, to the supervised and unsupervised learning examples; the log-likelihood criterion then becomes:

L(A; X) = -N log(| det(A)|)+ M i=1 S s=1 Ks k=1 l s ik log π s k ϕ (A -1 x i )s, µ s k , ν s k + N i=M+1 S s=1 log Ks k=1 π s k ϕ (A -1 x i )s, µ s k , ν s k . (15)
Note that l s ik ∈ {0, 1} Ks , l s ik = 1 if sample i comes from component c k of sources s and l s ik = 0 otherwise.

A generalised EM algorithm (GEM), (also noted here as Algorithm 1) can be designed to simultaneously maximise the likelihood function with respect to all model parameters. This algorithm is therefore similar to the EM algorithm used to estimate IFA parameters in an unsupervised setting, except for the E-step during which the posterior probabilities t s ik are only computed for the unlabelled samples. The score function g of each latent variable density is thus given by:

gs(z is ) =      -∂ log Ks k=1 l s ik π s k ϕ(zis;µ s k ,ν s k ) ∂zis , if i ≤ M -∂ log Ks k=1 π s k ϕ(zis;µ s k ,ν s k ) ∂zis , if i > M =    Ks k=1 l s ik (zis-µ s k ) ν s k , if i ≤ M Ks k=1 t s ik (zis-µ s k ) ν s k , if i > M (16)
Note that t s ik is the posterior probability that the sample i belongs to component k of the latent variable s, given z is = (W x i )s and the labels:

t s ik = π s k ϕ(z is ; µ s k , ν s k ) Ks k ′ =1 π s k ′ ϕ(z is ; µ s k ′ , ν s k ′ ) . ( 17
)
The proposed algorithm thus takes as inputs a centred observation matrix X (of size N × S), the available labels : l s ik , ∀i ∈ {1, . . . , M }, ∀s ∈ {1, . . . , S} and ∀k ∈ {1, . . . , Ks} and a constraints matrix (a binary matrix of size S × S) with zeros indicating the desired independence assumptions.

Algorithm 1 Pseudo-code for Semi-supervised IFA with GEM algorithm Inputs:Centred observation matrix X, cluster membership for the M labelled data l s ik , constraints matrix C # Random initialization of IFA parameter vector ψ (0) , q = 0 while convergence test do //Latent variable update

Z = X A (q) -1 t
//Update of the latent variable parameters (EM) for all s ∈ {1, . . . , S}, and k ∈ {1, . . . , Ks} do //E-Step t

s(q) ik = l s ik , ∀i ∈ {1, . . . , M } t s(q) ik = π s(q) k ϕ(z is ; µ s(q) k , ν s(q) k ) Ks k ′ =1 π s(q) k ′ ϕ(z is ; µ s(q) k ′ , ν s(q) k ′ ) , ∀i ∈ {M + 1, . . . , N }
end for for all s ∈ {1, . . . , S}, and k ∈ {1, . . . , Ks} do //M-step, Update of the parameter vector of each latent variable π

s(q+1) k = 1 N N i=1 t s(q) ik µ s(q+1) k = 1 N i=1 t s(q) ik N i=1 t s(q) ik z is ν s(q+1) k = 1 N i=1 t s(q) ik N i=1 t s(q) ik (z is -µ s(q+1) k ) 2
end for //Update of the score matrix G ( 16)

G = g (q+1) (Z) //Natural gradient (8) ∆A = A (q) 1 N N i=1 g z (q) i z (q) i t -I //Linear search τ (gradient step) τ * = Linearsearch(A (q) ,∆A) //mixing matrix update A (q+1) = A (q) + τ * .C • ∆A //Latent variable normalization to remove scale indetermination for all s ∈ {1, . . . , S} do σ 2 s = Ks k=1 π s(q+1) k (ν s(q+1) k + µ s(q+1) k 2 ) - Ks k=1 π s(q+1) k µ s(q+1) k 2 for all k ∈ {1, . . . , Ks do µ s(q+1) k = µ s(q+1) k /σs ν s(q+1) k = ν s(q+1) k /σ 2 s A (q+1) s. = A (q+1) s.
/σs end for end for q ← q + 1 end while Outputs:Estimated parameters : ψ ml , estimated latent variables : Z ml

Remark 1

The convergence of the algorithm can be checked in different ways. The increment in log likelihood can be monitored and a test such as:

L(ψ q ; X) -L(ψ (q-1) ; X) |L(ψ (q-1) ; X)| < ǫ, ( 18
)
where ǫ is a precision threshold set to a small value (10 -6 ), can be used to check the convergence. Another convergence test can be based on differences between successive estimates of the parameters.

In the following sections, we study the performance of semi-supervised noiseless IFA with mixing constraints using both simulated and real data sets. We compare our results with the results of classic unsupervised IFA. The experiments were designed to illustrate the capability of the proposed approach to enhance estimation accuracy and remove the indeterminacy commonly encountered in unsupervised IFA, such as indeterminacy related to permutation of sources.

Experiments using simulated and standard real data

To better understand our approach as compared to the unsupervised IFA model, different experiments were carried out to show the influence of incorporating prior information on estimation accuracy as well as on the complexity of the optimisation problem. The first set of experiments was designed to show the impact of mixing constraints on estimation accuracy. In the second set of experiments, we investigate the benefits of learning the IFA model when information regarding the component membership of some training samples is introduced. We show that such information can be exploited to enhance estimation accuracy and remove indeterminacy commonly encountered in unsupervised IFA, such as the indeterminacy due to permutation of sources.

Several simulated data sets were built as follows. Six independent latent variables were considered with the densities shown in Figure 2. The mixtures were then generated using the IFA model given in [START_REF] Amari | A New Learning Algorithm for Blind Signal Separation[END_REF], where each coefficient of the 6 × 6 mixing matrix was randomly generated according to a standard normal distribution. We aimed to recover six latent variables from six observed ones.

Influence of mixing constraints on estimation accuracy

We now present simulation results that highlight the practical usefulness of introducing mixing constraints when estimating the IFA model. The mixing coefficients were randomly generated according to a standard normal distribution, and zeros were randomly introduced in the mixing matrix; each coefficient has a probability equal to 0.3 of being null. The experiment aims to compare the results of the standard noiseless IFA and the constrained IFA in which prior information on the mixing process is incorporated. The comparison was made on data sets of different sizes N ∈ {100, . . . , 700} and assessed using the Amari performance index [START_REF] Cichocki | Adaptive Blind Signal and Image Processing[END_REF]. This index is a classical measure of performance in the context of ICA; it measures the closeness of the estimated un-mixing matrix to the inverse of the mixing matrix, which is assumed to be known:

Ap(Q) = 1 S S i=1     S j=1 |Q ij | max(|Q i. |) -1   +   S j=1 |Q ji | max(|Q .i |) -1     , (19)
Q = A W , (20)
Note that A is the exact mixing matrix and W = A -1 is the estimated un-mixing matrix. The lower the Amari index is, the more accurate the latent variable estimation is. The Amari indexes obtained for different sample sizes and for the two noiseless IFA models (with constraints or without constraints) are presented in Figure 3. The results were averaged over 30 different training sets, from which 25 random starting points were used for the GEM algorithm. Only the best solution according to likelihood was kept. The results show that constrained IFA (IFAc) outperforms IFA regardless of the sample size; the index improvement is indeed approximately 0.5. and the IFA with prior information on mixing matrix (IFAc).

Influence of labelling on estimation accuracy

The next experiment aims to illustrate the influence of the number of labelled samples on IFA model performance. At first, three different learning strategies were compared, namely, unsupervised IFA, semisupervised IFA using 25% of labelled samples for three sources and semi-supervised IFA using 25% of observations labelled over all sources. Figure 4 shows the performance measures obtained for these strategies when 500 simulated samples are used during the learning phase. The performances were quantified using the absolute correlation between the true sources and their estimates calculated based on a test set of 5, 000 samples. It can be seen that permutation is avoided for sources when labelled samples are provided; correlation is close to 1 in the diagonal terms. This result is because, in the semi-supervised setting, the likelihood is not invariant under permutations of the sources, as the numbering of the source is fixed by the known labels. This avoids the a posteriori analysis usually employed to identify the correspondences between data and sources. Visually, one also sees that the results are improved when semi-supervised IFA is employed, as a better contrast in the correlation matrix is achieved. Figure 5 shows the Amari performance index (1) and the Pearson's correlation coefficient r 2 between the true latent variables and their estimates when the number of labelled samples over all sources varies between 5% and 50%. These results were computed on 30 different learning data sets of 500 samples each. Fifty initialisations were performed for the GEM algorithm, and only the best solution according to likelihood was kept to avoid the problem of local maxima. From these figures, we can see the benefits of incorporating prior knowledge on labels in the estimation of the IFA model. As expected, when the number of labelled samples increases (i.e., to greater than 20%), the model behaves better, as both the mean correlation and the Amari performance index are significantly improved. When there is not enough prior information (i.e., less than 20%) provided to the model, the variability of the r 2 indicator seems to be important, which is mainly because the GEM algorithm converges to local maxima corresponding to the permutations of sources.

Influence of labelling on the likelihood landscape

The following experiment was designed to show the potential benefit of incorporating labels in terms of simplification of the optimisation problem. Figure 6(a) displays the computation time required for the GEM algorithm convergence as a function of the amount of labelled samples. These results were computed on 30 different learning data sets of 500 samples each. It can be seen that the computation time (or the number of iterations) decreases when the amount of labelled samples increases. The number of different local maxima found by the GEM algorithm was also counted for 100 random initialisations. Two local maxima were supposed to be distinct if the distance between the vectors containing all the parameters was larger than 0.001. This experiment was repeated 10 times with different data sets of 500 samples each, generated using the IFA model with the latent variable densities shown on Figure 2. The labelling concerns all the sources. Figure 6

Comparison on standard real datasets

In this subsection, the proposed method is compared with alternative approaches on standard data sets. We investigate the possibility of projecting original data onto a lower dimensional space, when using semisupervised IFA model with prior knowledge on the cluster membership of some training samples.

Experimental setup

Three real datasets of which the characteristics are given on Table 1 were considered. The first dataset is the Crabs dataset1 that contains different morphological measurements for crabs of different species and sex. The second dataset is the YaleB face database [START_REF] Georghiades | From Few to Many: Illumination Cone Models for Face Recognition under Variable Lighting and Pose[END_REF]. We chose the first 5 subjects from the dataset and got a total of 320 face samples captured under different illumination conditions (5 subjects × 64 illumination conditions). This dataset required a preprocessing step before using semi-supervised dimensionality reduction techniques. A model having pixel gray levels as input variables would have 32256 variables (one for each pixel). For this reason, we first processed the data by PCA and kept only the 40 leading principal components.The third dataset deals with the recognition of Thyroid disease using 5 biological measurements [START_REF] Coomans | Comparison of multivariate discriminant techniques for clinical data -application to the tyroid functional state[END_REF]. Table 1 Characteristics of real datasets.

In order to perform feature extraction with semi-supervised IFA, one has to specify the model structure including the number of significant independent factors and the number of Gaussian components modeling each latent variable density. We propose to use two strategies depending on the classification problem to be solved. If the labels can be decomposed into different classification problems, then each latent variable corresponds to a classification problem and the number of Gaussian components is given by the number of subclasses within the corresponding problem. For example, the crabs dataset can be categorized into four classes corresponding in fact to two binary classification problems concerning sex and species of the crabs. In such a situation it's natural to use two latent variables modeled as a mixture of two components (male and female for the sex variable and blue and orange for the species variable). When the classification problem cannot be decomposed in such way (as it is the case for the other two datasets), K latent variables can be used to solve a classification problem of K classes. In this case, each one of the K latent variables is modeled by a mixture of two components: one component encodes a sample's membership to a specific class, the second one encoding the remaining classes, in a one against all scenario. Finally, note that the noiseless IFA model requires recovering as many latent variables as observed ones, therefore the remaining latent variables are considered to be noise components and are modeled by a unique Gaussian distributions.

Results and discussion

Different learning strategies were compared, namely, principal component analysis (PCA), unsupervised IFA, semi-supervised discriminant principal component analysis DPCA [START_REF] Sun | A new discriminant principal component analysis method with partial supervision[END_REF], semi-supervised IFA using some labelled samples over all latent variables and the following supervised methods:LDA, QDA and independent component discriminant analysis ICDA [START_REF] Amato | Independent Component Discriminant Analysis[END_REF]. The DPCA algorithm optimizes a discriminant objective function composed of the PCA criterion and some discriminant criterion estimated on labelled data. However, in order to be used, one has to specify two regularization coefficients (λ, η), where λ controls the balancing of contributions of the two criterions (PCA and discriminant) while η is used to balance the contribution of between-class and generalized within-class scatter matrices. When PCA or DPCA was applied to reduce dimensionality, the nearest neighborhood (1-NN) classifier was employed for classification. 70% of data were randomly chosen for training, while the remaining data are used as a test set. The labelled data were selected at random from the training set. The performances were quantified using the correct classification rates according to each class calculated on a test set. The process was repeated for 30 times and the averaged results were recorded. Finally to avoid local minima, twenty random starting points were used for the GEM algorithm and only the best solution according to the likelihood was kept. Table 2 summarizes the results for supervised and unsupervised methods in order to supply baselines for comparisons with the two semi-supervised approaches compared here in Table 2.

With only 40% of labelled samples, the correct detection rate for the crabs dataset reaches 91.2%, which is close to the performance obtained by supervised approaches (94% on average). Note that unsupervised IFA correctly classifies only 37.5% of observations. PCA-1NN gives 63.1% of correct classification rate. For this dataset, DPCA does not significantly improve the performance. It gives 68.5% of correct classification rate with 40% of labelled samples. For the DPCA method, we set λ = 1 and η = 0.1 as suggested by the authors of the method. With 40% of labelled samples, classification performance for the YaleB dataset reaches 62.6% with our approach and 36.5% with DPCA (20.6% in the unsupervised IFA). For this dataset, we set λ = 1 and η = 0.1. Considering the thyroid dataset, semi-supervised IFA reaches 88.7% of correct classification with 40% of labelled samples. For this specific dataset, PCA and DPCA models outperforms semi-supervised IFA. For DPCA parameters, we set λ = 1 and η = 10. These experimental results show that efficient dimensionality reduction can be achieved in some problems for which unsupervised methods fail to capture the underlying structure of the data. The amount of labelled data required to reach a satisfactory accuracy depends on the discrimination problem complexity. For all the problems investigated here, we found that with only 10% of labelled data semi-supervised IFA gives much better results than standard IFA.

The aim of the next experiment is to highlight the practical interest of semi-supervised learning of the IFA model when labelling errors are present. Considering the crab dataset with 40% of labelled samples, the label noise was varied from 0 to 40%. Figure 7 shows the correct classification rate (%) averaged over 30 random test sets as a function of label noise. When the quantity of labelling errors provided to the model is small (i.e., less than 20%), the proposed approach does not suffer from label noise. We can also notice that the results are still better than those of unsupervised learning context even if the labelling error rate increases. --) and for semi-supervised IFA when the label noise varies between 0 to 40% with 40% of labelled samples. 70% of the samples were used for training, the remaining ones being used for testing.

Fault Diagnosis in Railway Track Circuits

The application considered in this paper concerns fault diagnosis in railway track circuits. This device is first described and the problem to be addressed is then explained. Finally, an overview of the proposed diagnosis method is presented.

The track circuit principle

The track circuit is an essential component of the automatic train control system (15; 28). Its main function is to detect the presence or absence of vehicle traffic within a specific section of railway track. On French high-speed lines, the track circuit is also a fundamental component of the track and vehicle transmission system. It uses a specific carrier frequency to transmit coded data to the train regarding, for example, the maximum authorised speed on a given section on the basis of safety constraints. The railway track is divided into different sections. Each one of them has a specific track circuit consisting of the following components (see Figure 8):

-A transmitter connected to one of the two section ends, which delivers a frequency modulated alternating current; -The two rails that can be considered as a transmission line; -At the other end of the track section, a receiver that essentially consists of a trap circuit used to avoid the transmission of information to the neighboring section; -Trimming capacitors connected between the two rails at constant spacing to compensate for the inductive behavior of the track. Electrical tuning is then performed to limit the attenuation of the transmitted current and improve the transmission level. The number of compensation points depends on the carrier frequency and the length of the track section. The rails themselves are part of the track circuit, and a train is detected when its wheels and axles short-circuit the track. The presence of a train in a given section induces the loss of the track circuit signal due to shorting by train wheels. The drop in the received signal below a preset threshold indicates that the section is occupied. In order to make the transmitted information specific to each track section and to minimize the influence of both longitudinal interference and transverse crosstalk, four frequencies are used for adjacent and parallel track circuits. Neighbouring track circuits are also isolated electrically by using tuned circuits (capacitor and inductance) on both the transmitter and the receiver.

The different parts of the system are subject to malfunctions that must be detected as soon as possible in order to maintain the system at the required safety and availability levels. In the most extreme case, a malfunction an unfortunate attenuation the signal that to the stopping of the train. The purpose of diagnosis is to inform maintainers about track circuit failures based on analysis of a current as recorded by an inspection vehicle. This paper will focus on trimming capacitor faults that affect their capacitance. Figure 9 shows an example of the inspection signal when the system is fault-free, while the other figures correspond to a defective 5 th capacitor and defective 5 th and 9 th capacitors. The aim of the diagnosis system is to detect faults in the track circuit and localise the defective capacitor by analysing the measurement signal.

Overview of the Diagnosis Method

The track circuit can be considered as a large-scale system made up of a series of spatially related subsystems that correspond to the trimming capacitors. A defect on one subsystem is represented by a continue value of the capacitance parameter. The proposed method is based on the following two observations (see Figure 9):

-The inspection signal has a specific structure, which is a succession of so many arches as capacitors; an arch can be approximated by a quadratic polynomial ax 2 + bx + c; -Each observed arch is influenced by the capacitors located upstream, that is, on the transmitter side.

The proposed method consists in extracting features from the measurement signal and building a generative model as shown in Figure 10, in which each observed variable X is corresponds to coefficients (b is , c is ) of the local polynomial approximating the arch located between two capacitors. Only two coefficients are used because of continuity constraints between each polynomial, as there exists a linear relationship between the third coefficient and the three coefficients of the previous polynomial. The continuous latent variable Z is of the generative model is linked the capacitance of the i th capacitor and the discrete latent variable Y is corresponds to the categorisation of the capacitor's state into one of the three states, namely, fault-free, minor defect, and major defect. As there is no influence between a trimming capacitor state and the inspection signal located upstream from it, some connections between latent and observed variables are omitted. This information can also introduced in the model estimation by using constraints on the mixing matrix. We can clearly see that this model is closely linked to the IFA model represented in Figure 1. Considering the diagnosis task as a blind source separation problem, the IFA model can be used to estimate the mixing matrix A and the MOG parameters for each source density. When these parameters are estimated it's possible to recover the latent components (i.e., capacitances of each capacitor) and the estimated capacitor's states from the observed variables alone. The capacitances are estimated by computing A -1 x i and the most probable capacitor state is estimated by setting it to the source component with maximum a posteriori probability.

Moreover, mixing constraints and component membership of some training samples may also be considered during the learning phase. The signals considered here are obtained for only one frequency (four frequencies are used for adjacent and parallel track circuits). The experiments given at this frequency will indeed show the number of observations that have to be labelled in order to obtain satisfactory performance. For the other frequencies, the idea is to label afterwards only this quantity of observations since this task is very time-consuming.

As already explained, a piecewise approach is adopted for signal representation; each arch is approximated by a second degree polynomial in which two coefficients are used as observed variables for each node in the model of Figure 1, resulting in 2 * S observed variables. Given an observation matrix, the aim is to recover S latent variables from 2 * S observed variables with the hope that they will be strongly correlated with the variable of interest, i.e., capacitance. As prior information on the mixing matrix is available, PCA cannot be used for preprocessing because the mixing structure would be lost. Note that 2 * S latent variables are therefore extracted; S latent variable densities corresponding to capacitances are assumed to be mixtures of three Gaussian components, one for each state of the capacitors, while the S other variables are assumed to be noise variables and are thus modelled by simple Gaussian distributions.

Remark 2 We note that with standard IFA model, S latent variables can be recovered from 2 × S observed ones. However, in this paper, we consider the noiseless IFA model, which is a straightforward way to incorporate prior information and to recover sources from the data. In this case, the number of latent variables must be equal to the number of observed variables.

Results and discussion

To assess the performance of this approach, we considered a track circuit of S = 18 subsystems (i.e., capacitors) and a database containing 2, 500 noised signals obtained for different values of the capacitance of each capacitor. We note that 500 signals were used in the training phase, while an additional 2, 000 signals were employed for the test phase. The experiments aim to illustrate the influence of both the number of labelled samples and the use of the mixing matrix constraints on the results.

The model provides two levels of interpretation corresponding to discrete and continuous latent variables. We first discuss the results for the continuous latent variables and then discuss results for the discrete latent variables.

Three different learning settings are compared in this section :

- Furthermore, the number of labelled data was varied between 0 and 500. Note that, when the number of labelled samples is equal to zero, semi-supervised IFA without constraints boils down to the standard IFA model (without any prior), which has indeed very poor performance in this task. Figure 10(a) shows the mean of the absolute value of the correlation between estimated latent variables and capacitances as a function of the number of labelled training samples when the mixing matrix is both constrained and not constrained. When more labelled samples are used, the performances reach a more satisfactory level. For all the settings, twenty random starting points were used for the GEM algorithm and only the best solution according to the likelihood was kept. This figure clearly highlights the benefit of using constraints when the amount of labelled samples is small. As expected, when the number of labelled data increases, the mean correlation also increases to reach a maximal value of 0.84 for the constrained IFA model with 250 labelled sampled and for the unconstrained one with 350 labelled samples. When a sufficient amount of labelled samples is provided to the model (> 350), the prior on the mixing process does not significantly improve the performances. It can also be noticed that unlabelled samples improve the performances of the approach, particularly when there are few labelled data available (< 300). In such a case, the unlabelled data are of great help in the estimation process as is it can be seen in Figure 11(a). Indeed, the contribution of unlabelled samples can be highlighted by comparing the performance of supervised IFA to those of semi-supervised IFA when the mixing matrix is not constrained. For example, with 250 labelled samples, the improvement in term of absolute correlation is around 0.5. Figure 11(b) shows the boxplot of the computation time as the proportion of labelled samples increases from 0 to 500. It can be seen that both the mean value and the variability of the computation time decrease when the amount of labelled samples increases. In the following section, the results are analysed according to the discrete latent variables. With regards to the three classes (i.e., fault-free ω 0 , minor defect ω 1 , and major defect ω 2 ), the confusion matrix between the true and the estimated classes for all capacitances observed in the data set is reported in Table 4 when 30% of samples are labelled. The results show that a good classification performance is obtained despite some confusions between neighbourhood classes (i.e., between ω 0 and ω 1 and ω 1 and ω 2 ). A detection matrix can be also computed by combining the two defect cases (i.e., ω 1 ∪ ω 2 ), which allows for an evaluation of the usual indicators such as false alarm (FA, F A = nb (ω0,D1∪D2) /nbω 0 = 1.55%), correct detection (CD, CD = nb (ω1∪ω2,D1∪D2) /nb (ω1∪ω2) = 92.43%) and false detection (FD, F D = nb (ω0,D1∪D2) /nb D1∪D2 = 24.22%) rates.

In order to accurately quantify the effect of the proportion of labelled samples, the evolution of both false alarm and correct detection rates is presented as a function of the number of labelled samples when the mixing matrix is either constrained or not constrained (Figure 12). When the amount of labelled samples provided to the model is small (i.e., less than 100), the false alarm rate is too high to exploit the results. However, when the proportion of labelled samples is sufficiently large to obtain a weak false alarm rate (i.e., 1.5%), the benefit of prior information on the mixing constraint and the cluster membership of some samples is clearly visible. With only 250 labelled samples, 92% correct detection is reached. Further improvement of the overall performance level would require a non-linear model. Correct detection rate False alarm rate # Fig. 12 Evolution of false alarm (a) and correct detection (b) rates when the number of labelled samples varies between 0 and 500 for semi-supervised IFA with (--) and without constraints (-). The rates have been evaluated on the test set.

Conclusions

In this paper, we have proposed a method for learning parameters of the noiseless IFA model that incorporates two kinds of prior information related to the mixing process, on the one hand, and cluster membership of some subset of training samples, on the other hand. The main assumption of this generative model is that the observed data have been generated by a linear mixture of latent variables assumed to be independent and non-Gaussian. Moreover, each individual latent variable has its own distribution, which is modelled semi-parametrically by a mixture of Gaussians (MOG). In this context, a criterion can be defined, and a GEM algorithm dedicated to its optimisation can be subsequently described. The main limitations of the approach concerns the hypotheses underlying the IFA model: if the mixing system is not linear or if the sources are not independent the method is likely to fail. However, the proposed method has been successfully applied to both simulated data and real data related to fault diagnosis in railway track circuits. The diagnosis system aims to recover the latent variables linked to track circuit defects using features extracted from the inspection signal and IFA hypothesis. A comparison between standard and proposed IFA models has been carried out to show that our approach takes advantage of prior information, thus significantly improving estimation accuracy and removing indeterminacy associated with unsupervised IFA, such as indeterminacy related to the permutation of sources. Further studies should be carried out to take into account imprecise and uncertain cluster memberships, such as those yielded subjectively by human experts. Indeed, soft labels seem to be appropriate when dealing with the imprecise knowledge that experts have regarding the monitoring data of complex systems and to deal with the label noise problem [START_REF] Côme | Learning from partially supervised data using mixture models and belief functions[END_REF].

Fig. 1

 1 Fig.1Graphical model for Independent Factor Analysis.

Fig. 2

 2 Fig. 2 Simulated sources densities.

Fig. 3

 3 Fig.3The Amari performance index as a function of the training sample size obtained using standard IFA (IFA) and the IFA with prior information on mixing matrix (IFAc).

Fig. 4

 4 Fig. 4 Absolute correlation between estimated and latent variables, computed on a test set for (a) unsupervised IFA, (b) semi-supervised IFA with 25% labelled samples for three sources (b) and (c) semi-supervised IFA with 25% of samples labelled over all sources (c).

Fig. 5

 5 Fig. 5 Influence of the proportion of labelled samples on the estimation of the semi-supervised IFA model: Boxplot of the Amari perfomance index and correlation coefficient function of the percentage of labelled samples.

Fig. 6

 6 Fig.6Influence of the proportion of labelled samples on the likelihood landscape: Boxplot of the computation time function of the percentage of labelled samples and number of local maxima detected over 100 random initializations for the GEM algorithm, as a function of labelled samples.

Fig. 7

 7 Fig.7Evolution of correct detection rate of the crab dataset for unsupervised IFA (--) and for semi-supervised IFA when the label noise varies between 0 to 40% with 40% of labelled samples. 70% of the samples were used for training, the remaining ones being used for testing.

Fig. 8

 8 Fig. 8 Diagram of a track circuit and example of a denoised inspection signal.

Fig. 9

 9 Fig. 9 Polynomial parametrization of inspection signal.

Fig. 10 A

 10 Fig. 10 A generative model for the diagnosis of track circuits represented by a graphical model including both continuous and discrete latent variables.

  supervised IFA without constraints : in this setting only labelled training data are used to learn the IFA model. -semi-supervised IFA without constraints : in this setting labelled data are completed with the remaining training samples (over the 500 training signals) this part of the training dataset is unlabelled. -semi-supervised IFA with constraints : in this settings both labelled and unlabelled training data are used and the constraints on the mixing matrix are also used.

Fig. 11

 11 Fig. 11 Influence of prior information : (a) Results of semi-supervised IFA with constraints (-), without constraints (--) and of supervised IFA without constraints (-. -. ) when the number of labelled samples varies between 0 and 500. (b) Boxplot of the computation time function of the number of labelled samples for semi-supervised IFA without constraints.

Table 2

 2 Correct detection rates (%) averaged over 30 random test sets, for five dimensionality reduction methods PCA, IFA, LDA, QDA and ICDA. For all the methods, 70% of the samples were used for training, the remaining ones being used for testing.

			Unsupervised	Supervised
			PCA	IFA	LDA	QDA	ICDA
		Crabs	63.1	37.56	94.67	93.83	92.92
		YaleB	25.8	20.59	96.77	87.53	82.63
	Thyroid	93.38	24.53	91.77	96.41	95.35
					Semi-supervised
	% labelled samples		10	20	30	40	50	60	70
	Crabs	IFA	85.28	89.5	88.5	91.22	90.28	91.56	91.5
		DPCA	67.56	67.61	67.65	68.49	68.45	69.63	69.27
	YaleB	IFA	32.47	47.47	60.07	62.57	68.19	74.10	79.27
		DPCA	33.25	34.18	35.06	36.5	39.43	38.53	42.48
	Thyroid	IFA	74.90	86.77	87.19	88.70	89.06	89.22	89.53
		DPCA	92	93.37	95.28	95.92	96.06	96.17	96.07

Table 3

 3 Correct detection rates (%) averaged over 30 random test sets, for two semi-supervised dimensionality reduction methods IFA and DPCA with different proportions of labelled samples (%). For all the methods, 70% of the samples were used for training, the remaining ones being used for testing.

Table 4

 4 Confusion matrix between true classes (ω 0 ,ω 1 ,ω 2 ) and their estimates (D 0 ,D 1 ,D 2 ), computed on a test set of 2000 track circuits in which 30% of samples have been labelled.

		ω 0	ω 1	ω 2
	D 0	98.45	11.24	1.73
	D 1	1.33	75.25	25.18
	D 2	0.22	13.51	73.09

http://rweb.stat.umn.edu/R/library/MASS/html/crabs.html