Julie Rivet

Guillaume Tochon

Serge Meimon

Michel Paques

Thierry Géraud

Michaël Atlan

Michael Atlan

Motion Compensation In Digital Holography For Retinal Imaging

Keywords: Motion Compensation, Cross-correlation, Digital Holography, Optical Coherence Tomography, Laser Doppler Holography

INTRODUCTION

Ophtalmology is an ancient science which has benefited from technological progresses in optics. The first images of eye fundus in vivo by holographic ophtalmoscopy [START_REF] Joseph | Holographic recording of a retina using a continuous wave laser[END_REF][START_REF] Wiggins | Fundus camera holography of retinal microvasculature[END_REF] paved the way to the development of time-averaged holography for retinal vessels imaging. This non-invasive and low-exposure optical interferometric imaging technology has made accessible the retinal imaging, leading to better diagnosis of several diseases like glaucoma, macular degeneration and diabetes. Later, retinal blood flow contrasts were also successfully rendered by narrowband [START_REF] Simonutti | Holographic laser doppler ophthalmoscopy[END_REF] and wideband [START_REF] Pellizzari | High speed optical holography of retinal blood flow[END_REF] digital holography.

Digital holography relies on the acquisition and the processing of interference images created by the superposition of two light waves from the same source. The first wave is a reference, whereas the second is back-scattered by the retina. Several facets of digital holography have been developped during the past few years to take maximum advantage of phase variations between both waves which contain numerous information about the retina. Optical Coherence Tomography (OCT) [START_REF] Fercher | Optical coherence tomography-principles and applications[END_REF][START_REF] Huang | Optical coherence tomography[END_REF] permits optical sectioning of retinal layers. On the other hand, laser Doppler imaging [START_REF] Magnain | Holographic laser doppler imaging of microvascular blood flow[END_REF] focuses on the detection of local dynamics.

Digital holography faces several issues altering images quality. Eye motion during image acquisition distorts the signal. In order to stabilize images of retinal features, image registration can be used. It consists in finding the geometric transformation to align two images. Medical image registration techniques have been developped during the past years [START_REF] Viergever | A survey of medical image registration -under review[END_REF]. Although most of the literature in this field concerns radiological modalities, several researches have been led in retinal imaging. The most widely used technique is featurebased registration of eye fundus images [START_REF] Liu | Retinal image registration via feature-guided gaussian mixture model[END_REF][START_REF] Ramli | Feature-based retinal image registration using d-saddle feature[END_REF]. However in our case, feature-based registration is not suitable because objects of interest can disappear or even be replaced by another structure if the movement of the eye is too large, which could bias the search of feature points. This article is focused on cross-correlation stabilization to compensate lateral and axial motion of the whole image structure in laser Doppler holography and swept-source holographic OCT, respectively.

HOLOGRAPHIC IMAGING

Setup

Optical interferograms are recorded with a laser Doppler instrument [START_REF] Simonutti | Holographic laser doppler ophthalmoscopy[END_REF] in human eyes, and with a holographic OCT instrument in synthetic phantoms. For the laser Doppler setup, the source is a laser diode emitting at the wavelength λ = 785 nm. The camera records 2048 × 2048-pixel images corresponding to the (x, y) plane, at a frame rate of ω S /(2π) = 80 Hz with 8 bit/pixel quantization. For OCT, the beam is emitted by a tunable laser and varies from ω 1 to ω 2 , which are linked to the wavelengths λ 1 = 870 nm and λ 2 = 820 nm by ω = 2πc/λ, where c is the speed of light, with a sweep time T = 0.5s. The camera records 1024 × 1024-pixel images at a frame rate of ω S /(2π) = 512 Hz with 16 bit/pixel quantization. Digitized interferograms from the camera are processed in real-time with the Holovibes1 software to compute and visualize holograms, and saved for offline processing.

Acquisition of interferograms

A Mach-Zehnder interferometer has been used to record interferogams. It consists in making interfere two beams from the same laser source. The source is split between reference and object arms. The light wave from the object arm is backscattered once the sample is reached and intereferes in the camera plane with the beam in the reference arm. The interferogram I acquired in the camera plane (x, y) at time t is expressed as:

I(x, y, t) = |E(x, y, t) + E LO (x, y, t)| 2 , (1)
where the object field is E and the reference field is E LO . The phase difference between both optical fields is noted φ.

The off-axis configuration of the setups allows a spatial separation between the constant term, and the cross-beating interferometric contribution, which will produce the image of interest, and its conjugate [START_REF] Verrier | Off-axis digital hologram reconstruction: some practical considerations[END_REF] (see Fig. 2a).

Processing of optically-acquired interferograms

Interferogram rescaling with wavelength. In order to avoid the distortion of the signal, the impact of the sweep in OCT needs to be considered to propagate the interferograms I from the camera plane to the retina or image plane. In fact the size of the pixels in the image plane depends on the wavelength of the beam and of the distance between camera and image planes [START_REF] Ferraro | Controlling image size as a function of distance and wavelength in fresnel-transform reconstruction of digital holograms[END_REF][START_REF] Lédl | Surface topography measurement by frequency sweeping digital holography[END_REF]: d = λz/(N d), where d and d are the lateral size of a pixel from camera plane and from image plane, respectively, z is the distance between both planes and N is the number of pixels in one lateral dimension. Because of the sweep of the source, the pixels of the image plane shrink with wavelength (arrow (1) of Fig. 1). To circumvent lateral field variation with wavelength, each interferogram is resampled by linear interpolation of the calculation grid with a different pitch (arrow (2) of Fig. 1). The rescaled interferogram is:

I (x, y, t) = I(xλ/λ 1 , yλ/λ 1 , t), ( 2
)
where λ is the current wavelength. In the case of laser Doppler imaging, the optical wavelength λ is kept constant, hence the interferogram does not need to be rescaled:

I ≡ I.
Spatial demodulation by Fresnel transform. The propagation of the fields from camera to image plane is carried out by a Fresnel transform [START_REF] Verrier | Off-axis digital hologram reconstruction: some practical considerations[END_REF] (arrow (3) of Fig. 1), which gives the hologram H(x, y, t):

H(x, y, t) = i λz exp (-ikz) I (x , y , t) × exp -iπ λz (x -x ) 2 + (y -y ) 2 dx dy . ( 3
)
The argument of the complex-valued hologram in the crosscontribution region is the phase difference φ between the optical fields E and E LO in the image plane.

Temporal demodulation by Fourier transform. In sweptsource OCT, the instantaneous beating frequency ∂φ/∂t = ω scales up linearly with axial depth z [START_REF] Fercher | Optical coherence tomography-principles and applications[END_REF], whereas in laser Doppler, it describes local velocities of the scatterers [START_REF] Magnain | Holographic laser doppler imaging of microvascular blood flow[END_REF].

Hence for both methods, temporal signal demodulation is performed by short-time Fourier transform (STFT) (arrow (4) of Fig. 1). Time-and-space-dependent spectrograms

)j (1) (2) (3) (4) (5)
Fig. 1: Sketch of the process of image propagation. The raw data in the camera plane are transferred in the image plane and are expressed with depth.

H(x, y, t, ω) are calculated from the stream of intermediate holograms H(x, y, t):

H(x, y, t, ω) = H(x, y, τ )g T (t -τ )e -iωτ dτ, (4)
where g T (t) is a time gate of width T at time t. Then, the envelope | H| of the STFT of H is formed (arrow (5) of Fig. 1).

The quantities | H(x, y, t, ω)| are the images in laser Doppler and OCT on which motion compensation will be performed.

MOTION COMPENSATION METHODS

Lateral Motion Compensation

We form Doppler-contrast images at the Nyquist frequency ω S /2 by calculating Eq. ( 4) with a time gate g T of two points (T = 25 ms). To compensate lateral motion (in (x, y) plane), we use a cross-correlation-based stabilization algorithm. The human eye is constantly subject to deviations, saccades and tremors. Added to heartbeat and respiration motion, those movements can shift or even remove the object of interest from the image during several frames. Thus, a processing chain has been built in order to keep important structures in the center of the en-face image. A reference image f 1 is made by averaging n = 10 consecutive images, and compared to a moving average f 2 of

f 1 (x, y) = n t=1 | H(x, y, t, ω S /2)| / n, f 2 (x, y) = n t=1 | H(x, y, t + τ, ω S /2)| / n. (5)
The quantitative comparison is performed with a normalized and centered cross-correlation [START_REF] Padfield | Masked object registration in the Fourier domain[END_REF]. The difference between the position of the maximum value and the center of the correlation matrix γ gives the displacement between the two images. Resulting shifts describe retina lateral motion.

Axial Motion Compensation

We form OCT images by calculating Eq. ( 4) with a time gate g T of 256 points (T = 0.5 s). Axial motion during the recording process has a negative impact on the reconstruction of OCT images: when an axial drift occurs, the depth information encoded in the beating frequency of the interferogram is inaccurate, adding an offset phase to the signal and decreasing in-depth accuracy. The principle of the compensation method is to identify where the phase shift due to the axial motion occurs by using STFT [START_REF] Hillmann | Common approach for compensation of axial motion artifacts in swept-source OCT and dispersion in Fourier-domain OCT[END_REF] on 20-point windows. Indeed, a global axial displacement of the sample during acquisition can be spotted in the time Fourier domain, where it corresponds to a frequency shift.

For y = y 0 , the reference image f 1 and the sliding image f 2 are defined as the modulus of the STFT of H(x, y 0 , t) whose sub-windows are fixed or moving, respectively:

f 1 (x, ω) = | H(x, y 0 , τ )g T (t 0 -τ )e -iωτ dτ |, f 2 (x, ω) = | H(x, y 0 , τ )g T (t 0 -τ )e -iωτ dτ |, (6)
where g T is a 20-point Hanning window, and t 0 and t 0 are the centers of the moving and the reference windows, respec- tively. The quantitative comparison between f 1 and f 2 is performed by calculating cross-correlation whose maximum position gives a frequency shift ∆ω(t 0 ) corresponding to the derivative of the phase caused by movement at time t 0 . Axial motion compensation is then performed by substracting the cumulated phase shift t 0 t 0 ∆ω(t)dt from the phase φ(t 0 ). This procedure is repeated for all times t 0 in the time gate g T used in Eq. ( 4) of width ω S T /(2π) = 256 points [START_REF] Hillmann | Common approach for compensation of axial motion artifacts in swept-source OCT and dispersion in Fourier-domain OCT[END_REF].

RESULTS

Fig. 2 shows retinal vessels acquired by laser Doppler imaging. Fig. 2a is the whole 2048 × 2048 image. The crossbeating interferometric contribution, which is the object of interest in the right bottom side of the image, is spatially separable from the other parts because of the off-axis configuration. Fig. 2b shows the focus on retinal vessels. Ves-sels are moving in the image and sometimes disappear during several frames. Then, vessels do not appear clearly on the time-averaged image of 324 consecutive frames (Fig. 2c). After lateral stabilization, vessels are more visible on the timeaveraged image: motion compensation is efficient. Although motion compensation is less efficient when the image changes too often, it improves overall image quality.

Fig. 3 and Fig. 4 are swept-source OCT images of samples composed of 1.5 mm diameter glass beads rolled in a single tape layer. Fig. 3 shows the in-depth profile (x, ω) of this sample corresponding to the y-axis represented by dotted lines in Fig. 4. Fig. 4 shows the en-face (x, y) images at two different depths corresponding to the dotted lines in Fig. 4. Fig. 4a and 4b correspond to the red dotted line (depth of 1.1 mm) in Fig. 3, and Fig. 4c and 4d correspond to the green dotted line (depth of 1.5 mm) in Fig. 3. In Fig. 4a and4b, the target layer is located between tape and beads: the beads on the left are at the same elevation as the tape layer on the right which starts to be sectioned.

A difference of lateral resolution is visible on images of Fig. 4. The interferogram rescaling allows the contours of the beads to be cleaner. Comparing Fig. 3a and Fig. 3b shows an improvement in axial resolution, which is also noticeable on en-face images: in Fig. 4c, beads on the right seem almost at the same elevation than beads on the left, because the depth accuracy is low, while in Fig. 4d, right beads clearly belong to a different layer. Axial motion correction can be seen in Fig. 3c: the accuracy is improving and the different layers are better separated. Besides, the optical bench will be strenghtened to reduce mechanical noise, in order to further improve the axial accuracy of depth images.

CONCLUSION AND PROSPECTS

We have demonstrated lateral and axial motion compensation in laser Doppler holography and holographic OCT by crosscorrelation stabilization, respectively. This method is suited to our images and cancels efficiently the effects of motion. Besides, the reported results are non-iterative and compatible with real-time processing at high throughput on graphics processing units. A real-time implementation of lateral motion compensation has been implemented2 . The reported results in motion compensation pave the way towards the design of high resolution computational imaging for the retina in realtime by digital holography.

Aside from motion, retinal imaging suffers from phase distortions caused by the eye: the light backscattered by the retina crosses lenses in the eye, which can cause aberrations. As a future work, we will try to correct optical aberrations.

Fig. 2 :

 2 Fig. 2: Retinal image stabilization by lateral (x, y) motion compensation of intensity holograms (Section 3.1). The white scale bars indicate 200 µm.

Fig. 3 :Fig. 4 :

 34 Fig. 3: Axial motion compensation in swept-source holographic OCT of glass beads rolled in a tape layer. The registered holograms are slices in the (x, ω) plane, for the ycut in Fig. 4. Figures represent the sample: (a) without any correction, (b) with interferogram rescaling (Eq. (2)), and (c) with interferogram rescaling and axial motion compensation (Section 3.2). The white scale bars indicate 0.2 mm.

http://holovibes.com

https://youtu.be/RhVPXBnPhXc

Funding: This work has been supported the European Research Council (ERC Synergy HELMHOLTZ, grant agreement #610110).


Agriculture is not only a food producing system but a mechanism that needs to be well driven for sustainability in landscape: Some teachings from French and Taiwanese situation

Agriculture is not only a food producing system but a mechanism that needs to be well driven for sustainability in landscape: Some teachings from French and Taiwanese situation

Yves Michelin, B. Cicilien, Kuang Yu Wang

To cite this version:

Yves Michelin, B. Cicilien, Kuang Yu Wang. Agriculture is not only a food producing system but a mechanism that needs to be well driven for sustainability in landscape: Some teachings from French and Taiwanese situation. PECSRL Conference :Mountains, Upland, Lowlands -European landscapes from an altitudional perspective, Sep 2016, Innsbruck Seefeld, Austria. hal-02183125

Agriculture is not only a food producing system but a mechanism that needs to be well driven for sustainability in landscape: some teachings from French and Taiwanese situations Our study analyses the difficuhies they face and the scracegies ch ey use, wich a view co draw up guidelines for becter management of suscainable local developmenc and landscape.

HAL

  is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come

  farming practices are nor only producing food but are also closely connecced to local culture and environment. Nowadays th ey have been transformed inco serving the re gional and global networks for food production, often valorized by competitiveness of productivity and price. ln the dominant idea of the free market and the global di vision of food production, agricultural landscapes gradually get detached from indi vidualiry and local contexts. At the same time, urbanization has a major inBuence on the change of agriculture and tural landscapes. Many rural landscapes are exploited and developed for leisure purposes and tourist attractions in order to serve the needs of an urban population for rural fancasy. But the meaning of rural agricultural landscapes can be more than that of food producing machines or a 'Disneyland'. Biodiversity, cultural heritage, social net works, etc., could also be values related to agricultu.ral practices while local commu nities seek new developments to keep the cultural heritage alive. Our study focuses on the comparison of three contrasted situations: ln Yilan County (Taiwan), the spring onion production has a special quality de rived from its good match with the local environmental conditions. However, growing rourism and rural housing have encouraged land speculation and in creased the land prices, resulting in the decrease of farming. Sorne efforts were made to support local farming practices, conserving landscape characters and promoting sustainable rural developments. ln the French Cevennes terraces, sweet onion production has increased since it has been protected by a specific regulation based on the guarantee of origin with the support of specific laws inside a national park. As a result, terraces have been renovated and new producers develop their activities. ln the suburb of Clermont-Ferrand (France), a traditional garlic production is disappearing due to the lack of specific regulation or a d yn arnic network of farm ers in a context of competition with intensive crop production and development demands.

Wang, Kuang-Yu 1 ; Michelin, Yves' & Cicilien, Benoit' 1 Chuan Yuan Universfr)I Taiwan 1 VetA gr oSup, France

	K eyw ords: local food, agricultural landscape, public policies, Taiwan, French Mas sif Central
	Rural agricultural landscapes have a long history in Europe and East Asia where


Clement Leroy
email: clement.leroy@irisa.fr

Eric Anquetil
email: eric.anquetil@irisa.fr

Nathalie Girard
email: nathalie.girard@irisa.fr

ParaFIS:A new online fuzzy inference system based on parallel drift anticipation

Keywords: Evolving Fuzzy System EFS, Non Stationary Environment NSE, brutal drift adaptation, long-life learning

This paper proposes a new architecture of incremental fuzzy inference system (also called Evolving Fuzzy System -EFS). In the context of classifying data stream in non stationary environment, concept drifts problems must be addressed. Several studies have shown that EFS can deal with such environment thanks to their high structural flexibility. These EFS perform well with smooth drift (or incremental drift). The new architecture we propose is focused on improving the processing of brutal changes in the data distribution (often called brutal concept drift). More precisely, a generalized EFS is paired with a module of anticipation to improve the adaptation of new rules after a brutal drift. The proposed architecture is evaluated on three datasets from UCI repository where artificial brutal drifts have been applied. A fit model is also proposed to get a "reactivity time" needed to converge to the steady-state and the score at end. Both characteristics are compared between the same system with and without anticipation and with a similar EFS from stateof-the-art. The experiments demonstrates improvements in both cases.

I. INTRODUCTION

More and more data are generated by stream over a long period leading to the necessity for classification methods to include long-life learning in their algorithm to stay reliable in time. As an example, a command-gesture system [START_REF] Bouillon | Online active supervision of an evolving classifier for customized-gesture-command learning[END_REF], where the user draws a gesture to apply a command, produces this kind of data stream. The set of gestures is chosen by the user and can change at any time. For example a novice user often draw slowly a gesture, and once he is used to make it, he can speed up the gesture leading to a shift of the target concept in the feature space. This modification is called concept drift. Such environment is said to be non stationary. It means that the probability distribution that generates data (here the user) depends on time. To maintain high performance in gestures recognition over time, the classification system must adapt its parameters and/or structure to the concept drift. Thus, the learning process has to be incremental. For each new data, the system extends its knowledge (i.e., it learns and adapts its parameters). Moreover, to get a constant complexity in time and avoid a memory overflow, along a long-life learning, data should not be saved. Recent researches have shown that Evolving Fuzzy Systems -EFS (also called Incremental Fuzzy Inference System) could deal with such environment thanks to their high structural flexibility [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF], [START_REF] Angelov | Simplified fuzzy rule-based systems using nonparametric antecedents and relative data density[END_REF]. The structure (number of rules, and antecedents/consequent parameters) of these fuzzy rule-based systems evolve with the stream of data. Especially, we distinguish between two scales of adaptation to concept drift. One concerns the adaptation of the rules' parameters and tackles incremental drifts (or smooth drifts). The other relates to the adaptation of the structure, by the addition or deletion of rules, and tackles brutal drifts (brutal shifts of the data distribution). Both adaptions use two main algorithmic approaches. The first one relies on temporal (adaptive) sliding windows where only the most recent data are taking into account. As an example, the fuzzy windows concept drift adaptation method (FW-DA) [START_REF] Liu | Fuzzy time windowing for gradual concept drift adaptation[END_REF] relies on the method and shows promising results outperforming state-of-the-art. The second one weights the data according to their age and interest. An example is [START_REF] Bouillon | Decremental learning of evolving fuzzy inference systems using a sliding window[END_REF] which uses a decremental learning on the premise and conclusion parts with the Differed Directional Forgetting (DDF). However, both approaches have to figure out a "forgetting" parameter which controls the relevance of the old data. If the forgetting parameter is too high, the system will be less stable leading to lower performance. If the parameter is too low, the system will not be reactive to the change in the environment leading to bad performance. This issue is often called the plasticity-stability dilemma. This dilemma can result in instability, particularly with EFS where the number of rules, their size and their reactivity greatly depend on parameters given a priori. To prevent instability, [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF] recently proposed a new rule splitting method. The idea is to split into two, a rule which makes too many mistakes or which is too large. However, it takes time after the splitting for the rules to readjust to the local distribution of points, introducing inertia in the learning. On the contrary, [START_REF] Song | A self-adaptive fuzzy network for prediction in non-stationary environments[END_REF] proposed a method to accelerate the learning of new rules after a drift based on the generation of data from new drifted concept using a GAN [START_REF] Goodfellow | Generative adversarial nets[END_REF]. Similarly, we propose a complementary method to prevent instability of the rules while improving reactivity. To do so, we introduce a new original architecture of EFS, called ParaFIS. In ParaFIS, a generalized evolving fuzzy system (principal system) is learned synchronously with an anticipation module. For each rule of the principal system, two rules are anticipated in the anticipation module. Thus, the anticipation module enforces the system to locally adapt the distribution of points with two rules rather than just one, in order to anticipate a concept drift before it occurs. The paper is organized as follows. Section II introduces the generalized EFS and its learning model with a discussion on its drawbacks. Section III introduces our contribution, the ParaFIS evolving system. The experiments showing that the anticipation improves the plasticity of the system while keeping stability, are detailed in section IV. In this section, we propose to measure performance of the system on artificial brutal drifts with fitted parameters of a handcraft model. This evaluation protocol allows to quantify the time of reactivity of the system and its stability in the steady-state.

II. GENERALIZED EVOLVING FUZZY SYSTEMS

In this paper, we focus on the generalized evolving fuzzy system that uses a generalized version of Takagy-Sugeno fuzzy systems, already used [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF], [START_REF] Almaksour | Improving premise structure in evolving takagi-sugeno neuro-fuzzy classifiers[END_REF], [START_REF] Lemos | Multivariable gaussian evolving fuzzy modeling system[END_REF].

A. Model Architecture

A Takagy-Sugeno (TS) fuzzy system is a set of fuzzy inference rules R = {r i , 0 ≤ i ≤ N } with an antecedent part (also called premise), and a consequent part. Each rule's antecedent is defined with a prototype that is set by a cluster with a center µ i . The structure of a rule r i , is as follows:

IF x is close to µ i THEN y 1 i = l 1 i (x) .. y c i = l c i (x) (1)
With l j i a polynomial function for r i of class j; c the number of class and N the number of rules. The degree of the polynomial function is set to 1 with π j ik the polynomial coefficients (see Eq. (2)).

y j i = l j i (x) = π j i0 + π j i1 x 1 + π j i2 x 2 + .. + π j in x n = Π j i x (2)
The membership of x to a rule r i , denoted β i (x), is given by a normalized Radial Basis Function K (RBF), of the distance from x to µ i (see Eq.( 3)). The RBF is often a multivariate Gaussian or Cauchy function [START_REF] Almaksour | Improving premise structure in evolving takagi-sugeno neuro-fuzzy classifiers[END_REF], [START_REF] Lemos | Multivariable gaussian evolving fuzzy modeling system[END_REF].

β i (x) = K(||x -µ 1 || 2 ) (3)
In the generalized version of TS, the Mahalanobis distance is used to get rotated hyper-ellipsoid clusters as follows:

||x -µ i || 2 = (x -µ i )A -1 (x -µ i ) T (4)
With A the covariance matrix. Finally, the predicted class for x is given by Eq. ( 5), [START_REF] Song | A self-adaptive fuzzy network for prediction in non-stationary environments[END_REF].

class(x) = y = argmax j y j (x) (5)
Where

y j (x) = N i=1 β i (x)y j i (6)

B. Rule's adaptation

Each new incoming data x t is used to adapt the model parameters. In the premise part, only the most activated rule adapts its center and covariance matrix according to Eq. ( 7), [START_REF] Almaksour | Improving premise structure in evolving takagi-sugeno neuro-fuzzy classifiers[END_REF] where t is the number of samples of the rule.

µ t = t -1 t µ t-1 + 1 t x t (7)
A t = t -1 t A t-1 + 1 t (x t -µ t )(x t -µ t ) T (8)
The forgetting capacity is put in the equation, by setting t = min(k, tmax) (see [START_REF] Bouillon | Decremental learning of evolving fuzzy inference systems using a sliding window[END_REF]) with tmax a threshold that defined the forgetting capacity, and k the number of samples that activated the most the rules. tmax is often written with a forgetting factor α = tmax-1 tmax ∈ [0, 1] where α = 1 when there is no forgetting capacity. The consequent part is learned using a Weighted Recursive Least Square method (WRLS). The membership functions β is assumed to be almost constant to converge to the optimal solution. To reduce computation time, the local learning of the consequent part is often preferred. And, the rules are assumed to be independent to apply RLS on each one. The conclusion matrix

Π i(t) = [Π 1 i(t) , .., Π c i(t)
] of the rule r i at time t (i.e. after t data points) is recursively computed according to:

Π i(t) = Π i(t-1) + C i(t) β i (x)C i x(Y t -xΠ i(t-1) ) (9) Where C i(t) = C i(t-1) - β i (x)C i(t-1) xx T C i 1 + β i (x)x T C i x (10)
With C i a correlation matrix initialized by C i(t=0) = ΩId where Id is the identity matrix and Ω a constant often fixed to 100 (see [START_REF] Almaksour | Improving premise structure in evolving takagi-sugeno neuro-fuzzy classifiers[END_REF], [START_REF] Angelov | An approach to online identification of takagi-sugeno fuzzy models[END_REF]).

C. Rule creation condition

The adaptation of the parameters is not relevant to handle brutal drifts, as when a class must be represented by several clusters or when the target concept shift in the feature space. Dealing with brutal drifts requires the adaptation of the fuzzy system structure, like rule addition. In a context of online learning from scratch, all classes start with one prototype (i.e. one rule). The structure adaptation relies on specific conditions observed on data, via the existing rules. Several criteria are defined to detect brutal drifts. Most EFS uses a distancebased criteria [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF], [START_REF] Almaksour | Improving premise structure in evolving takagi-sugeno neuro-fuzzy classifiers[END_REF] that compares a certain threshold T (P ) (depending on parameters P ) with the distance between the prototype of the closest rule (its center µ) and the new incoming points x t . If T (P ) < dist(µ -x t ), then the rule creation criteria is met, a new rule is created over the last incoming point according to Eq. [START_REF] Gama | On evaluating stream learning algorithms[END_REF]. Two scales of adaptation co-exist in EFS, the adaptation of rule parameters and the structure adaptation. Smooth drifts are tackled by introducing forgetting capacity in the parameter adaptation whereas brutal drifts are tackled by the creation of new rules. However, the system is degraded by the parameter adaptation when a brutal drift occurs. Indeed, all rule creation conditions lead to a trade-off between speed of detection and sensitivity to noise. But, for all, it exists a time ∆T > 0 between the true occurrence of the drift, and the detection time. As shown in Figure 1 and 2, during this ∆T time, one (ore more) rule tries to adapt to the drift and changes its parameters making it less fit to its previous concept. But, when the rule creation is triggered, this previous adaptation is not cancelled making the old rule perhaps unstable. Moreover, the new rule is created over one single point (the one that triggers the rule creation) although all points during ∆T could have been used to initialize the new rule. This results in a longer time τ for the new rule adaptation to the new concept.

µ new = x t Π new = 0 cov kl = 1 100 δ kl ∀k, l ∈ [1, n] (11)

III. CONTRIBUTION: PARAFIS SYSTEM

In order to attend stability and plasticity, we present in this section the ParaFIS system.

A. Model's architecture

As shown in Figure 3, the proposed model is based on the generalized EFS similar to this described section II. This part is dedicated to smooth drifts to keep stability (no structural change). And, for each rule r i , a module of anticipation is added to deal with brutal drifts. This anticipation module is composed of two sub-rules r i1 and r i2 that have an antecedent part (a center, a covariance matrix, a Cauchy membership function) and a consequent part with c hyperplanes (c the number of classes). The system classifies the data at any time

B. Model's learning

Each new point coming into the system will be used to learn both the principal system and the anticipation module. As in the generalized EFS, the principal system will adapt the antecedent part of the most activated rule r m by updating its center and its covariance matrix using Eq. ( 7),( 8) with a factor α = 1 (no forgetting capacity). Then, the sub-rules r m1 and r m2 are updated using the same equations with α1, α2 ≤ 1. The two sub-rules have two different forgetting factors leading to two temporal scales of learning. r m1 has a low forgetting factor and is learned on the most recent data whereas r m2 has a high forgetting factor and is learned on a long history. In this way, r m1 quickly reacts to a change in the distribution of points whereas r m2 preserves the old concept with a slow adaptation. The consequent part in the principal system is learned as usual (Eq. ( 9),( 10)). In the anticipation module, r i1 , r i2 have the same consequent part as r i (the c hyperplanes y j i with j ∈ [1, c]).

C. Detection of brutal concept drifts

Contrary to current rule creation criterion which use no information of neighborhood [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF], [START_REF] Almaksour | Improving premise structure in evolving takagi-sugeno neuro-fuzzy classifiers[END_REF], we propose here to integrate a brutal drift detector based on a clustering separability criteria. The idea is to assume that if the two subrules in the anticipation module are enough separated, then a brutal drift occurred. Then, r m2 learned over the large history matches the old concept and r m1 learned over the few last points matches the new drifted concept. Eq. ( 12) presents the proposed separability criteria that is based on the covariance of both clusters.

Condition 1 ||µ i -µ j || > σ i + σ j (12)
Where, as depicetd in Figure 4, σ i (resp. σ j ) is the distance between µ i (resp. µ j ) and the hyper-ellipsoid's envelop of cluster i (resp. j), along (µ i , µ j ) axis. Besides, to force the rules to take into account a certain number of points n min (20 Fig. 4: Separability criteria using covariance of both clusters by default) before deciding to create a new rule, the following inertia criteria is added:

Condition 2 k i > n min (13)

D. Initialization of new rules based on the anticipation module

If the rule creation conditions are met for the most activated rule r m from the principal system, then r m is replaced by r m1 and r m2 in the principal system. Then for each new r mi of the principal system, two new sub-rules r mi 1 , r mi 2 are initialized in the anticipation module, as follows:

µ rm i 1 = µ rm i 2 = µ rm i A rm i 1 = A rm i 2 = A rm i Π rm i 1 = Π rm i Π rm i 2 = 0 k rm i 1 = k rm i 2 = 0 (14)
The idea is to keep the learned information of r mi in the subrule r mi 1 , and to initialize the second r mi 2 from scratch. All steps of the learning are summarized in Algorithm 1. 7)-( 8), update the antecedent part of rm with α = 1, and rm 1 , rm 2 with α1 = α2 According to Eq. ( 9)-( 10), update the consequent part of all rules ri, ri 1 and ri 2 else

Algorithm 1 Learn ParaFIS model

Replace rm by rm 1 ,rm 2 in the principal system Initialize 4 new sub-rules in the anticipation module Eq. ( 14) end if else

Create a new rule Eq. ( 11) Initialized two sub-rules in the anticipation module with (11) end if

IV. EXPERIMENTAL VALIDATION

A. Evaluation protocol 1) Prequential test with artificial brutal drift: To evaluate the performance of an online classifier, it is current to use the prequential test [START_REF] Gama | On evaluating stream learning algorithms[END_REF]. In this test, data are given one by one to the system. The system first tests the new data to get a score (1 if the class is well classified, 0 otherwise) and then learns TABLE I: Information on the datasets used in the experiments on it. In this way, all the data are used to test the system and then, to train it while maintaining independence between each phase. Scores are then averaged over a certain window (of size n) to get a smooth curve of the performance over time. This test simulates a real usecase as for many online application where the system must adapt to the behavior of a user along time. In the following experiments, n = 5 (except to plot the figures n = 100 to smooth the score).

There is no existing dataset with annotation of the nature of the drift or with the occurrence time of the drift making the evaluation and comparison of online classifiers complex. To make the task easier, the paper proposes to generate artificial brutal drifts in real data at a chosen time. To do so, each dataset is split into three sub-datasets A, B, C to create a specific data stream. As illustrated in Figure 5, A is composed of the first T 1 data belonging to n 1 classes. B (resp. C) is composed of the data in the stream between T 1 and T 2 (resp. between T 2 and T 3), this data belongs to n 2 ≤ n 1 (resp. n 3 ≤ n 1 ) different classes and are relabeled with the n 1 labels to produce the brutal drift. In this ways, n 2 brutal drifts are done at the time T 1 and n 3 others at time T 2. Thereafter, we call this approach the protocol P . In the context of a commandgesture system, the protocol is equivalent to a user changing the gesture of a command while keeping the possibility to make the old gesture.

2) Benchmark dataset: Regarding the context of commandgestures system, three datasets from the handwritten pattern recognition community have been chosen to assess the performance of the system: PenDigits [START_REF] Alimoglu | Combining Multiple Classifiers for Pen-Based Handwritten Digit Recognition[END_REF], Letters [START_REF] Frey | Letter recognition using Holland-style adaptive classifiers[END_REF] and LaViola [START_REF] Laviola | A practical approach for writerdependent symbol recognition using a writer-independent symbol recognizers[END_REF]. They are all available in the UCI machine learning repository [START_REF] Dheeru | UCI machine learning repository[END_REF]. All three datasets contain different features extracted from the handwritten patterns (digits, letters or symbols handwritten by different writers). Description of each dataset is given in Table I. These datasets are built in a static context meaning that there is no order between data. Thus, the prequential test with artificial brutal drifts can be done several times (m) by shuffling the dataset. All results given thereafter are averaged over m prequential tests with m = 100. The experimental parameters are given in Table 6: Fit of the three phases score using Eq. ( 15)

II. The n1, n2, n3 classes of each dataset follow the order of occurrence in the data file from UCI.

3) Characterization of the plasticity and stability: In order to measure the plasticity and stability of the system, we propose to fit the prequential score with an handcraft model given by Eq. [START_REF] Dheeru | UCI machine learning repository[END_REF].

y(t) = S(1 -e -t τ ) + s min (15)
Where we define a characteristic time τ , which represents the reactivity time (or the plasticity). In particular after a time t = τ , 63% of the score have been reached until the steady state. At the steady state, the score is given by S + s min . A least square method is used to determine the parameters which best fit the curve with the model. The example of fits on the Letters dataset (experimented with the protocol P) given in Figure 6 shows that this model fits well the prequential score.

B. Evaluation of the anticipation 1) Importance of covariance matrix initialization: No study has yet tackled the issue of optimizing the initialization parameters of the antecedent part. Yet, the initialization of the prototype (the antecedent part) greatly influence the score. To show it, we propose to compare three methods currently used in incremental fuzzy system to initialize the covariance of a new rule: Method I1 (Eq. ( 16), [START_REF] Bouillon | Decremental learning of evolving fuzzy inference systems using a sliding window[END_REF]), Method I2 (Eq. ( 17), [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF], [START_REF] Lemos | Multivariable gaussian evolving fuzzy modeling system[END_REF]), Method I3 (Eq. ( 16)). The center of prototypes is initialized just on the last point for all methods.

∀k, l ∈ [1, n], ∀j ∈ [1, N ] I1:
cov kl = min j (cov j kl )δ kl (16)

I2: cov kl = 1 100 * δ kl (17)
I3:

cov kl = mean j (cov j kl ) 10 (18)
To compare them, we use the generalized evolving fuzzy system as described section II with Condition 1 and Condition 2 to create rules (section III). Results are obtained from the three datasets and a fit using Eq. ( 15) is done. Results are [START_REF] Dheeru | UCI machine learning repository[END_REF] displayed in Figure 7 column A, and the fitted parameters are shown in Table III. Results show that the best choice of initialization depends on the dataset. However it is clear that the initialization of the rules widely impacts the results.

For instance, there is 5.6% of difference for the S + s min parameters between methods I3 and I2 for letters dataset at the phase B. This highlights the importance of well initializing the covariance matrix of the new rule. This can be explained by the fact the covariance matrix plays a great role in the learning. Indeed, only the antecedent part of the most activated rule is learned on the new incoming point and the activation, computed with a Mahalanobis distance, crucially depends on the shape of the covariance matrix. If it is too small, few data from the new concept will activate the new rule. If it is too big, data from other concepts may also activate the new rule.

2) Comparison with our proposition: Now, the impact of the anticipation will be evaluated.The antecedent part initialization is compared when using the best method among {I1,I2,I3} for each dataset, and using the anticipation module Eq.( 14). However the drift detector can depend on the parameters of the system as it is in ours, Eq. ( 12) (this depends on the covariance matrices). Thus different initialization of prototype change the times of detection. To avoid this bias, the drift detector is first used on the system with anticipation and the times t, at which the drifts are detected, are saved. Then, rather than using the detector, we use the saved files with all times t to create rules. In such ways, the system can be compared fairly without bias of the detection. Two sets of parameters are used for the anticipation, P ara 1 = (α 1 = 1, α 2 = 0.9) and P ara 2 = (α 1 = 1, α 2 = 0.95) with α 1 (resp. α 2 ) the forgetting factor of the sub-rules i1 (resp. i2) of the secondary system for i ∈ [1, r]. Moreover, a comparison is done with our own implementation of an EFS similar to Gen-Smart EFS [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF] called GEFS*. In this last one, the rule creation criteria (given by Eq. ( 19)-(20)) is the same than Gen-smart EFS [START_REF] Lughofer | Incremental rule splitting in generalized evolving fuzzy systems for autonomous drift compensation[END_REF] with the same rule initialization method (I2). However, there are no rule merging or rule splitting methods as there are in Gen-Smart EFS.

(x -µ i )A -1 (x -µ i ) T > r 2 i (19) r i = κp 1 √ 2 1.0 (1 -1 ki+1 ) m (20)

Fig. 1 :Fig. 2 :

 12 Fig. 1: Problems met in generalized EFS when brutal drift occurs. The antecedent structure of each rule is represented by the ellipse (the covariance matrix) and the star (the center)

Fig. 3 :

 3 Fig. 3: Architecture of our proposition ParaFIS

  Require: xnew, current fuzzy rule system containing set of rules R and set of classes C if class(xnew) ∈ C then Find max activated rule rm ∈ R according to (3) Check the condition 1 condition 2 Eq. (12-13) with rm if Condition 1 & Condition 2 == FALSE then According to Eq.(

Fig. 5 :

 5 Fig. 5: Protocol P -Generation of data stream with brutal drifts Dataset Classes Features Samples Scriptwriters Letters 26 16 20000 20 LaViola 48 50 16891 34 PenDigits 10 16 10992 44

Fig. 7 :

 7 Fig. 7: Prequential score (y-axis) with respect to data (x-axis) for the three datasets. The first (resp. the second) column refers to the first (resp. second) experiment

TABLE II :

 II Experimental parametersFig.

TABLE III :

 III Fitted parameters of the handcraft model

TABLE IV :

 IV Mean score over the three phases A,B,C