Nisrine Ghadban

Paul Honeine

Farah Mourad-Chehade

Clovis Francis

Joumana Farah

Mobility using first and second derivatives for kernel-based regression in wireless sensor networks

This paper deals with the problem of tracking and monitoring physical phenomena using wireless sensor networks. It proposes an original mobility scheme that aims at improving the tracking process. To this end, a model is defined using kernelbased methods and a learning process. The sensors are given the ability to move in a manner that minimizes the approximation error, and thus improves the efficiency of the model. First and second derivatives of the approximation error are used to define the new positions of the nodes. The performance of the proposed method is illustrated in the context of monitoring gas diffusion with wireless sensor networks.

I. INTRODUCTION

Wireless Sensor Networks (WSNs) constitute a multidisciplinary research field at the crossroads of communities such as microelectronics, wireless communication, and signal processing. Originally developed for military purposes, their scope of application continues to expand by including the monitoring in civil and industrial domains [START_REF] Conner | Workplace applications of sensor networks[END_REF], [START_REF] Czubak | On applications of wireless sensor networks[END_REF]. A WSN is composed of a large number of nodes, also called microsensors, which have communication, computing, and sensing functions [START_REF] Akyildiz | A survey on sensor networks[END_REF]. However, these tiny devices usually have limited amounts of memory, reduced processing capabilities, limited power resources, and low range communication capacities.

Challenges in WSNs gave rise to an active research field in machine learning, with kernel-based methods for regression and pattern recognition, including problems such as localization [START_REF] Nguyen | A kernel-based learning approach to ad hoc sensor network localization[END_REF], detection [START_REF] Nguyen | Nonparametric decentralized detection using kernel methods[END_REF], and regression [START_REF] Honeine | Functional estimation in hilbert space for distributed learning in wireless sensor networks[END_REF], [START_REF] Honeine | A decentralized approach for non-linear prediction of time series data in sensor networks[END_REF] in WSNs. The framework presented in [START_REF] Honeine | Distributed regression in sensor networks with a reduced-order kernel model[END_REF] allows the representation of the investigated field by a global model linking information from all sensors. These sensors are often densely and randomly deployed to insure a good coverage of the region under scrutiny. In order to overcome this problem, one may consider mobile sensors, which move in a manner that optimizes coverage [START_REF] Shorey | Mobile, Wireless, and Sensor Networks: Technology, Applications, and Future Directions[END_REF], [START_REF] Song | A mobile sensor network system for monitoring of unfriendly environments[END_REF], [START_REF] Mourad | Controlled mobility sensor networks for target tracking using ant colony optimization[END_REF]. In this paper, some sensors are assumed to be robots having a controlled mobility, and thus they move in such a way to collect more relevant information.

More specifically, this study tackles the problem of sensors mobility using a kernel-based regression in a wireless sensor network. For this purpose, we propose a kernel-based model and we update it using information from robot sensors. We show that the proposed framework allows us to derive an efficient scheme for sensor mobility, by minimizing the approximation error. It turns out that the resulting problem has tight connections with a well-known problem in machine learning, the so-called pre-image problem [START_REF] Honeine | A closed-form solution for the preimage problem in kernel-based machines[END_REF], [START_REF] Honeine | Preimage problem in kernel-based machine learning[END_REF]. We take advantage of recent developments in this area to derive two optimization schemes for the mobility problem.

The rest of this paper is organized as follows: In the next section, we derive the framework of kernel-based methods for learning in WSNs. In Section III, we present the optimization problem and describe several optimization schemes in Section IV. Section V provides experimentation results and discussions, whereas Section VI concludes the paper.

II. KERNEL-BASED METHODS IN WSNS

We aim to estimate a diffusion field using a WSN by constructing a model based on information collected from sensors. This paper considers the framework of kernel-based methods. Let N be the number of sensors deployed in a region X, where X ⊂ R 2 or X ⊂ R 3 , for a two-or threedimensional space. Let x i ∈ X be the position of sensor i, 1 ≤ i ≤ N , y i ∈ R be the measurement made by sensor i of the studied physical quantity such as a temperature measure or a gas concentration. The field is modeled using the information of the N sensors. For this purpose, their measurements are collected and the field is estimated by the real-valued function ψ(•) defined on X such that ψ(x j ) ≈ y j for all j.

In order to define the function ψ(•), the reproducing kernel formalism is considered [START_REF] Schölkopf | A Generalized Representer Theorem[END_REF]. Let H be a reproducing kernel Hilbert space constituted by functions from X to R. We denote by •, • H the associated scalar product, and by • H the corresponding norm. A (positive definite) kernel κ(•, •) denotes a real-valued function defined on X × X, and verifying

∀x j ∈ X, ∀ψ(•) ∈ H, ψ(x j ) = ψ(•), κ(•, x j ) H .
For some arbitrary reproducing kernel Hilbert space H given by the choice of the corresponding kernel κ(•, •), the optimization problem is given by minimizing the mean-squared-error between the model output ψ(x j ) and the desired output y j :

ψ(•) = arg min ψ∈H N j=1 |y j -ψ(x j )| 2 + ν ψ 2 H , (1)
where ν is a regularization parameter that controls the tradeoff between fitting the available data and the smoothness of the solution. Note that the appropriate choice of a neighborhood also allows such control, thereby reducing the contribution of ν in the optimization problem. Due to the Representer Theorem [START_REF] Schölkopf | A Generalized Representer Theorem[END_REF], the solution to this regularized optimization problem is given by a linear combination of the κ(x j , •), for j ∈ {1, 2, • • • , N }, namely:

ψ(•) = N j=1 α j κ(x j , •).
(

) 2
Since the model is linear with the parameters α j , the latter can be easily estimated as follows. Let α be the vector containing the α j for j ∈ {1, 2, • • • , N }. We can show that

α = (K + νI) -1 y, ( 3
)
where y is the vector containing y j and K is the Gram matrix whose entries are κ(

x j , x l ), for j, l ∈ {1, 2, • • • , N }.
In this paper, we use radial kernels of the form:

κ(x i , x j ) = g( x i -x j 2 ), (4)
where g(•) is a positive real function. Radial kernels are very natural for the studied problem, since they depend only on the distance and are translation invariant. Let g ′ (•) denote the first derivative of the function g(•) with respect to its argument, i.e., g ′ (z) = ∂g(z)/∂z. The most used radial kernel is the Gaussian kernel, where g(•) = exp -1 2σ 2 (•) , with σ the tunable bandwidth parameter. In this case, we have g ′ (z) = -1 2σ 2 g(z).

III. MOBILITY

The aim of the proposed method is to adapt iteratively the learning set in such a way to improve the regression model. Let A t = {(x 1,t , y 1,t ), (x 2,t , y 2,t ), • • • , (x N,t , y N,t )} be the learning set at iteration t, where A 0 is the set with the initial learning sensors, and let ψ t (•) be the model computed according to Section II using A t , namely

ψ t (•) = N j=1 α j,t κ(x j,t , •), (5)
where the coefficients α j,t are obtained using expression (3) with the set A t . The quadratic validation error of the model at any given (x + , y + ) is given by the squared difference between the desired output y + and the estimated model output ψ t (x + ):

ε t (x + ) = |y + -ψ t (x + )| 2 . ( 6
)
Optimizing jointly all the positions x j,t and the coefficients α j,t , for j = 1, 2, . . . , N , is an intractable problem. This is mainly due to the nonlinearity of ψ t (•) with respect to each of the samples x j,t , while it remains linear with respect to the coefficients. Moreover, the estimated coefficients α j,t can cancel the model error at any entry of the learning set A t .

To overcome these difficulties, we propose to use M robot sensors. These robot sensors are denoted by the set

B t = {(x * 1,t , y * 1,t ), (x * 2,t , y * 2,t ), • • • , (x * M,t , y * M,t )}, where x * i ∈ X and y * i ∈ R for i ∈ {1, 2, • • • , M }.
Only a single robot is moved at each iteration. Consider the mobility of the i-th robot at iteration t. Therefore, the entry

(x * i,t-1 , y * i,t-1 ) of the set B t-1 is updated into (x * i,t , y * i,t ), namely B t = B t-1 \{(x * i,t-1 , y * i,t-1 )} ∪ {(x * i,t , y * i,t )}. (7)
Several strategies can be adopted to choose the robot to move at each iteration: a cyclic selection, a randomization, or a selection by considering the largest prediction error. It is worth noting that the latter outperforms other selection criteria. It is therefore considered in our method.

The proposed optimization strategy operates in the following steps at each iteration t:

1) Estimate the model ψ t (•) using the learning subset A t , according to Section II, with

ψ t (•) = N j=1 α j,t κ(x j,t-1 , •).
2) Compute, for i = 1, 2, . . . , M , the quadratic prediction error for each robot (x * i,t-1 , y * i,t-1 ), according to (6) with

ε t (x * i,t-1 ) = |y * i,t-1 -ψ t (x * i,t-1 )| 2 .
3) Select, from all M indices, the index i whose robot is yielding the largest error ε t (x * i,t-1 ) and then move the robot from x * i,t-1 to x * i,t in order to reduce the prediction error, namely by minimizing the function

ε t (•), x * i,t = arg min x∈X ε t (x). (8)
All other robots remain fixed, i.e., x * j,t = x * j,t-1 for all j = 1, 2, . . . , M and j = i. 4) Update the set B t-1 to B t as given in [START_REF] Honeine | A decentralized approach for non-linear prediction of time series data in sensor networks[END_REF]. 5) Update the set A t-1 to A t by adding the previous information (x * i,t-1 , y * i,t-1 ) to the set A t-1 , and increasing by one the number N : A t = A t-1 ∪ {(x * i,t-1 , y * i,t-1 )}. The main difficulty is the optimization problem (8), i.e., the minimization with respect to x of the quadratic cost function

ε t (x) = y * i,t-1 - N j=1 α j,t κ(x j,t-1 , x) 2 . (9)
The exact solution corresponds to solving ∇ε i,t (x) = 0 for which H(x) is positive semi-definite, where ∇ε t (x) is the gradient of the cost function ε t (•) with respect to its argument and H(x) is its Hessian matrix applied on x. However, this is a hard problem, since it is a non-linear and nonconvex optimization problem, mainly due to the nature of the considered kernel. Having a first guess which is the initial position x * i,t-1 , one is able to optimize it within its vicinity using iterative optimization schemes. In this case, x is within the vicinity of x * i,t-1 where the desired output is assumed to be constant and equal to y * i,t-1 . The rest of this paper is devoted to the minimization process using such schemes.

IV. OPTIMIZATION SCHEMES

In this section, we derive two first-order and one secondorder iterative optimization schemes for adaptive sampling.

To this end, we propose three possible schemes: the fixedpoint method, the gradient descent method, and the Newton's method. While the proposed methods can be applied for any type of kernel, the formulation is given here for radial kernels, i.e., of the form (4). In this case, the gradient of the error function ( 9) is given by

∇ε t (x) = -4 y * i,t-1 - N j=1 α j,t κ(x j,t-1 , x) × N j=1 α j,t g ′ ( x -x j,t- 1
2 )(xx j,t-1 ).

In this expression, we assume that x is within the vicinity of x * i,t-1 where the desired output remains constant and equal to y * i,t-1 .

A. Fixed-point method

The minimum of the quadratic cost function ( 9) is obtained when its gradient goes to zero. Therefore, at the minimum ∇ε t (x) = 0, namely

N j=1 α j,t-1 g ′ ( x -x j,t-1 2 ) (x -x j,t-1 ) = 0.
This equation leads to the so-called fixed-point expression, with the following update rule of sample x * i,t-1 to the new sample x * i,t :

x * i,t = N j=1 α j,t-1 g ′ ( x * i,t-1 -x j,t-1 2 ) x j,t-1 N j=1 α j,t-1 g ′ ( x * i,t-1 -x j,t-1 2 ) . (10)

B. Gradient descent method

The gradient descent scheme adapts in the direction opposed to the gradient, with the update rule

x * i,t = x * i,t-1 -η t ∇ε t (x * i,t-1 ), (11)
where η t > 0 is the stepsize parameter, which depends on the iteration t. The use of the stepsize here allows controlling the convergence, as opposed to the fixed-point method, the price to pay being the appropriate choice of its value. For appropriate convergence, the stepsize η t should decrease at each iteration [START_REF] Robbins | A Stochastic Approximation Method[END_REF]. A common choice of the stepsize is η t = η 0 /t, where η 0 is a positive constant parameter. Another possibility is a "search then converge" approach [START_REF] Honeine | Online kernel principal component analysis: a reducedorder model[END_REF], with η t = η0 1+t/τ . In this case, the delay parameter τ determines the duration of the initial search phase, with η t ≃ η 0 when t ≪ τ , before a converge phase where η t decreases as η 0 /t when t ≫ τ . This is well suited to the hypothesis of the sensors' limited energy.

C. Newton's method

The step size η is taken based on the Hessian matrix H. Thus, we obtain:

x * i,t = x * i,t-1 -H -1 ∇ε t (x * i.t-1 ). ( 12
)
This method allows to follow the curvature of ε t (x * i,t-1 ) by considering its second partial derivatives. By this way, the step size is imposed by the curvature of the cost function. However, this calculation has a greater computational complexity compared to the gradient descent method.

V. EXPERIMENTS

In order to illustrate the results of this work, we consider the diffusion of a gas in a two-dimensional space X = [-0.5, 0.5] × [-0.5, 0.5]. This distribution is governed by the following differential equation:

∂G(x, t) ∂t -c ∇ x 2 G(x, t) = Q(x, t), (13)
where G(x, t) is the density of gas depending on the position and time, ∇ x 2 G(x, t) is the Laplace operator, Q(x, t) corresponds to the added quantity of gas, and c is the conductivity of the medium. A gas source placed at position S(0, 0) is activated. We consider N = 36 sensors deployed with a uniform distribution, and M = 5 mobile (robot) sensors randomly deployed in the region under scrutiny X. The volume of the sensors is assumed to be very small compared to the dimension of the space. The simulations, realized using Matlab, aim at determining the quantity of gas at any spot of the space. Once it has reached a certain equilibrium, the distribution of gas in the space is assumed to be invariant with respect to time. In the following, two settings are investigated.

In the first one, the sensors are assumed to be fixed, i.e., static, and a global model is estimated ψ s (subscript s for static).

In the second one, robot sensors are mobile, and the model denoted by ψ m (subscript m for mobile) is computed and robot sensors are able to move according to either update rules [START_REF] Song | A mobile sensor network system for monitoring of unfriendly environments[END_REF], [START_REF] Mourad | Controlled mobility sensor networks for target tracking using ant colony optimization[END_REF] or [START_REF] Honeine | A closed-form solution for the preimage problem in kernel-based machines[END_REF]. At each time step, the learning set is updated by moving a single robot sensor, selected from the 5 available mobile sensors by considering the largest predictive error. The global model is then estimated using the whole updated learning set. The relevance of the resulting model is evaluated with the root mean-squared error estimated using 16 test sensors. These test sensors are uniformly deployed in the region X under scrutiny and measure at these positions the corresponding gas density. Called test sensors, they have fixed positions and the corresponding information is not used in the learning process. In the following, we consider 100 iterations, i.e., t = 1, 2, . . . , 100. At each iteration, a single sensor is moved using a mobility scheme according to one of the methods given in Section IV. Moreover, since the mobility must be performed within the space X, if the calculation sends the sensor outside X, the sensor remains in its place. The Gaussian kernel is considered in all simulations.

In order to illustrate the proposed methods, we need to tune up some parameters. In kernel-based machine learning, the kernel parameter need to be adjusted. The bandwidth of the Gaussian kernel can be fine-tuned in the static setting, by using any cross-validation scheme. The bandwidth of the Gaussian kernel is set to σ = 0.5. The fixed-point optimization method and Newton's method do not require the tuning of any additional parameter, while the gradient descent method depends on the stepsize parameter η t , which is taken of the form η0 1+t/τ as recommended in Section IV-B. We set η 0 to 2 -12 and τ to 70 for the gradient descent optimization scheme.

Figure 1 shows the root-mean-squared errors at each time step for the static settings, i.e., without mobility, and the mobility governed by the three proposed optimization methods. The plot shows that the proposed methods are able to reduce significantly the validation error and thus improve the estimated model. Moreover, the Newton's method performs better than the gradient descent and the fixed-point one, since the error decreases continuously until getting to convergence.

Beyond the approximation error, one is often interested in the traveled distance, which is directly related to the energy consumption in a mobile WSN. Keeping limited node movements is crucial, because of the limited energy constraint of sensors. That is why we are interested in small paths for the sensors. Table I shows the average of the traveled distances of all sensors, as well as the largest distance. It is easy to see that the fixed-point method provides large variations, mainly due to the absence of a stepsize parameter to control the mobility. The Newton's method has, on average, a better efficiency. However, this method has a more important computational cost compared to the others. Therefore, we can conclude that, depending on the application case and the constraints imposed, one can choose one or the other of these methods.

VI. CONCLUSION

In this paper, we model a gas field measured by a network of mobile wireless sensors. We optimize this model by moving five robot sensors deployed randomly in the area to be monitored. Models built using both fixed sensors and mobile sensors are compared. In both cases, the objective is to define the global model that best matches the measured neighboring amount of gas. These models are determined by means of a non-linear learning process using Gaussian kernels. According to the realized simulations, we demonstrated a great deal of error reduction by the proposed methods, compared to the static case.

  Figure 1. Root mean squared error over the test sensors.

		0.46
		0.44
	Root mean square error on test sensor	0.32 0.34 0.36 0.38 0.4 0.42		without mobility mobility according to gradient mobility according to Newton mobility according to fixed-point
		0.3
		0	50	100	150	200
				Number of steps
				Table I
	AVERAGE DISTANCE AND MAXIMUM DISTANCE FOR EACH
				MOBILITY SCHEME.
	Mobility scheme: Fixed-point Gradient Newton
	Average distance:	0.1660	0.0664	0.0768
	Largest distance:	0.8299	0.1875	0.1854


Sandy Mahfouz

Paul Honeine

Farah Mourad-Chehade

Joumana Farah

Hichem Snoussi

Combining a physical model with a nonlinear fluctuation for signal propagation modeling in WSNs

Keywords: Distance estimation, kernel functions, multikernel learning, RSSI, semiparametric regression

In this paper, we propose a semiparametric regression model that relates the received signal strength indicators (RSSIs) to the distances separating stationary sensors and moving sensors in a wireless sensor network. This model combines the well-known log-distance theoretical propagation model with a nonlinear fluctuation term, estimated within the framework of kernel-based machines. This leads to a more robust propagation model. A fully comprehensive study of the choices of parameters is provided, and a comparison to state-of-the-art models using real and simulated data is given as well.

I. INTRODUCTION

Received signal strength indicator (RSSI), a standard feature in most radios, has shown attractive properties that could be exploited for localization in WSNs [START_REF] Lau | Enhanced RSSI-based real-time user location tracking system for indoor and outdoor environments[END_REF]. In fact, one can take advantage of the attenuation of the signal strength with the increase of the traveled distance, to estimate the distances between sensors, using only RSSI information, without the need of additional hardware. However, estimating the exact distances using RSSI turns out to be really challenging, since the measurements of signals' powers could be significantly altered by the presence of additive noise, multipath fading, shadowing, and other interferences. Therefore, it becomes more and more important to find a mathematical model which can accurately describe the relationship between the RSSI values and the distances.

Several models have been proposed in the literature to characterize the relationship between the RSSIs and the distances. A popular one is the Okumura-Hata model, also known as the log-distance propagation model [START_REF] Medeisis | On the use of the universal Okumura-Hata propagation prediction model in rural areas[END_REF], [START_REF] Patwari | Locating the nodes: cooperative localization in wireless sensor networks[END_REF], [START_REF] Zanella | Rss-based ranging by multichannel rss averaging[END_REF]; even though this model has many limitations, it is still widely used because of its simplicity. However, this model is basically for outdoors, since it predicts the signal strength without taking the surrounding environment into account, such as the walls and the floors. Therefore, it becomes inaccurate in cases where there is no line of sight between sensors. Different models have been proposed to overcome this problem, such as modified versions of the log-distance model [START_REF] Seidel | 914 mhz path loss prediction models for indoor wireless communications in multifloored buildings[END_REF], [START_REF] Sandeep | Wireless network visualization and indoor empirical propagation model for a campus wi-fi network[END_REF], in which the attenuations due to floors and walls are explicitly included. In other ones [START_REF] Wang | An indoors wireless positioning system based on wireless local area network infrastructure[END_REF], [START_REF] Yang | Indoor localization using improved rss-based lateration methods[END_REF], the authors determine a mathematical relationship between RSSIs and distances, without taking physical properties into consideration. To this end, they estimate an empirical model based on polynomial regression. Depending on the applications and the considered environment, several other models can also be found in the literature [START_REF] Neskovic | Modern approaches in modeling of mobile radio systems propagation environment[END_REF], [START_REF] Sarkar | A survey of various propagation models for mobile communication[END_REF], [START_REF] Andrade | IEEE 802.11 WLANs: A comparison on indoor coverage models[END_REF].

Lately, kernel machine methods have been greatly used in the context of function approximation and model fitting, since they simplify the specification of a nonparametric model, especially in the case of nonlinear data [START_REF] Mahfouz | Kernel-based localization using fingerprinting in wireless sensor networks[END_REF], [START_REF]Decentralized localization using fingerprinting and kernel methods in wireless sensor networks[END_REF], [START_REF] Mahfouz | Target tracking using machine learning and kalman filter in wireless sensor networks[END_REF]. In addition, kernel machine methods are being widely used for semiparametric regression, to characterize relationships between data in the medical domain as in [START_REF] Liu | Semiparametric regression of multidimensional genetic pathway data: Least-squares kernel machines and linear mixed models[END_REF], or in hyperspectral imaging as in [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF]. In this paper, we propose a semiparametric propagation model to define the relationship between the RSSIs and distances. The proposed model combines the well-known log-distance propagation model with a nonlinear fluctuation term, defined in a reproducing kernel Hilbert space. The nonlinear term added in the proposed model compensates for the missing term in the log-distance model, i.e., all environment's sources of noises, therefore allowing a better modeling of the RSSI/distance relationship. To our knowledge, such a model has not yet been proposed in the literature.

The rest of the paper is organized as follows. Section II gives a brief description of two well-known state-of-theart models, along with an estimation of their parameters. Section III describes the proposed semiparametric model. Section IV provides a comparison of the proposed model to the two state-of-the-art models using real and simulated data, while Section V concludes the paper.

II. STATE-OF-THE-ART MODELS

As previously explained, many signal propagation models, that characterize the relationship between the RSSIs and the distances, are proposed in the literature. Some models take into account the physical characteristics of the environment, such as the path-loss in the well-known theoretical log-distance model [START_REF] Medeisis | On the use of the universal Okumura-Hata propagation prediction model in rural areas[END_REF], [START_REF] Patwari | Locating the nodes: cooperative localization in wireless sensor networks[END_REF]. Others avoid the physical properties by finding a mathematical relationship based on polynomial regression [START_REF] Wang | An indoors wireless positioning system based on wireless local area network infrastructure[END_REF], [START_REF] Yang | Indoor localization using improved rss-based lateration methods[END_REF].

A training phase is necessary for both types of models in order to determine their parameters. To this end, N reference positions, denoted by p ℓ , ℓ ∈ {1, . . . , N }, are generated uniformly or randomly in the studied environment. A stationary sensor with known position s continuously broadcasts signals in the network at a fixed initial power, and a sensor is placed consecutively at the reference positions to detect the broadcasted signals and measure their RSSIs. Let ρ ℓ be the power received from the stationary sensor at position s by the sensor at position p ℓ . Then, the distances d ℓ between each reference position p ℓ and the stationary sensor are computed. In this way, a training set of N pairs (ρ ℓ , d ℓ ) is obtained, that can be used for the estimation of the model's parameters.

After determining the model's parameters, one can estimate the distance between the stationary sensor at position s and a sensor at an unknown position x i using the estimated model and the RSSI value ρ i measured by the sensor.

A brief description of the log-distance model and the estimation of its parameters is provided in the first subsection. In the second subsection, a description of the mathematical model based on polynomial regression is provided, along with its parameters.

A. Log-distance propagation model

Both theoretical and measurements-based propagation models indicate that the received signal strength decreases logarithmically with distance. Therefore, the log-distance propagation model, that characterizes the relationship between RSSIs and distances, is expressed as follows:

ρ i = ρ 0 -10 n P log 10 d i d 0 , (1)
where ρ i is the power received from the stationary sensor at position s by the sensor at position x i , ρ 0 is the power at the reference distance d 0 , d i = sx i is the Euclidian distance between the position x i of the considered node and the position s of a stationary sensor, and n P is the path-loss exponent. The value of n P depends on the specific propagation environment, i.e. type of construction material, architecture, location within building, dimensions of the coverage area, etc. For low values of n P , the signal loss towards distance is low. Developing equation (1) leads to the following:

log 10 d i = ρ 0 10 n P + log 10 d 0 - ρ i 10 n P . (2)
Now let a = -1 10 nP and b = ρ0 10 nP + log 10 d 0 . Then equation ( 2) can be written as follows:

log 10 d i = a ρ i + b. (3)
The parameters a and b depend on the characteristics of the environment and can be estimated from the N reference measurements using the Least-Squares (LS) method. They must be chosen in a way to minimize the sum of squared residuals, i.e. the mean squared error on the training set, given by:

1 N N ℓ=1 ε 2 ℓ . (4)
Here, ε ℓ = log 10 d ℓa ρ ℓb, and ℓ ∈ {1, . . . , N }. By taking the partial derivatives of ( 4) with respect to a and b and setting them to zero, one gets the following in matrix form:

1 ⊤ D D = 1 ⊤ ρ ρ N 1 a b , (5)
where 1 is the N × 1 vector of ones, D is the N × 1 vector whose ℓ-th entry is equal to log 10 d ℓ , and ρ is the N × 1 vector whose ℓ-th entry is equal to ρ ℓ , with ℓ ∈ {1, . . . , N }.

Having computed the values of a and b using ( 5), one can estimate the distance di between the considered stationary sensor at position s and a sensor at an unknown position x i in the network using the following equation: di = 10 a ρi+b , where ρ i is the power received from the stationary sensor at position s by the sensor at position x i .

B. Polynomial model

Having described the log-distance propagation model and its parameters, we introduce in this subsection the polynomial model. The objective of using a polynomial regression is to determine a mathematical relation between RSSIs and distances, without having to take physical properties into consideration. In such case, the signal propagation model is then the n-th degree polynomial given by the following:

d i = a 0 + a 1 ρ i + a 2 ρ 2 i + • • • + a n ρ n i , (6)
where a j , j ∈ {0, . . . , n}, are the polynomial's coefficients to be determined. The parameters a j should be chosen in a way to fit the n-th degree polynomial through the training set. The LS method is used here as well; therefore, we need to find the parameters a j that allow us to minimize the mean squared error on the training set as given in (4), with

ε ℓ = d ℓ - a 0 -a 1 ρ ℓ -a 2 ρ 2 ℓ + • • • -a n ρ n ℓ
, and ℓ ∈ {1, . . . , N }. The partial derivatives of (4) with respect to a 0 , a 1 , . . . , a n are then set to zero and rearranged in order to obtain the following matrix form of the solution:

     d 1 d 2 . . . d N      =      ρ 0 1 ρ 0 2 . . . ρ 0 N . . . . . . . . . . . . ρ n 1 ρ n 2 . . . ρ n N           a 0 a 1 . . . a n      . ( 7
)
Having estimated the parameters a j , j ∈ {1, . . . , n} using [START_REF] Wang | An indoors wireless positioning system based on wireless local area network infrastructure[END_REF], one can then estimate the distance di separating the stationary sensor at position s from the sensor at unknown position x i using (6).

Remark 1. It is interesting to notice that replacing log 10 d i by d i in (3) leads to a particular case of [START_REF] Sandeep | Wireless network visualization and indoor empirical propagation model for a campus wi-fi network[END_REF], where the parameters a j , j ∈ {2, . . . , n}, are null. In fact, as already mentioned, it is proven that the received signal strength decreases logarithmically with distance; in other words, the distance in logarithm scale is a linear function of the received signal strength as shown in [START_REF] Patwari | Locating the nodes: cooperative localization in wireless sensor networks[END_REF].

III. PROPOSED SEMIPARAMETRIC MODEL As we already stated, the received signal strength decreases logarithmically with distance. However, the logdistance signal propagation model of (1) alone is insufficient to characterize the RSSI/distance relationship, since it assumes that the latter is fully parametric, neglecting then all other factors in the environment, such as physical obstacles, multipath propagation, additive noises, etc. In this section, we propose a new model to characterize the relationship between the RSSIs and the distances. Our model is a semi-parametric one that combines the theoretical log-distance signal propagation model with a nonlinear fluctuation term, that represents a combination of all unknown factors affecting the RSSI measures. This term should provide the log-distance physical model with flexibility, resulting in a more accurate model.

Let ψ(•) denote the model that associates to each RSSI ρ i the logarithm of the distance log 10 (d i ). According to the definition given in the beginning of this section, ψ(•) can be decomposed into twos terms: a linear term and a nonlinear fluctuation term, as follows:

ψ(•) = ψ lin (•) + ψ nlin (•) (8)
To determine ψ(•), we consider the log-distance model of (1), and add a noise term ϕ. Then the relationship between the RSSIs and the distances is given by:

ρ i = ρ 0 -10 n P log 10 d i d 0 + ϕ. (9)
This model can be written as follows:

log 10 d i = ρ 0 10 n P + log 10 d 0 - ρ i 10 n P + ϕ 10 n P . (10)
One can see that [START_REF] Sarkar | A survey of various propagation models for mobile communication[END_REF] is a combination of a linear model in terms of ρ i and a nonlinear model. Thus, one can conclude the following:

ψ lin (ρ i ) = α ρ i + β, ψ nlin (ρ i ) = ϕ 10 nP , (11)
where α = -1 10 nP and β = ρ0 10 nP + log 10 d 0 . As for the nonlinear term ψ nlin (•), we assume that it lies in a reproducing kernel Hilbert space denoted by H nlin , and generated by a positive definite kernel function κ nlin (•, •). From the representer's theorem [START_REF] Kimeldorf | Some results on tchebycheffian spline functions[END_REF], [START_REF] Schölkopf | A generalized representer theorem[END_REF], ψ nlin (•) can be written as a linear combination of kernels:

ψ nlin (•) = N ℓ=1 γ ℓ κ nlin (ρ ℓ , •), (12)
where κ nlin : IR × IR → IR, and γ ℓ , ℓ ∈ {1, . . . , N }, are parameters to be estimated. According to [START_REF] Haussler | Convolution kernels on discrete structures[END_REF], it can be shown that the direct sum of H lin and H nlin of the RKHS of kernels κ lin (ρ u , ρ v ) = ρ ⊤ u ρ v and κ nlin is also a RKHS with the following kernel function:

κ(ρ u , ρ v ) = ρ ⊤ u ρ v + κ nlin (ρ u , ρ v ), = κ lin (ρ u , ρ v ) + κ nlin (ρ u , ρ v ).
= c + ρ ⊤ u ρv q Exponential κ nlin (ρu, ρv) = exp 1 σ ρ ⊤ u ρv Gaussian κ nlin (ρu, ρv) = exp -1 2σ 2 ρu -ρv 2
Table I shows some of the most commonly used kernel functions κ nlin . Now the N + 2 parameters, α, β, and γ = (γ 1 . . . γ N ) ⊤ , to be estimated can be found using the training set and the Least-Squares (LS) method. They must be chosen in a way to minimize the mean squared error on the training set, given by:

1 N N ℓ=1 ε 2 ℓ + η ψ nlin 2 H nlin . (13)
Here, ε ℓ = log 10 d ℓα ρ iβ -N j=1 γ j κ(ρ j , ρ ℓ ), ℓ ∈ {1, . . . , N }, and the quantity η is a regularization parameter that controls the tradeoff between the training error and the complexity of the solution.

By multiplying (13) by N for convenience and writing it in matrix form, one gets the following:

D ⊤ D + α 2 ρ ⊤ ρ + N β 2 + γ ⊤ K ⊤ nlin K nlin γ -2αD ⊤ ρ -2β1 ⊤ D -2D ⊤ K nlin γ + 2αβ1 ⊤ ρ + 2αρ ⊤ K nlin γ + 2β1 ⊤ γ + η N γ ⊤ K nlin γ, (14)
where K nlin is the N ×N matrix whose (u, v)-th entry is κ nlin (ρ u , ρ v ), for u, v ∈ {1, ..., N }.

Then, the partial derivatives in matrix form of ( 14) with respect to α, β and γ are given by the following:

     1 2 ∂ ∂α = ρ ⊤ ρ α + 1 ⊤ ρ β + ρ ⊤ K nlin γ -D ⊤ ρ, 1 2 ∂ ∂β = 1 ⊤ ρ α + N β + 1 ⊤ K nlin γ -1 ⊤ D, 1 2 ∂ ∂γ = ρ α + 1 β + (K nlin + ηN I)γ -D, ( 15
)
where I is the N ×N identity matrix.

Setting the derivatives in [START_REF] Liu | Semiparametric regression of multidimensional genetic pathway data: Least-squares kernel machines and linear mixed models[END_REF] to zero will lead to a linear system having the form B = AX, where

X =   α β γ   , B =   1 ⊤ D D D ⊤ ρ   , A =   1 ⊤ ρ ρ ρ ⊤ ρ N 1 1 ⊤ ρ 1 ⊤ K nlin K nlin + ηN I ρ ⊤ K nlin   ( 16
)
The solution is then given by:

X = (A ⊤ A) -1 A ⊤ B. (17)
After computing the model's parameters using [START_REF] Kimeldorf | Some results on tchebycheffian spline functions[END_REF], one can find the logarithm of any distance separating a stationary sensor and a moving sensor in the network using only the RSSI ρ i , as follows:

log 10 d i = ψ(ρ i ) = αρ i + β + N ℓ=1 γ ℓ κ(ρ ℓ , ρ i ). ( 18
)
The distance is then given by d i = 10 ψ(ρi) .

Remark 2. It is interesting to see that the solution obtained in ( 5) is actually included in the final solution [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF] of the proposed semiparametric model. Indeed, the main difference between the proposed model and the log-distance model lies in the nonlinear fluctuation term, whose inclusion induces the introduction of γ, thus the incorporation of N new equations to the system.

Remark 3. In this paper, we use the gaussian kernel of Table I. Nevertheless, let us consider the polynomial kernel and study the effect of such choice on our model. From [START_REF] Schölkopf | A generalized representer theorem[END_REF] and Table I, one can write the following:

log 10 d i = αρ i + β + N ℓ=1 γ ℓ c + ρ ⊤ ℓ ρ i q , ( 19
)
where c > 0, and q ∈ IN + . For q = 1, the terms of ( 19) can be rearranged as follows:

log 10 d i = β + c N ℓ=1 γ ℓ + (α + N ℓ=1 γ ℓ ρ ℓ )ρ i .
This equation takes the same form as [START_REF] Patwari | Locating the nodes: cooperative localization in wireless sensor networks[END_REF]. Therefore, when choosing the polynomial kernel with q = 1, the obtained model is equivalent to the log-distance propagation model. Now let q = 2; equation ( 19) can then be written as follows:

log 10 d i = β+c 2 N ℓ=1 γ ℓ +(α+2c N ℓ=1 γ ℓ ρ ℓ )ρ i +( N ℓ=1 γ ℓ ρ 2 ℓ )ρ 2 i .
This equation is of the following form:

log 10 d i = a 0 + a 1 ρ i + a 2 ρ 2 i .
For larger values of q, one will eventually be able to rearrange the terms and write the following:

log 10 d i = a 0 + a 1 ρ i + a 2 ρ 2 i + • • • + a q ρ q i , (20)
where a j , j ∈ {0, . . . , q}, are coefficients that are functions of α, β, c, and γ ℓ , ℓ ∈ {1 . . . N }. Finally, one can conclude that choosing the polynomial kernel in our model leads to [START_REF] Zanca | RSSI measurements in indoor wireless sensor networks[END_REF], which is the same model described in Remark 1.

IV. EXPERIMENTAL RESULTS

We now propose to evaluate the accuracy of the proposed model, when used for distance estimation, in the case of real data and simulated data. In the first subsection, a set of collected measurements available from [START_REF] Zanca | RSSI measurements in indoor wireless sensor networks[END_REF] is used for the evaluation, and the results are then compared to ones obtained using the two state-of-the-art models of Section II. Then, in the second subsection, the model is evaluated on simulated data; the results are also compared to ones obtained using the models of Section II.

A. Evaluation of the model on real data

In this subsection, we use the set of collected measurements available from [START_REF] Zanca | RSSI measurements in indoor wireless sensor networks[END_REF] for the evaluation of the proposed model. The measurements are carried out in a room of approximately 10m × 10m, where 48 uniformly distributed EyesIFX sensor nodes are deployed. Furniture and people in the room cause multi-path interferences affecting the collected RSSI values. We consider that there are 4 fixed stationary sensors at known positions, and 44 other sensors with known positions for the training and test phases. The left plot of Fig. 1 shows the topology of the testbed.

It is important to note that the average values over time of the RSSIs are used in this section. In fact, the RSSIs vary a lot with respect to time and movements, as one can see in the right plot of Fig. 1. These variations are known as short-term or multi-path fading. On the other hand, the local average of the signal varies slowly. These slow fluctuations depend mostly on environmental characteristics, and they are known as long-term fading. Therefore, it is more suitable to use the average values of the RSSIs than to use all the collected values [START_REF] Neskovic | Modern approaches in modeling of mobile radio systems propagation environment[END_REF].

Based on the described scenario, one can see that there are 4 signal propagation models to be determined, i.e., one model per stationary sensor. Each model has a different set of training data to be used, and different parameters that need to be found. In fact, the RSSIs of the signals exchanged between each stationary sensor and the 44 other sensors are used along with the distances separating this stationary sensor from the other sensors. This information is then used in the training phase as described in III to compute the model's parameters. We use the Gaussian kernel from Table I, whose parameters are chosen in a way to minimize the error on the training set. The value of this error for the 4 computed models is given in Table II, along with the error on the training set for the physical log-distance propagation model used in [START_REF] Medeisis | On the use of the universal Okumura-Hata propagation prediction model in rural areas[END_REF], [START_REF] Patwari | Locating the nodes: cooperative localization in wireless sensor networks[END_REF], [START_REF] Zanella | Rss-based ranging by multichannel rss averaging[END_REF] and for the polynomial model used in [START_REF] Wang | An indoors wireless positioning system based on wireless local area network infrastructure[END_REF], [START_REF] Yang | Indoor localization using improved rss-based lateration methods[END_REF]. We also use the leave-one-out (LOO) technique in order to evaluate the performance of the proposed model in the case of data that are not part of the training set. The LOO technique involves using a single observation from the collected data as the validation data, and the remaining observations as the training data. This is repeated 44 times, such that each observation is used once as the validation data. This technique is interesting because it allows us to compare our proposed model to the log-distance model and the polynomial model, even though the set of collected data is not really large. Finally, the mean estimation error averaged over the 44 simulations is stored in Table II for the 4 computed models, in the case of the proposed model, as well as the two already-described ones. One can see from Table II that the proposed model yields better results than the state-of-the-art models, when comparing the mean estimation error.

B. Evaluation of the model on simulated data

In this subsection, we evaluate the proposed model using simulated data. To this end, we consider the average walls model described in [START_REF] Andrade | IEEE 802.11 WLANs: A comparison on indoor coverage models[END_REF] to generate the RSSI measures. This model is a modified version of the log-distance model that explicitly takes into account the attenuations due to walls. The received signal strength indicator is then given by the following:

ρ i = ρ 0 -10 n P log 10 d i d 0 -N wi L wi + ǫ i , (21)
where ρ i (in dB) is the power received from the stationary sensor at position s by the sensor at position x i , ρ 0 is the power at the reference distance d 0 set to 1dB, while the path-loss exponent n P is set to 4. As for the quantities L wi and N wi , they denote respectively the loss due to walls and the number of penetrated walls. The quantity L wi is taken equal to 6.9dB, since we consider the case of heavy thick walls [START_REF] Andrade | IEEE 802.11 WLANs: A comparison on indoor coverage models[END_REF]. Finally, the quantity ǫ i is a zero mean additive white noise, whose standard deviation is taken equal to 1% of the standard deviation of the RSSIs. Now consider the 25m × 5m area given in Fig. 2, where 2 fixed stationary sensors and 45 known positions for the training phase are considered. This figure shows that there are 5 rooms in the area; therefore, the signal penetrates a maximum of 4 walls during its propagation, i.e., N wi ∈ {0 . . . 4}. As for the test phase, 100 positions are randomly generated in the studied area, and their RSSIs are obtained using (21). Fig. 3 shows the RSSIs measured at the 45 training positions along with the distances separating these positions from the stationary sensor 1 of Fig. 2. One can see that the signal decreases with the traveled distances, and that the loss gets higher every time the signal penetrates a wall.

We now have to determine the two signal propagation models' parameters. As we previously explained, the RSSIs of the signals and the distances are used to estimate these parameters. The Gaussian kernel from Table I is used here as well. Table III shows the errors on the distances that are estimated using the proposed model, the physical logdistance propagation model, and the polynomial model. Compared to the estimation error obtained using the stateof-the-art models, the estimation error is reduced by half or more when using the proposed model, proving that the latter outperforms the other models in terms of accuracy.

V. CONCLUSION

In this paper, we proposed a semiparametric regression model that combines the well-known log-distance propagation model with a nonlinear fluctuation term, estimated within the framework of kernel-based machines. Evaluation on real and simulated data showed that the proposed model outperforms state-of-the-art models in terms of accuracy. Future works will handle improvements of the quality of the RSSIs; for instance, a filtering process can be considered to reduce the noise on the RSSIs in the case of real data before finding the signal propagation model.

Fig. 1 :

 1 Fig. 1: Left plot: Topology of the testbed, where ∆ represents the stationary sensors and + represents the training positions. Right plot: RSSIs measured at the training positions as a function of the distances separating these positions from the stationary sensor 1.

Fig. 2 :

 2 Fig. 2: Topology of the simulated area, where ∆ represents the stationary sensors and + represents the training positions.

Fig. 3 :

 3 Fig. 3: RSSIs measured at the training positions as a function of the distances separating these positions from the stationary sensor 1 of Fig. 2 -Simulated data.

TABLE I :

 I Some commonly used reproducing kernels, with parameters c, σ > 0, and q ∈ IN + .

	Type	General form
	Polynomial	κ nlin (ρu, ρv)

TABLE II :

 II Comparison between models (error in meters) for real data

			Physical	Polynomial	Polynomial	Polynomial	Proposed
			model	model (n = 2) model (n = 3) model (n = 4)	model
	Model 1	Training error	1.14	1.13	1.09	1.10	0.44
		LOO error	1.19	1.24	1.21	1.76	1.07
	Model 2	Training error	1.55	1.41	1.40	1.40	1.42
		LOO error	1.62	1.48	1.51	1.54	1.51
	Model 3	Training error	1.47	1.44	1.44	1.41	1.32
		LOO error	1.52	1.55	1.88	2.53	1.46
	Model 4	Training error	1.66	1.65	1.65	1.62	1.47
		LOO error	1.72	1.78	1.78	1.83	1.71
	Mean error	Training error	1.46	1.41	1.39	1.38	1.16
		LOO error	1.51	1.51	1.60	1.91	1.43

TABLE III :

 III Comparison between models (error in meters) for simulated data

			Physical	Polynomial	Polynomial	Polynomial	Proposed
			model	model (n = 2) model (n = 3) model (n = 4)	model
	Model 1	Training error	1.36	0.55	0.51	0.52	0.21
		Test error	1.25	0.60	0.57	0.57	0.23
	Model 2	Training error	1.13	0.57	0.55	0.55	0.17
		Test error	1.17	0.58	0.57	0.58	0.22
	Mean error	Training error	1.25	0.56	0.53	0.54	0.19
		Test error	1.21	0.59	0.57	0.58	0.23


Jie Chen
email: jie.chen@unice.fr

Paul Honeine
email: paul.honeine@utt.fr

ESTIMATING ABUNDANCE FRACTIONS OF MATERIALS IN HYPERSPECTRAL IMAGES BY FITTING A POST-NONLINEAR MIXING MODEL

Keywords: Nonlinear unmixing, post-nonlinear mixing model, kernel methods, hyperspectral data processing

Within the area of hyperspectral data processing, nonlinear unmixing techniques have emerged as promising alternatives for overcoming the limitations of linear methods. In this paper, we consider the class of post-nonlinear mixing models of the partially linear form. More precisely, these composite models consist of a linear mixing part and a nonlinear fluctuation term defined in a reproducing kernel Hilbert space, both terms being parameterized by the endmember spectral signatures and their respective abundances. These models consider that the reproducing kernel may also depend advantageously on the fractional abundances. An iterative algorithm is then derived to jointly estimate the fractional abundances and to infer the nonlinear functional term.

INTRODUCTION

There has been tremendous effort in the last decade to solve the linear unmixing problem, in the supervised, semisupervised and unsupervised learning setting, motivated by its simplicity and acceptable approximation error in many scenarios. However, in the very recent years, researchers have started developing nonlinear unmixing techniques, as nonlinear models may provide more accurate estimations of mixtures, in particular when there are significant physical interactions of photons scattered by multiple materials [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: geometrical, statistical, and sparse regression-based approaches[END_REF].

Several algorithms have been proposed to address the nonlinear unmixing problem. For instance, bilinear models were considered to handle complex scenarios such as multilayered scenes [START_REF] Raksuntorn | Nonlinear spectral mixture analysis for hyperspectral imagery in an unknown environment[END_REF][START_REF] Halimi | Nonlinear unmixing of hyperspectral images using a generalized bilinear model[END_REF], by introducing additional interaction terms in the linear model. An unmixing algorithm based on a manifold learning process was investigated in [START_REF] Bachmann | Exploiting manifold geometry in hyperspectral imagery[END_REF], under the assumption that hyperspectral data may be embedded into a low-dimensional manifold. Kernel methods aim to avoid high computational complexity by using more simple physics-inspired models. Some kernels were designed to be sufficiently flexible for matching various nonlinearity de-grees [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF]. Finally, artificial neural networks were considered to conduct unmixing with training data [START_REF] Plaza | Nonlinear neural network mixture models for fractional abundance estimation in AVIRIS hyperspectral images[END_REF].

More recently, it has been argued that both linear and nonlinear mixtures can be present in a single pixel [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF][START_REF] Chen | A novel kernel-based nonlinear unmixing scheme of hyperspectral images[END_REF][START_REF] Close | Estimating the percentage of linear and nonlinear mixing at the subpixel level in hyperspectral imaging[END_REF], meaning that partially linear models may be advantageously used. In [START_REF] Chen | A novel kernel-based nonlinear unmixing scheme of hyperspectral images[END_REF], we proposed a kernel-based nonlinear unmixing technique based on the partially linear assumption. In [START_REF] Close | Estimating the percentage of linear and nonlinear mixing at the subpixel level in hyperspectral imaging[END_REF], the authors combined the linear mixture model with intimate mixtures. In [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF], we investigated this problem under the framework of functional learning theory, where the tradeoff between the linear and the nonlinear components is automatically adjusted by solving a convex multi-kernel learning problem. However, in [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF][START_REF] Chen | A novel kernel-based nonlinear unmixing scheme of hyperspectral images[END_REF], in order to make the problem convex and tractable, the nonlinear component of the model was assumed to be only a function of the material spectral signatures, and thus to be independent of their abundances.

In this paper, we consider a more general problem where the nonlinear mixing term also depends on the fractional abundances. This means that the reproducing kernel depends on the abundances, and is directly learned from data. A specific kernel is proposed, and an efficient iterative algorithm based on the augmented Lagrangian principle is derived to jointly estimate the abundances and to infer the nonlinear functional term.

NOTATIONS AND THE MIXTURE MODEL

Notations

Suppose that in a scene, there exist R significant endmembers with spectral signature m i ∈ R L , for i = 1, . . . , R, where L denotes the number of spectral bands. Let r ∈ R L be an observed pixel, and α ∈ R R be the associated vector of abundances. Let us denote by M = [m 1 , . . . , m R ] ∈ R L×R the matrix of the endmember spectra. For the sake of convenience, the -th row of M will be denoted by m λ ∈ R L , that is, the column vector m λ contains the R endmember signatures at the -th wavelength band.

Post-nonlinear mixture model

In order to estimate the mixing ratios of the endmembers, we modeled in [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF] the nonlinear mixing mechanism ψ by a linear trend parameterized by the abundance vector α, combined with a nonlinear fluctuation term represented by ψ nlin , namely

ψ(m λ , α) = ψ lin (m λ , α) + ψ nlin (m λ ) with ψ lin (m λ , α) = α m λ
where the nonlinear fluctuation ψ nlin (m λ ) is independent of the abundance vector α. It can be expected that incorporating the abundance into the nonlinear function should be advantageous for a more accurate modeling. In this paper, we consider the following post-nonlinear model for

ψ nlin ψ(m λ , α) = ψ lin (m λ , α) + ψ nlin (α m λ ) with ψ lin (m λ , α) = α m λ . (1)
This model can characterize many different nonlinearities. For example, with the parabolic transform ψ nlin (x) = x 2 , model ( 1) is a polynomial post-nonlinear mixture model that mimics second-order interactions between material signatures, including cross-terms of the bilinear mixture model discussed in [START_REF] Halimi | Nonlinear unmixing of hyperspectral images using a generalized bilinear model[END_REF].

KERNEL-BASED NONLINEAR UNMIXING WITH THE POST-NONLINEAR MODEL

Problem formulation

Considering the problem (1), we intend to estimate the abundance vector α and to infer the nonlinear function ψ nlin by solving the following functional optimization problem1

ψ * = arg min ψ 1 2 ψ lin 2 Hlin + ψ nlin 2 Hnlin + 1 2µ L =1 e 2
where ψ = ψ lin + ψ nlin with

ψ lin (m λ ) = α m λ subject to e = r -(α m λ + ψ nlin (α m λ )) α 0 and 1 α = 1. (2)
The functional ψ nlin is restricted to be an element of the reproducing kernel Hilbert space H nlin , with the reproducing property

ψ nlin (x) = ψ nlin , κ nlin (•, x) ,
where κ nlin is the reproducing kernel associated to the Hilbert space H nlin . It should be noticed that solving problem ( 2) is challenging because it is no longer convex, and the dual variables cannot be expressed in closed-form as in [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF]. Fortunately, we can take advantage of its dual form for seeking a local optimum in an efficient way, via a variable splitting technique.

Let us introduce a new variable u that replaces α in the nonlinear fluctuation ψ nlin . An equality constraint that connects u and α must also be introduced to ensure the equivalence between the transformed problem and the original one. The problem ( 2) is now formulated as Thanks to this strategy, the derivatives with respect to each of the two variables α and u can now be evaluated separately.

α * = argmin
On the one hand, the variable α is estimated using an algorithm similar to the one derived in [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF]. On the other hand, the variable u is updated by using a gradient descent approach.

Problem solving with an ADMM-like algorithm

The augmented Lagrange function of the problem (3) is

J o a (ψ, α, u, e, ν) = J(ψ, α, u, e) + ν (u -α) + 1 2ζ u -α 2
where ζ > 0 is a positive parameter and ν ∈ R R is the dual variable associated to the constraint u = α. For simplicity, the function J o a can also be written in an equivalent scaled form J a with d the scaled dual variable [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF] J a (ψ, α, u, e, d) = J(ψ, α, u, e)

+ 1 2ζ u -α -d 2
Minimizing the above cost function can be performed by iterating the following three steps:

• Step 1: Update the variables α, ψ and e by α (k+1) , ψ (k+1) , e (k+1) = argmin α,ψ,e J a (ψ, α, u (k) , e, d (k) )

• Step 2: Update the variable u via a gradient descent

u (k+1) = u (k) -η ∇ u J a (ψ (k+1) , α (k+1) , u, e (k+1) , d (k) )
with η the step size.

• Step 3: Update the scaled dual variable d by

d (k+1) = d (k) + (u (k+1) -α (k+1) ).
This three-step iterative strategy corresponds to the ADMM algorithm in the convex case, which is known to converge to the global optimum. Unfortunately, the problem (3) is not convex, and the algorithm may lead to a local minimum depending on the initial values for the variables [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF]. However, the algorithm can take advantage of two properties of the postnonlinear mixture model. Firstly, the non-negativity and sumto-one constraints over α severely limit the space of acceptable solutions. Secondly, it is possible to initialize the algorithm with an appropriate solution obtained with conventional algorithms, e.g., a linear unmixing algorithm. During the extensive numerical simulations we carried out, we observed that the proposed algorithm always converged to reasonable solutions and usually provided better estimation results than state-of-the-art algorithms.

After presenting the principle of the algorithm, we shall now provide some extra details on Steps 1 and 2 calculation. In Step 1, variables α, ψ nlin and e can be updated via the Lagrange duality formulation, by fixing the other variables. By introducing the Lagrange multipliers β , γ and λ, where the superscript (k) has been omitted for simplicity of notation, the Lagrange function of the optimization problem described in Step 1 writes

G a = 1 2 α 2 + ψ nlin 2 Hnlin + 1 2µ L =1 e 2 - L =1 β (e -r + ψ nlin (u (k) m λ )) - R r=1 γ r α r + λ(1 α -1) + 1 2ζ u (k) -α -d (k) 2 . ( 4
)
The conditions for optimality of G a with respect to the primal variables are given by

     α (k+1) = ζ ζ+1 (M β * + γ * -λ * 1 + 1 ζ (u (k) -d (k) )) ψ (k+1) nlin = L =1 κ nlin (•, u (k) m λ ) β * e (k+1) = µ β * .
By substituting these conditions into equation ( 4), we get the dual problem (5) (see next page), where K nlin is the Gram matrix with ( , p)-th entry

[K nlin ] p = κ nlin (α m λ , α m λp ).
This allows us to estimate the optimum dual variables and thereby to obtain the optimal values of the primal variables.

In Step 2, due to the highly nonlinear kernel function, we suggest to update u using a gradient descent step

∇ u J a = - 1 2 ∂β * K nlin β * ∂u + 1 ζ (u -α (k+1) -d (k) ). (6)

Proposed kernel and gradient update

We propose the following kernel function for the model ( 2)

κ nlin (α m λ , α m λp ) = exp(α m λ m λp α).
It is easy to see that this kernel corresponds to an infiniteorder polynomial function of the linear mixture. This leads to the following expression

∂β K nlin β ∂u = M diag{β}K nlin diag{β}M u,
which allows to complete the evaluation of the gradient in (6).

EXPERIMENTAL RESULTS

Experiments on synthetic data

We first report experimental results on synthetic images. Three and five material spectra with 420 bands from ENVI software library were selected respectively to generate two scenes. For each scene, 2500 pixels were generated with different mixture models to evaluate and compare the performance of several algorithms. These models were the linear model, the bilinear mixture model, and the post-nonlinear mixing model (PNMM) defined by r = (M α) 0.7 + n. The abundance vectors α n , with n = 1, . . . , 2500, were uniformly generated in the simplex defined by non-negative and sum-to-one constraints. Finally, the images were corrupted with an additive white Gaussian noise n with SNR of 30 dB. Before comparing the unmixing performance of different methods, we illustrate the convergence of our algorithm when applied to a pixel mixed by the post-nonlinear model. Typical evolutions of the abundance vector α and its counterpart u are shown in Figure 1. Our approach was compared with the fully constrained least square method (FCLS) [START_REF] Heinz | Fully constrained least squares linear mixture analysis for material quantification in hyperspectral imagery[END_REF], the extended endmember matrix method (ExtM) [START_REF] Raksuntorn | Nonlinear spectral mixture analysis for hyperspectral imagery in an unknown environment[END_REF], and our previously proposed K-Hype method [START_REF] Chen | Nonlinear unmixing of hyperspectral data based on a linear-mixture/nonlinear-fluctuation model[END_REF]. The root mean square error of the estimated abundances was used to compare these algorithms. The comparative results are reported in Table 1.

Experiments on real data

In order to overcome the absence of ground truth reference for evaluating unmixing algorithms, unmixing was performed with various algorithms to estimate abundances in each pixel. The abundance vectors were used as features with an SVM classifier, based on a Gaussian kernel. Note that other classification algorithms such as [START_REF] Abdallah | An improved training algorithm for nonlinear kernel discriminants[END_REF] could have been used. This

β,γ,λ Ga(β, γ, λ) = - ρ 2   β γ λ     1 ρ (K + µI) M -M 1 M I -1 -1 M -1 R     β γ λ   +    r -ρ ζ M (u (k) -d (k) ) -ρ ζ (u (k) -d (k) ) ρ ζ 1 (u (k) -d (k) ) -1      β γ λ   subject to γ 0 (

2

 2 ψ,α,u,e J(ψ, α, u, e) subject to u = α where ψ = ψ lin + ψ nlin with ψ lin (m λ ) = α m λ subject to e = r -(α m λ + ψ nlin (u m λ )) α 0 and 1 α = 1.

Fig. 1 .

 1 Fig. 1. Convergence of the proposed algorithms.

  max

Fig. 2 .

 2 Fig. 2. Classification map with 10% training data. From left to right: ground-truth, FCLS (61.36%), K-Hype (71.39%), Proposed (73.38%)

  .00002 0.0758±0.0019 0.0604±0.0017 0.0134±0.0002 0.1137±0.0032 0.1428±0.0039 ExtM 0.0079±0.0001 0.0312±0.0013 0.0601±0.0016 0.0157±0.0003 0.0575±0.0024 0.1427±0.0040 K-Hype 0.0208±0.0004 0.0349±0.0013 0.0446±0.0020 0.0231±0.0004 0.0307±0.0008 0.0398±0.0012 Proposed 0.0123±0.0002 0.0288±0.0009 0.0159±0.0003 0.0225±0.0005 0.0245±0.0006 0.0299±0.0008

Table 1 .

 1 RMSE comparison

Table 2 .

 2 Classification accuracies with abundance features allowed us to evaluate the performance of unmixing algorithms via an alternative way that benefits from the ground truth reference of classification tasks. The scene used in our experiment is the celebrated Indian Pines region captured by AVIRIS. In this experiment, we used FCLS, K-Hype and the proposed algorithm for unmixing-based feature extraction. Training set was constructed by randomly selecting 5%, 10% and 15% of pixels. The comparative results are illustrated in Table2and Figure2. These classification results clearly highlight the advantage of nonlinear unmixing algorithms, and in particular the proposed algorithm.

	Training percentage	5%	10%	15%
	FCLS	56.41 61.36 62.32
	K-Hype	67.67 71.39 74.68
	Proposed	68.64 73.38 74.18

The sum-to-one constraint, which may not be necessary in some contexts, can be removed without effecting the structure of the algorithm.


Sandy Mahfouz
email: sandy.mahfouz@utt.fr

Farah Mourad-Chehade

Paul Honeine
email: paul.honeine@utt.fr

Hichem Snoussi
email: hichem.snoussi@utt.fr

Joumana Farah
email: farah.mourad@utt.fr

Kernel-based localization using fingerprinting in wireless sensor networks

Keywords: fingerprinting, kernel methods, localization, wireless sensor networks

Indoor localization is an important issue in wireless sensor networks for a very large number of applications. Recently, localization techniques based on the received signal strength indicator (RSSI) have been widely used due to their simple and low cost implementation. In this paper, we propose an algorithm for localization in wireless sensor networks based on radiolocation fingerprinting and kernel methods. The proposed method is compared to another well-known localization algorithm in the case of real data collected in an indoor environment where RSSI measures are affected by noise and other interferences.

I. INTRODUCTION

Wireless sensor networks (WSNs) have become a major research field during the last years. They are composed of a large number of tiny autonomous wireless sensors that are spatially distributed and can communicate and exchange data [START_REF] Akyildiz | A survey on sensor networks[END_REF], [START_REF] Sohraby | Wireless Sensor Networks: Technology, Protocols, and Applications[END_REF]. The sensors are mostly battery powered; therefore, they have a limited lifetime, and hence all built-in algorithms must focus on reducing energy consumption. Nowadays, WSNs are rapidly gaining importance in many fields of applications, especially in military, environmental, and healthcare domains [START_REF]Special issue on self-organizing distributed collaborative sensor networks[END_REF]. Indeed, WSNs can be efficiently used for monitoring or tracking of enemies in a battlefield [START_REF] Lee | Wireless sensor network design for tactical military applications : Remote large-scale environments[END_REF], and for real-time forest fire detection [START_REF] Yu | Real-time forest fire detection with wireless sensor networks[END_REF]. Sensor networks can also be used to detect foreign chemical agents in the environment [START_REF] De | Wireless chemical sensor networks for air quality monitoring[END_REF].

In all applications, awareness of location information is fundamental since collected data are meaningless without any geographical context. A direct way to locate the nodes is the use of GPS (Global Position System) devices; however, this solution is not practical because of the high cost of integrating GPS receivers and also of the limited efficiency of GPS signals in indoor environments. An alternative solution consists of considering two types of sensors, anchors and nonanchor nodes. Anchors have known fixed locations, whereas non-anchors nodes, or simply nodes, have unknown positions and thus need to be localized. Many localization algorithms using anchors have been proposed. They are mainly based on estimating the distances between the anchors and the nodes, using the Time Of Arrival (TOA) [START_REF] Shen | Passive location using toa measurements[END_REF], the Time Difference Of Arrival (TDoA) [START_REF] Eickhoff | 3d-accuracy improvements for tdoa based wireless local positioning systems[END_REF], the Angle of Arrival (AoA) [START_REF] Rong | Angle of arrival localization for wireless sensor networks[END_REF], or the Received Signal Strength Indicator (RSSI) [START_REF] Parker | Vehicular node localization using receivedsignal-strength indicator[END_REF], [START_REF] Zanca | Experimental comparison of rssi-based localization algorithms for indoor wireless sensor networks[END_REF]. The RSSI methods are based on the attenuation of the signal strength over distance. They have become popular due to their low-power consumption and cost competitiveness. For instance, in [START_REF] Mourad | Model-free interval-based localization in manetes[END_REF] and [START_REF]A robust localization algorithm for mobile sensors using belief functions and interval theory[END_REF], interval-based methods using RSSI measurements are proposed; these methods perform an outer approximation of the solution area, leading to boxes guaranteed to include the actual location. In [START_REF] Essoloh | Distributed localization in wireless sensor networks as a pre-image problem in a reproducing kernel hilbert space[END_REF] and [START_REF] Honeine | Preimage problem in kernel-based machine learning[END_REF], the authors introduce a distributed strategy for localization in wireless sensor networks by solving the pre-image problem in a reproducing kernel Hilbert space.

Another technique for RSSI-based localization is fingerprinting [START_REF] Lin | Performance comparison of indoor positioning techniques based on location fingerprinting in wireless networks[END_REF], [START_REF] Robles | 3d fingerprint-based localization for wireless sensor networks[END_REF], [START_REF] Meng | Secure and robust wi-fi fingerprinting indoor localization[END_REF]. The advantage of location fingerprinting technique is that it takes into account the stationary characteristics of the environment such as multipath propagation, i.e. wall attenuation. This technique consists of two phases: the offline phase and the online phase. The offline phase is a training phase where the goal is to build a training database for each reference location; a radio map is then created where each reference node has RSSI measures received from the anchors associated to its position. To determine the position of a node in the online phase, RSSI values collected from the node are processed and the estimated location coordinates are calculated based on the knowledge built during the offline phase. A well-known algorithm based on location fingerprinting is the Weighted K-Nearest Neighbor (WKNN) algorithm [START_REF] Koyuncu | A 2d positioning system using wsns in indoor environment[END_REF]. In the basic form of this algorithm, the RSSI values collected by the node whose position is unknown are compared to the samples in the database and the output is the weighted average of the K nearest samples; the nearness indicator for this method is based on the Euclidean distance between RSSIs.

In this paper, we propose a localization technique combining location fingerprinting and kernel methods. The proposed method employs fingerprinting information to define a model that associates to the measured data the positions where they are made. To this end, a fingerprinting database containing RSSI measures for each reference node associated with its position is constructed. The model is then defined using kernel methods ( [START_REF] Aronszajn | Theory of reproducing kernels[END_REF], [START_REF] Poggio | The mathematics of learning: Dealing with data[END_REF]) and the constructed database. Then, in order to localize the nodes, each one collects, while moving, RSSI measures at every time-step. These RSSI measures are used with the defined model to estimate the position of the considered node. The rest of the paper is organized as follows.

Section II outlines the proposed localization algorithm. Section III describes the use of kernel methods in our algorithm. Section IV illustrates the effectiveness of the method by comparing it to KNN-based methods. Finally, Section V concludes the paper.

II. DESCRIPTION OF THE PROPOSED ALGORITHM

We consider two types of sensors: anchors and non-anchor nodes. Anchors have known fixed locations, denoted by a i , i ∈ {1, . . . , N a }, where N a is their number. Non-anchor nodes or simply nodes are mobile having unknown positions, and hence needing to be localized regularly. Nodes positions at a given time t are denoted by x j (t), j ∈ {1, . . . , N x }, where N x is their number.

The proposed algorithm consists of two phases: offline and online phases. In the offline phase, N offline positions, denoted by p ℓ , ℓ ∈ 1, ..., N , are generated in a uniformly distributed manner in the studied environment. A node is then placed at each of these offline positions to perform measurements. Indeed, each anchor broadcasts signals in the network with the initial power ρ 0 . By placing the node at a position p ℓ , it becomes possible to measure the Received Signal Strength Indicator (RSSI) of anchors signals at this location. The RSSI is assumed to follow the Okumura-Hata model [START_REF] Medeisis | On the use of the universal okumurahata propagation prediction model in rural areas[END_REF] given by:

ρ a i ,p ℓ = ρ 0 -10n P log 10 a i -p ℓ + ε i,ℓ (1)
In this expression, ρ a i ,p ℓ (in dBm) is the power received from the anchor a i by the node at position p ℓ , ρ 0 is the initial power (in dBm), n P is the path-loss exponent, a ip ℓ is the Euclidian distance between the considered node at position p ℓ and the anchor a i , and ε i,ℓ is the noise affecting the RSS measures. Therefore, a database of n couples (ρ ℓ , p ℓ ) is obtained, where ℓ ∈ {1, . . . , N }, and where ρ ℓ = (ρ a 1 ,p ℓ , . . . , ρ a Na ,p ℓ ) ⊤ is the vector of received powers at position p ℓ .

Based on the information collected in the database, the objective is to define a set of functions ψ d (•) : IR Na → IR, that associate to each entry ρ ℓ the corresponding output p ℓ,d , where d ∈ {1, . . . , D}, D being the space's dimension and p ℓ,d is an element of p ℓ = (p ℓ,1 , . . . , p ℓ,D ) ⊤ . As we will show in Section III, kernel methods provide an elegant framework to determine these functions. Once ψ d (•) has been defined, a node's position, at a given time t, can be directly estimated in the online phase. Indeed, consider a node j to be localized at a given time t. This node receives signals from the N a anchors in the environment, measures their RSSI values, and stores them in a vector ρ j (t). The estimated coordinates of the node j are then given by xj,d (t) = ψ d (ρ j (t)), where d ∈ {1, . . . , D}.

III. DEFINING ψ (•) USING KERNEL METHODS

Let us consider a reproducing kernel κ : IR Na ×IR Na → IR, and let us denote by H its reproducing kernel Hilbert space (RKHS) with inner product •, • H . From the reproducing property, we can say that every ψ d (•) of H can be evaluated at any

ρ ℓ ∈ IR Na by ψ d (ρ ℓ ) = ψ d (•) , κ(•, ρ ℓ ) H .
((p ℓ,d -ψ d (ρ ℓ )) 2 + η ψ d 2 H , (2)
where η is a regularization parameter that controls the tradeoff between the training error and the complexity of the solution, as measured by the norm in the RKHS. To minimize the functional in (2), we take the functional derivative with respect to ψ d (•) and set it to zero [START_REF] Poggio | The mathematics of learning: Dealing with data[END_REF], namely

-2 N ℓ=1 (p ℓ,d -ψ d (•) , κ(•, ρ ℓ ) H ) κ(•, ρ ℓ ) + 2ηψ d (•) = 0.
This optimality condition leads to the expression

ψ d (•) = N ℓ=1 (p ℓ,d -ψ d (•) , κ(•, ρ ℓ ) H ) η κ(•, ρ ℓ ). (3)
Having

ψ d (ρ ℓ ) = ψ d (•) , κ(•, ρ ℓ )
H from the reproducing property, we can write (3) as follows:

ψ d (•) = N ℓ=1 (p ℓ,d -ψ d (ρ ℓ )) η κ(•, ρ ℓ )
Finally, we get the following model :

ψ d (•) = N ℓ=1 α ℓ,d κ(ρ ℓ , •), (4)
where

α ℓ,d = p ℓ,d -ψ d (ρ ℓ ) η . (5)
By injecting expression (4) in the cost function ( 2), we get a least squares optimization problem that can be written in matrix form as follows :

(K + ηI N )α d = P d , (6)
where

α d = (α 1,d , . . . , α N,d ) ⊤ , P d = (p 1,d , . . . , p N,d )
⊤ , I N is the N -by-N identity matrix, K is a N -by-N matrix whose (i, j)-th entry is given by κ(ρ i , ρ j ).

In this paper, we consider the case of the Gaussian kernel given by:

κ(ρ i , ρ j ) = exp -ρ i -ρ j 2 2σ 2 ,
where σ is the bandwidth of the Gaussian kernel. In this case, equation ( 4) approximates the unknown function by a weighted superposition of Gaussian functions. The weight α ℓ,d of each Gaussian is taken such as to minimize the error on the training set. The bandwidth of the Gaussian kernel, together with the regularization parameter, control the degree of smoothness, noise tolerance, and generalization of the solution.

IV. SIMULATION AND RESULTS

In order to evaluate the performance of our method, we tested it using various scenarios considering a 2D environment. In the first part of this section, results are shown in the case of Matlab-generated data where RSSI values are obtained using the Okumura-Hata model [START_REF] Medeisis | On the use of the universal okumurahata propagation prediction model in rural areas[END_REF]. However, experimental studies have shown that most of the state of the art localization algorithms, once deployed in real testbeds, achieve much worse performance than what was expected. Therefore, it is essential to test localization algorithms in the presence of real data. In the second part of this section, the proposed method is tested with real data gathered in a 10m × 10m real indoor environment [START_REF] Zanca | Rssi measurements in indoor wireless sensor networks[END_REF]. In the third section, the results using the proposed method are compared to ones obtained with the well-known K-Nearest Neighbor localization algorithm, with different weighting functions.

A. Evaluation of the proposed method on simulated data

For the training phase, 100 offline positions and 16 anchors are generated in a uniformly distributed manner over a 1000m × 1000m area (Fig. 1). Anchors are assumed to be static. The RSSI values needed to construct the database of N couples (ρ ℓ , p ℓ ), where ℓ ∈ {1, . . . , N } (N = 100), are generated using the Okumura-Hata model given by (1) with parameters set to n P = 4 and ρ 0 = 150dBm.

We generated a trajectory and calculated the RSSI values received by the moving node using the same RSSI model without noise; we then estimated the position of the node during its movement using our localization algorithm. Fig. 2 shows the generated trajectory and the estimated one. The two trajectories are superimposed. Let the estimation error be the root mean squared distances between the exact positions and the estimated ones. Then the estimation error is equal to 1.063m for η = 2 -20 and σ = 2 7 . Results for different values of η and σ are shown in Table I. In Fig. 3, a zero-mean Gaussian white noise of standard deviation σ b = 1dB is added to the RSSI values. The estimated trajectory follows the actual one and the estimation error is equal to 10.656m for η = 2 -12 and σ = 2 6 ; results for different parameter values are shown in Table II.

B. Evaluation of the proposed method using real data

The set of collected measurements used in this section are available from [START_REF] Zanca | Rssi measurements in indoor wireless sensor networks[END_REF]. The authors deployed 48 EyesIFX nodes in a room measuring approximately 10m × 10m. The nodes are uniformly distributed over the area and suspended at approximately 75 centimeters from the ceiling. Furniture and people present in the room cause multi-path interferences affecting the collected RSSI values. Fig. 4 shows the topology of the testbed. We chose 5 nodes to be anchors in our simulation, leaving us with 43 nodes to use as offline positions. However, to get better results, we used the present data to generate 57 new offline positions in order to get a total of 100 offline positions (Fig. 5). To evaluate our algorithm in the case of a moving node, we also generated a trajectory using the data. Fig. 6 shows the generated trajectory and the estimated one. The estimation error is equal to 0.195m for η = 2 -19 and σ = 16. Fig. 7 shows how η and σ affect the estimation error, where η = 2 s with s ∈ {-20, -19, • • • , -1} and σ = 2 t with t ∈ {1, 2, • • • , 10}.

(σ b = 1dB) σ η = 2 -20 η = 2 -18 η = 2 -16 η = 2 -14 η = 2 -12

C. Comparison to Weighted K-Nearest Neighbor algorithm

The Weighted K-Nearest Neighbor (WKNN) algorithm is a well-known localization algorithm based on the fingerprinting technique and the Euclidean distances between the RSSI values received by the node to be localized and the fingerprints in the database. Let us denote by D j,ℓ = ρ j (t)ρ ℓ the Euclidean distance between the RSSI vector ρ j (t) of the node j we want to localize at time t and ρ ℓ of the database where ℓ ∈ {1, . . . , n}. Let I j be the set of indices of ρ ℓ yielding the K smallest distances D j,ℓ at time t. I j of cardinality K corresponds to the set of the K nearest neighbors (offline positions) to the node j at time t. Then node j's position is estimated by:

xj (t) = k∈Ij w k p k ,
where p k is one of the nearest offline positions neighboring the node, and w k is the normalized weight factor depending on the k th minimal distance. The weight factor is important in contributing to position accuracy; in fact, neighbors with high weight contribute more to position coordinates than elements with low weight, since weights and distances vary in inverse proportion to one another. Weight values used in simulations are listed in Table III.

Table IV shows the estimation error obtained while estimating the same trajectory of Fig. 2 (without adding noise to the RSSI values) using WKNN algorithm for K ∈ {1, • • • , 15} and for five different weights (Table III). The lowest error obtained using WKNN in this case is equal to 16.863m when

"A" 1/K "B" 1/D j,k v∈I j 1/D j,v "C" 1/D 2 j,k v∈I j 1/D 2 j,v "D" 1/D 3 j,k v∈I j 1/D 3 j,v "E" exp(-D j,k ) v∈I j exp(-D j,v )
using weight "C" and K = 4. Using our proposed method, we found an estimation error of 1.063m for η = 2 -20 and σ = 2 7 (see Table I).

After adding a zero-mean Gaussian white noise of standard deviation σ b = 1dB to the RSSI values, we get the results in Table V, where the estimation error is presented in the case of WKNN algorithm. The lowest error obtained using WKNN is equal to 17.644m when using weight "C" and K = 5. Our proposed method gives an estimation error of 10.656m for η = 2 -12 and σ = 6 (see Table II) which is also smaller than the best result obtained with WKNN.

We now consider the same scenario as in Section IV.B (real RSSI values are used). Table VI shows the estimation error obtained while estimating the trajectory of Fig. 6 using WKNN algorithm. The lowest error obtained is equal to 0.340m when using weight "E" and K ∈ {10, 11, 12, 13}. Our proposed method gives an estimation error of 0.195m for η = 2 -19 and σ = 16 which is a better result than the one obtained with WKNN.

V. CONCLUSION

In this paper, we presented a new localization algorithm based on radio-location fingerprinting and kernel methods. The simulation results show that the proposed method allows a more accurate localization than the WKNN algorithm, for different weight models and parameters. We also analyzed the behavior of our algorithm in the presence of real data, and the results were satisfying. However, we cannot but notice that the results obtained in the case of Matlab simulated data without the additive noise are more than excellent, leaving us with the conclusion that, in the case of an accurate radio channel model, the proposed model technique overcome many existing localization techniques in terms of accuracy.

  The function ψ d (•) must minimize the mean quadratic error between the model's outputs ψ d (ρ ℓ ) and the desired outputs p ℓ,d : N ℓ=1

Fig. 2 .

 2 Fig. 1. Generated topology

Fig. 3 .

 3 Fig. 3. Estimation of the trajectory in the presence of noise

Fig. 5 .

 5 Fig. 4. Original testbed

Fig. 6 .Fig. 7 .

 67 Fig.6. Estimation of the trajectory (real data)

TABLE I .

 I ESTIMATION ERROR (IN METERS) FOR DIFFERENT VALUES OF η AND σ -MATLAB-SIMULATED DATA WITHOUT NOISE

	σ	η = 2 -20	η = 2 -18	η = 2 -16	η = 2 -14	η = 2 -12
	2 4	32.764	32.764	32.763	32.757	32.735
	2 5	6.352	6.353	6.371	6.505	6.875
	2 6	3.172	2.974	2.670	2.487	4.119
	2 7	1.063	2.455	5.235	7.066	7.953
	2 8	6.846	7.442	7.693	7.930	10.961
	2 9	7.666	7.806	10.671	25.648	41.617
	2 10	10.620	25.664	41.902	49.256	51.238

TABLE II .

 II ESTIMATION ERROR (IN METERS) FOR DIFFERENT VALUES OF η AND σ -MATLAB-SIMULATED DATA WITH NOISE

TABLE III .

 III LIST OF SELECTED WEIGHT FACTORS w k

	Type	Expression of w k

TABLE IV .

 IV

		ESTIMATION ERROR (IN METERS) BETWEEN REAL AND
	ESTIMATED TRAJECTORY -MATLAB-SIMULATED DATA WITHOUT NOISE
	K	"A"	"B"	"C"	"D"	"E"
	1	48.376	48.376	48.376	48.376	48.376
	2	36.538	27.603	27.647	29.801	38.935
	3	28.882	20.535	20.108	23.030	37.641
	4	29.114	19.865	16.863	19.516	37.245
	5	28.969	20.062	16.966	19.242	37.230
	6	35.913	24.720	18.424	18.959	37.228
	7	35.705	25.655	19.472	19.592	37.227
	8	38.163	27.503	20.287	19.618	37.227
	9	39.979	29.327	21.120	19.636	37.227
	10	37.017	27.591	19.970	18.965	37.227
	11	36.504	27.353	19.829	18.728	37.227
	12	36.870	27.059	19.468	18.456	37.227
	13	38.415	27.958	19.664	18.337	37.227
	14	39.871	28.105	19.491	18.087	37.227
	15	39.939	28.295	19.825	18.304	37.227
	TABLE V.	ESTIMATION ERROR (IN METERS) BETWEEN REAL AND
	ESTIMATED TRAJECTORY -MATLAB-SIMULATED DATA WITH NOISE
			(σ b = 1dB)
	K	"A"	"B"	"C"	"D"	"E"
	1	49.205	49.205	49.205	49.205	49.205
	2	36.639	28.390	27.536	29.491	40.328
	3	30.285	22.605	21.218	23.506	39.121
	4	27.907	20.414	17.805	20.167	38.768
	5	28.792	20.814	17.644	19.474	38.698
	6	36.218	26.035	19.792	19.737	38.691
	7	36.134	26.708	20.363	19.878	38.690
	8	38.966	28.550	20.931	19.734	38.690
	9	40.597	30.016	21.945	20.046	38.690
	10	37.429	28.559	20.986	19.413	38.690
	11	37.629	28.677	20.838	19.126	38.690
	12	37.398	27.831	20.041	18.618	38.690
	13	39.135	28.912	20.384	18.541	38.690
	14	40.254	29.231	20.512	18.443	38.690
	15	40.334	29.079	20.599	18.471	38.690

TABLE VI .

 VI ESTIMATION ERROR (IN METERS) BETWEEN REAL AND ESTIMATED TRAJECTORY -REAL DATA (SECTION IV.B)

	K	"A"	"B"	"C"	"D"	"E"
	1	0.500	0.500	0.500	0.500	0.500
	2	0.543	0.486	0.475	0.475	0.467
	3	0.480	0.431	0.426	0.433	0.406
	4	0.550	0.432	0.410	0.417	0.400
	5	0.564	0.438	0.404	0.408	0.383
	6	0.606	0.453	0.405	0.408	0.378
	7	0.596	0.448	0.391	0.400	0.362
	8	0.618	0.460	0.382	0.393	0.348
	9	0.644	0.480	0.384	0.389	0.341
	10	0.636	0.485	0.388	0.389	0.340
	11	0.615	0.483	0.394	0.393	0.340
	12	0.620	0.483	0.394	0.393	0.340
	13	0.649	0.496	0.398	0.396	0.340
	14	0.714	0.527	0.407	0.397	0.341
	15	0.758	0.547	0.413	0.399	0.341

ACKNOWLEDGMENT

The authors would like to thank the Champagne-Ardenne Region in France for funding this project.


Sandy Mahfouz
email: sandy.mahfouz@utt.fr

Farah Mourad-Chehade

Paul Honeine
email: paul.honeine@utt.fr

Joumana Farah
email: farah.mourad@utt.fr

Hichem Snoussi
email: hichem.snoussi@utt.fr

Localisation par fingerprinting et méthodes à noyaux dans les réseaux de capteurs sans fil

Keywords: localisation, réseaux de capteurs, fingerprinting, méthodes à noyaux

Introduction

Les réseaux de capteurs sans fil sont des réseaux composés d'un grand nombre de capteurs intelligents collaborant ensemble. Ces derniers sont munis d'un capteur physique mesurant une certaine grandeur physique, d'unités de traitement et de stockage de données et d'un module de communication sans fil [START_REF] Akyildiz | A survey on sensor networks[END_REF]. L'une des problématiques principales des réseaux de capteurs sans fil est la localisation géographique des capteurs qui s'avère primordiale à la fois pour les protocoles de routage que pour la plupart des applications. Plusieurs approches ont été considérées pour la résolution du problème de localisation [START_REF] Honeine | Preimage problem in kernel-based machine learning[END_REF][START_REF] Mourad | Anchor-based localization via interval analysis for mobile ad-hoc sensor networks[END_REF]. Dans ce papier, nous proposons une approche basée sur le fingerprinting [START_REF] Meng | Secure and robust wi-fi fingerprinting indoor localization[END_REF], permettant de construire une base de données par rapport à l'état du réseau. Une modélisation à l'aide des méthodes à noyaux [START_REF] Honeine | Preimage problem in kernel-based machine learning[END_REF] est ensuite réalisée, celleci ayant recours à la base de données construite au préalable. Une fonction est ainsi définie qui, pour un ensemble de mesure en entrée, donne en sortie une estimation de la position du capteur.

Description de l'approche

Le réseau est formé de deux types de capteurs : les capteurs fixes à positions connues appelés ancres, notés a i , i ∈ {1, . . . , N a }, et les capteurs à positions inconnues appelés noeuds, notés x m (t), m ∈ {1, . . . , N x }. Les ancres sont supposées uniformément distribuées par rapport à la région de surveillance. L'approche proposée est composée de deux étapes, la première permettant de construire une base de données à l'aide du fingerprinting [START_REF] Meng | Secure and robust wi-fi fingerprinting indoor localization[END_REF], la deuxième permettant de réaliser un apprentissage du système. A la première étape, n positions p ℓ , ℓ ∈ {1, . . . , n}, dites positions offline, sont générées d'une manière uniformément espacée le long de la zone de surveillance. Les ancres émettent par la suite des signaux dans le réseau avec une même puissance initiale ρ 0 . En plaçant les noeuds aux positions offline, il devient possible de mesurer les puissances des signaux reçus, celles-ci décroissant d'une manière monotone avec l'accroissement de la distance parcourue par le signal. Notons ρ ai,p ℓ la puissance mesurée à la réception du signal émis par l'ancre a i et reçu par un noeud placé à la position p ℓ . Une base de données est ainsi obtenue, celle-ci est formée de n couples de mesures (ρ ℓ , p ℓ ), où ℓ ∈ {1, . . . , n}, ρ ℓ = (ρ a1,p ℓ , . . . , ρ a Na ,p ℓ ) T est le vecteur puissance relevé en p ℓ .

Une fois la base de données offline est construite, l'objectif est de définir une fonction ψ (•) qui associe à chaque entrée ρ ℓ la sortie correspondante p ℓ , avec ψ : IR Na → IR d où d est la dimension de l'espace. Comme préconisé dans la littérature des méthodes dites à noyaux [START_REF] Honeine | Preimage problem in kernel-based machine learning[END_REF], la fonction recherchée prend la forme suivante :

ψ(•) = n ℓ=1 α ℓ k(ρ ℓ , •),
K (α 1 , . . . , α n ) T = (p 1 , . . . , p n ) T ,
où K étant la matrice de Gram de taille n × n dont le (i, j) ième élément est k(ρ i , ρ j ). Une fois cette fonction identifiée, un noeud mobile j mesure les puissances des signaux qu'il reçoit des ancres à l'instant t et les stocke dans un vecteur ρ m (t). Sa position estimée est alors donnée par ψ (ρ m (t)).

Dans ce papier, nous considérons le cas d'un noyau Gaussien donné par :

k(ρ i , ρ j ) = exp -ρ i -ρ j 2 2σ 2

Conclusions et perspectives

La méthode proposée permet une localisation des capteurs dans le réseau en utilisant la technique de fingerprinting, suvie d'un apprentissage à l'aide des méthodes à noyaux. Les résultats de simulation illustrent clairement les performances de la méthode proposée. Des améliorations vis à vis des performances du modèle en présence du bruit seront à considérer pour permettre une localisation encore plus précise.

FIG. 1 -

 1 FIG.1-Estimation de la trajectoire sans bruit.

FIG. 2 -

 2 FIG.2-Estimation de la trajectoire avec bruit.

0 ,

 0 où • désigne la norme Euclidienne. Afin d'illustrer les performances de la méthode proposée, nous considérons un noeud mobile se déplaçant dans une zone carrée de 1000m × 1000m. Nous considérons également 16 ancres et 100 positions offline uniformément espacées le long de la zone de surveillance. La Figure 1 illustre la trajectoire estimée à l'aide de données non bruitées. Elle montre également la trajectoire réelle du noeud ainsi que les positions offline et celles des ancres. L'erreur quadratique entre les positions réelles et les positions estimées vaut 3, 39m. La Figure 2 illustre la trajectoire estimée en utilisant des données bruitées. En ajoutant du bruit, l'erreur augmente à 15, 2m, avec toujours une bonne estimation.