Une méthode statistique de détection d’anomalie pour les modèles à espace d’état non linéaires

05/09/2017
Publication e-STA e-STA 2008-2
OAI : oai:www.see.asso.fr:545:2008-2:19826
DOI :

Résumé

Une méthode statistique de détection d’anomalie pour les modèles à espace d’état non linéaires

Métriques

25
8
136.16 Ko
 application/pdf
bitcache://f12feb7f1916065e365fa0ae6b2b365db525fe2f

Licence

Creative Commons Aucune (Tous droits réservés)
<resource  xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
                xmlns="http://datacite.org/schema/kernel-4"
                xsi:schemaLocation="http://datacite.org/schema/kernel-4 http://schema.datacite.org/meta/kernel-4/metadata.xsd">
        <identifier identifierType="DOI">10.23723/545:2008-2/19826</identifier><creators><creator><creatorName>Ghislain Verdier</creatorName></creator><creator><creatorName>Nadine Hilgert</creatorName></creator><creator><creatorName>Jean-Pierre Vila</creatorName></creator></creators><titles>
            <title>Une méthode statistique de détection d’anomalie pour les modèles à espace d’état non linéaires</title></titles>
        <publisher>SEE</publisher>
        <publicationYear>2017</publicationYear>
        <resourceType resourceTypeGeneral="Text">Text</resourceType><dates>
	    <date dateType="Created">Tue 5 Sep 2017</date>
	    <date dateType="Updated">Tue 5 Sep 2017</date>
            <date dateType="Submitted">Fri 20 Apr 2018</date>
	</dates>
        <alternateIdentifiers>
	    <alternateIdentifier alternateIdentifierType="bitstream">f12feb7f1916065e365fa0ae6b2b365db525fe2f</alternateIdentifier>
	</alternateIdentifiers>
        <formats>
	    <format>application/pdf</format>
	</formats>
	<version>33646</version>
        <descriptions>
            <description descriptionType="Abstract"></description>
        </descriptions>
    </resource>
.

1 Une m´ethode statistique de d´etection d’anomalie pour les mod`eles `a espace d’´etat non lin´eaires Ghislain Verdier1,2, Nadine Hilgert1 et Jean-Pierre Vila1 1 UMR Analyse des Syst`emes et Biom´etrie - SUPAGRO/INRA - 2 place Viala, 34060 Montpellier, France 2 I3M, UMR CNRS 5149, Universit´e Montpellier II, cc 51, Place Eug`ene Bataillon, 34095 Montpellier, France {ghislain.verdier,nadine.hilgert,jean-pierre.vila}@supagro.inra.fr R´esum´e : Le test CUSUM est un des tests statistiques les plus efficaces pour la d´etection d’anomalie dans les syst`emes dynamiques. Ce test repose sur la connaissance de la forme analytique de la densit´e conditionnelle des observations par rapport `a leurs valeurs pass´ees. Cette forme n’est g´en´eralement pas connue pour toute une classe de mod`eles largement rencontr´es dans l’industrie : les mod`eles `a espace d’´etat non lin´eaires. Nous proposons ici une r`egle de d´ecision de type CUSUM construite `a partir d’estimations des vraisemblances conditionnelles inconnues. Ces estimations sont obtenues en utilisant un filtre particulaire `a noyaux de convolution, plus efficace qu’un filtre de Kalman ´etendu. La r`egle obtenue est appliqu´ee sur un proc´ed´e de d´epollution biologique en vue de d´etecter un biais sur un capteur de concentration en biomasse. Mots cl´es : D´etection d’anomalie, R`egle du CUSUM, Filtre particulaire, Proc´ed´e de d´epollution. I. INTRODUCTION DEPUIS de nombreuses ann´ees, la d´etection d’anomalie, ou d´etection de panne, occupe une grande place dans la supervision de processus, dans des applications de type industrie agro-alimentaire, sismologie, syst`eme de guidage ou encore pour la surveillance de proc´ed´es de d´epollution biologique. Les grandes disciplines traitant ce probl`eme sont l’Automatique et l’Intelligence Artificielle, mais depuis quelques ann´ees, l’int´erˆet pour l’utilisation de m´ethodes sta- tistiques, compl´ementaires des m´ethodes traditionnelles, n’a cess´e d’augmenter (voir Basseville et Nikiforov [1]). Ces approches stochastiques permettent de prendre en compte une partie de l’incertitude sur le syst`eme par des bruits sur le mod`ele. Elles permettent de plus, de construire des r`egles de d´ecision statistiques en fonction de contraintes fix´ees par l’exp´erimentateur : temps moyen entre deux fausses alarmes fix´e, probabilit´e de faux diagnostic fix´ee... Une des toutes premi`eres r`egles de d´ecision statistiques mises au point est la r`egle des sommes cumul´ees (r`egle du CUSUM). Elle fut propos´ee par Page [3] en 1954, pour d´etecter un changement de valeur d’un param`etre dans la densit´e d’une suite d’observations ind´ependantes. Elle a ´et´e par la suite adapt´ee `a des donn´ees d´ependantes, de sorte `a d´etecter un mˆeme changement de valeur dans un mod`ele dynamique. Cette r`egle n´ecessite pour sa construction, la connaissance de la densit´e conditionnelle des observations par rapport `a leurs valeurs pass´ees. La forme analytique de cette densit´e conditionnelle n’est g´en´eralement pas accessible pour toute une classe de mod`eles souvent rencontr´es en pratique : les mod`eles dynamiques non lin´eaires `a espace d’´etat. L’objet de cette communication est de pr´esenter une m´ethode de d´etection de changement de valeur de param`etre dans des mod`eles `a espace d’´etat non lin´eaires, correspondant `a une anomalie dans le processus mod´elis´e. Notre approche consiste `a utiliser conjointement la th´eorie du filtrage par- ticulaire, et plus particuli`erement une m´ethode de filtrage originale bas´ee sur une estimation non param´etrique du filtre optimal, avec des int´egrations de Monte Carlo, pour estimer la vraisemblance conditionnelle des observations. Nous sommes alors en mesure de construire une r`egle de d´ecision statistique construite sous la forme du CUSUM. Dans ce travail, la proc´edure de d´etection obtenue est appliqu´ee sur un probl`eme r´eel, `a savoir la d´etection d’anomalie sur un proc´ed´e de digestion ana´erobie. Dans un premier temps, nous pr´esentons le contexte de la d´etection de changement d’un point de vue statistique et nous rappelons l’´ecriture du CUSUM. Le cas des mod`eles dynamiques `a espace d’´etat est ensuite ´evoqu´e et l’algorithme de filtrage-d´etection est introduit. Enfin, dans une derni`ere partie, des r´esultats sur des donn´ees r´eelles sont pr´esent´es. II. D´ETECTION DE CHANGEMENT ET R`EGLE DU CUSUM L’objectif de la d´etection de changement de valeur de param`etre, ou de panne, est de pouvoir rep´erer, le plus rapide- ment possible, le passage d’un ´etat de fonctionnement normal, appel´e r´egime H0, `a un ´etat de fonctionnement anormal, ou panne, appel´e r´egime H1, pour le syst`eme sous surveillance. On peut caract´eriser ce changement de r´egime par un chan- gement de valeur d’un param`etre θ dans la mod´elisation du syst`eme, qui va passer de la valeur θ0 sous H0 `a la valeur θ1 sous H1. Plus pr´ecis´ement, consid´erons le cadre suivant : sous r´egime de fonctionnement normal, la densit´e condition- nelle des observations Yn sachant le pass´e Y1, ..., Yn−1 est pθ0 (Yn|Y1, ..., Yn−1). A un instant tp inconnu, le param`etre θ prend la valeur θ1 et les observations ont alors pour densit´e e-STA copyright © 2008 by see Volume 5 (2008), N°2 pp 13-16 2 conditionnelle pθ0,θ1,tp (Yn|Y1, ..., Yn−1) pour n ≥ tp, densit´e qui d´epend de l’instant de changement tp. On suppose pour l’instant que les valeurs θ0 et θ1 sont connues. Le temps d’arrˆet de la r`egle du CUSUM est alors d´efini comme le premier instant o`u la statistique de test gn, en moyenne proche de 0 sous H0, positive et croissante sous H1, franchit un seuil h fix´e par l’exp´erimentateur : ta = inf{n : gn ≥ h}, avec gn = max 1≤j≤n n i=j log pθ0,θ1,j(Yi|Y1, ..., Yi−1) pθ0 (Yi|Y1, ..., Yi−1) . Cette r`egle de d´ecision est assez efficace et poss`ede de plus des propri´et´es d’optimalit´e (Lai [2]). Cependant, elle n´ecessite pour son utilisation, la connaissance analytique des vraisemblances conditionnelles des observations par rapport `a leurs valeurs pass´ees sous H0 et sous H1. III. D´ETECTION DANS LES MOD `ELES DYNAMIQUES `A ESPACE D’´ETAT A. Mod`eles dynamiques `a espace d’´etat Un mod`ele dynamique `a espace d’´etat est constitu´e de va- riables non observ´ees, Xn, les variables d’´etat, et de variables d’observation, Yn, qui permettent d’obtenir de l’information sur l’´etat du syst`eme. Dans un cadre stochastique, un tel mod`ele s’´ecrit sous la forme : Xn = f(Xn−1, θ, vn) Yn = g(Xn, θ, wn) , (1) o`u v et w sont respectivement les bruits des ´equations d’´etat et d’observation. Nous supposons que f et g sont des fonctions connues et θ est le param`etre caract´erisant la panne. Nous supposons par ailleurs que la valeur θ0 est connue et que la valeur θ1, qui est fixe, est inconnue mais appartient `a un intervalle Θ1. Il s’agit du cadre d’application g´en´eralement rencontr´e en pratique. On peut en effet supposer que dans une phase pr´eliminaire `a la surveillance, le param`etre caract´erisant le fonctionnement normal du syst`eme a ´et´e bien estim´e. L’intensit´e d’une panne ´etant souvent inconnue au pr´ealable, la valeur du param`etre caract´erisant ce r´egime de panne est la plupart du temps inconnue. Pour ce type de mod`ele, la densit´e conditionnelle des observations par rapport `a leurs valeurs pass´ees n’´etant g´en´eralement pas accessible mˆeme quand les lois des bruits sont connues, il est impossible d’appliquer la r`egle du CUSUM originelle. B. La r`egle de Willsky et Jones Lorsque le mod`ele `a espace d’´etat est lin´eaire et la panne ad- ditive, une des approches les plus efficaces consiste `a appliquer un filtre de Kalman. Cette approche a ´et´e ´etudi´ee par Willsky et Jones [8] qui d´efinissent un processus d’innovation gaussien de moyenne nulle sous H0 et de moyenne non nulle sous H1. Le probl`eme de la d´etection de changement de param`etre dans ce mod`ele revient donc `a d´etecter un changement de moyenne dans une suite de variables gaussiennes. Il s’agit d’un probl`eme classique qui peut se traiter `a l’aide de la r`egle du CUSUM. Lorsque le mod`ele est non lin´eaire, on peut envisager d’utiliser la mˆeme approche `a l’aide d’un filtre de Kalman ´etendu. Mais cette m´ethode n’a plus aucune justification th´eorique et se montre alors inefficace si la non lin´earit´e est trop forte. C. L’algorithme de filtrage-d´etection La th´eorie du filtrage, qui permet de reconstruire l’´etat Xn du processus en fonction des observations Y1, ..., Yn dans le mod`ele (1), est largement utilis´ee dans l’industrie, notamment les outils de type filtre de Kalman. R´ecemment, de nouvelles m´ethodes de filtrage ont ´et´e mises au point pour les mod`eles non lin´eaires, il s’agit des filtres `a particules dont le principe est le suivant : un grand nombre de trajectoires ( ˜X, ˜Y ) sont simul´ees suivant le mod`ele (1) et des poids sont associ´es `a ces trajectoires `a partir des observations Y . A l’instant n, l’ensemble des particules pond´er´ees ˜Xn forment une mesure empirique de la distribution de l’´etat conditionnellement aux observations. Nous consid´erons ici un filtre particulaire original, le filtre `a noyaux de convolution (Rossi et Vila [4]), qui permet, non pas d’obtenir une mesure empirique de l’´etat du syst`eme, mais une estimation fonctionnelle de la densit´e pr´edictive pθ(Xn|Y1, ..., Yn−1), not´ee ˆpN θ (Xn|Y1, ..., Yn−1), o`u N est le nombre de particules utilis´es pour l’estimation non pa- ram´etrique de pθ. Cette estimation, combin´ee `a une m´ethode d’int´egration de Monte Carlo classique, permet d’estimer la densit´e conditionnelle des observations pθ(Yn|Y1, ..., Yn−1). Il suffit donc de consid´erer l’ensemble des mod`eles qui vont caract´eriser tous les ´etats de fonctionnement possibles du syst`eme : d’une part un mod`ele d´ecrivant le fonctionnement normal du proc´ed´e, ∀ n ≥ 1 : Xn = f(Xn−1, θ0, vn) Yn = g(Xn, θ0, wn) , (2) et d’autre part `a l’instant n, n mod`eles caract´erisant le fonc- tionnement sous r´egime de panne, chacun reli´e `a un instant de panne possible j, 1 ≤ j ≤ n : mod`ele j : Xi = f(Xi−1, θ0, vi) Yi = g(Xi, θ0, wi) , pour i < j (3) et    Xi θi = f(Xi−1, θ1, vi) θi−1 Yi = g(Xi, θ1, wi) , pour i ≥ j. (4) Remarque : En th´eorie du filtrage, une approche classique pour traiter les param`etres inconnus est de les consid´erer comme des variables al´eatoires, que l’on rajoute `a l’´etat du syst`eme. On estime alors ces param`etres `a partir du filtre `a particules. C’est l’approche suivie ici en posant θi = θi−1 (puisque θ1 est suppos´e constant) dans l’´etat du syst`eme (mod`ele (4)). Le mod`ele (2) va nous permettre d’estimer la vraisemblance conditionnelle des observations pθ0 (Yi|Y1, ..., Yi−1) pour tout e-STA copyright © 2008 by see Volume 5 (2008), N°2 pp 13-16 3 i ≥ 1. En effet, puisque, pθ0 (Yi|Y1, ..., Yi−1) = pθ0 (Yi|xi).pθ0 (xi|Y1, ..., Yi−1) dxi, la vraisemblance conditionnelle sera donc estim´ee par int´egration de Monte Carlo : ˆlN,m 0,i = 1 m m k=1 pθ0 (Yi|xi,0(k)) pθ0 (Yi|Y1, ..., Yn−1), quand N et m sont grands, o`u xi,0(k) sera un ´echantillon al´eatoire issu de l’estimation ˆpN θ0 (Xi|Y1, ..., Yi−1). Cette es- timation est obtenue grˆace au filtre `a noyau de convolution appliqu´e au mod`ele (2). N correspond au nombre de parti- cules du filtre et m repr´esente la taille de l’´echantillon pour l’int´egration de Monte Carlo. De la mˆeme fac¸on, tous les autres mod`eles (3)-(4) ca- ract´erisant les comportements possibles sous r´egime de panne, permettent d’obtenir les estimations des vraisemblances condi- tionnelles : pour tout j ≥ 1 et i ≥ j : pH1,j(Yi|Y1, ..., Yi−1) = pH1 (Yi|xi, θ).pH1,j(xi, θ|Y1, ..., Yi−1) dxidθ. et alors, ˆlN,m 1,i,j = 1 m m k=1 pH1,j(Yi|xi,j(k), θi,j(k)) pH1,j(Yi|Y1, ..., Yi−1) o`u (xi,j(k), θi,j(k))k=1,...,m est un ´echantillon al´eatoire de ˆpN H1,j(xi, θ|Y1, ..., Yi−1). Nous proposons alors une r`egle de d´ecision construite de la mˆeme fac¸on que celle du CUSUM : ˆta = inf    n : max 1≤j≤n n i=j log ˆlN,m 1,i,j ˆlN,m 0,i ≥ h    . (5) o`u h est choisi tel que le temps moyen entre deux fausses alarmes soit ´egal `a une constante fix´ee par l’exp´erimentateur. Nous avons d´emontr´e les propri´et´es optimales de cette r`egle de d´ecision (Verdier et al. [7]) dans le cas simplifi´e o`u le param`etre θ1 est connu. IV. APPLICATIONS Nous avons ´etudi´e le comportement de la r`egle ˆta sur des donn´ees r´eelles provenant d’un proc´ed´e de digestion ana´erobie. Le syst`eme dynamique que nous avons consid´er´e mod´elise le fonctionnement d’un bior´eacteur de traitement des boues de vinification apr`es les vendanges. Ce mod`ele a ´et´e mis au point au Laboratoire de Biotechnologie de l’Environnement (LBE) de l’INRA `a Narbonne au cours des dix derni`eres ann´ees. Le mod`ele original (Steyer et Bernard [5]) met en jeu deux biomasses (X1 et X2) et deux substrats (S1 et S2). Nous nous sommes int´eress´es `a la premi`ere r´eaction qui consiste en la d´egradation du substrat S1 par la biomasse X1, selon le mod`ele :    X1(i + 1) = X1(i) + T (µ1 − αD(i))X1(i) + v1 i+1 S1(i + 1) = S1(i) + T (D(i)(Sin 1 − S1(i)) −k1µ1X1(i)) + v2 i+1, (6) 0 500 1000 1500 2000 2500 3000 3500 4000 −5 0 5 10 15 20 25 30 35 40 temps, t Qin (L/h) Fig. 1. Trac´e du d´ebit th´eorique (en fonc´e) et du d´ebit mesur´e (en clair) Qin o`u X1 et S1 sont les concentrations en bact´eries acidog`enes et en substrat organique respectivement. v1 et v2 sont des bruits blancs gaussiens de variances respectives 10−6 et 10−4 . D est le taux de dilution : D = Qin V o`u Qin et V sont respectivement le d´ebit d’alimentation (voir figure 1) et le volume de travail du bior´eacteur (constant et ´egal `a 350L). D est la variable de contrˆole du syst`eme (6). µ1 est le taux de croissance de la biomasse et suit une loi de Monod : µ1(S1) = µmax S1 KS1 + S1 , avec µmax = 1.2 et KS1 = 8.875. Sin 1 est la concentration initiale en substrat et vaut 9g/L. T est le pas de discr´etisation, l’intervalle de temps entre deux mesures ´etant de deux mi- nutes. Le terme α repr´esente la proportion de biomasse non soumise `a l’effet de dilution et vaut 0.5. k1 est un rendement de conversion et est fix´e `a 42.12. Nous nous sommes int´eress´es `a des dysfonctionnements sur le capteur mesurant la concentration en substrat. Le mod`ele (6) repr´esente ainsi l’´equation d’´etat du mod`ele `a espace d’´etat et l’´equation d’observation est donn´ee par : CS1 (i + 1) = S1(i + 1) − θ + wi+1 (7) o`u CS1 repr´esente la mesure du capteur, w est le bruit sur le capteur suppos´e gaussien et de variance 10−3 , et θ caract´erise le dysfonctionnement, `a savoir ici l’apparition d’un biais entre la concentration r´eelle S1 et la mesure fournie par le capteur. On suppose que sous H0, θ = θ0 = 0, le biais est nul. Et sous H1, θ = θ1 qui est inconnue mais appartient `a l’intervalle [0.6; 1.5]. Le premier graphe de la figure 2 repr´esente les donn´ees mesur´ees de la concentration S1 (en clair) ainsi que des donn´ees simul´ees suivant le mod`ele `a espace d’´etat (6)-(7) sous l’hypoth`ese qu’il n’y a pas de dysfonctionnement du capteur (en fonc´e). On voit donc apparaˆıtre un biais entre le capteur e-STA copyright © 2008 by see Volume 5 (2008), N°2 pp 13-16 4 2600 2800 3000 3200 3400 3600 3800 4000 4200 4400 0.5 1 1.5 2 2.5 3 temps, t CS 1 Concentration en substrat S1 2600 2800 3000 3200 3400 3600 3800 4000 4200 4400 0 500 1000 1500 temps, t gt Statistique de test Fig. 2. D´etection de l’apparition d’un biais sur le capteur `a l’aide de l’algorithme de filtrage-d´etection et le mod`ele vers le pas de temps t = 4000 (voir aussi Steyer et Bernard [5]). Le deuxi`eme graphe de la figure 2 repr´esente l’´evolution de la statistique de test de la r`egle de d´ecision propos´ee. On voit que la statistique de test augmente de fac¸on significative `a partir de t = 3950, ce qui permet de d´etecter la panne tr`es rapidement, `a l’aide d’un seuil obtenu par exemple, de mani`ere adaptative (Verdier et al. [6]). V. CONCLUSION Nous proposons une m´ethode de d´etection d’anomalie pour les mod`eles dynamiques `a espace d’´etat non lin´eaires ins- pir´ee de la r`egle du CUSUM. Cette approche utilise une m´ethode de filtrage particulaire originale qui peut s’appliquer dans des conditions plus g´en´erales que les filtres usuels. La proc´edure propos´ee poss`ede des propri´et´es d’optimalit´e et son comportement en simulation montre qu’elle est globalement plus efficace que la r`egle classique de Willsky et Jones dans le cas de mod`eles non lin´eaires (voir Verdier et al. [7]). La mise en place de cette r`egle de d´ecision sur un proc´ed´e de d´epollution donne des r´esultats encourageants puisque la panne consid´er´ee, l’apparition d’un biais sur le capteur, est d´etect´ee tr`es rapidement. L’adaptation de cette r`egle de d´ecision au cas du diagnostic, c’est `a dire lorsqu’on consid`ere un ensemble de K > 1 pannes possibles sur le proc´ed´e, est directe. L’objectif est alors de d´etecter une panne et de la localiser parmi les K. Il suffit de consid´erer autant de mod`eles (3)-(4) qu’il y a de pannes possibles, chacun ´etant reli´e `a un type de panne. REFERENCES [1] M. Basseville and I. Nikiforov, Detection of Abrupt Changes. Theory and Application. Prentice-Hall, 1993. [2] T. L. Lai, “Information bounds and quick detection of parameter changes in stochastic systems,” IEEE Trans. Inform. Theory, vol. 44, pp. 2917– 2929, Nov. 1998. [3] E. S. Page, “Continuous inspection schemes,” Biometrika, vol. 41, pp. 100–115, 1954. [4] V. Rossi and J.-P. Vila, “Nonlinear filtering in discrete time : a particle convolution approach,” An. Inst. Stat. Univ. Paris, vol. 3, pp. 71–102, 2006. [5] J.-P. Steyer and O. Bernard, “An exemple of the benefits obtained from the long term use of mathematical models in wastewater biological treatment,” in Proc. of the 4th MATHMOD International Symposium on Mathematical Modelling, Vienna, Austria, 2003, pp. 245–251. [6] G. Verdier, N. Hilgert, and J.-P. Vila, “Adaptive threshold computation for CUSUM-type procedures in change detection and isolation problems,” Computational Statistics and Data Analysis, En r´evision. [7] ——, “Optimality of CUSUM rule approximations in change-point de- tection problems - applications to nonlinear state-space systems,” IEEE Trans. Inform. Theory, Soumis. [8] A. S. Willsky and H. L. Jones, “A generalized likelihood ratio approach to detection and estimation of jumps in linear systems,” IEEE Trans. Automat. Contr., vol. AC-21, pp. 108–112, Feb. 1976. e-STA copyright © 2008 by see Volume 5 (2008), N°2 pp 13-16