![\"\"](\"..\/..\/lise.vaudor\/Rfigures\/Ondelettes_discretes\/draw_headband-1.png\")
<\/p>\n<\/p>\n
A\u00efe a\u00efe a\u00efe!! Quelle prise de t\u00eate que les transform\u00e9es d\u2019ondelettes quand on n\u2019est pas tr\u00e8s matheux\u2026 Comprendre de A \u00e0 Z comment elles fonctionnent, en d\u00e9roulant tout le raisonnement de la d\u00e9finition d\u2019une ondelette \u00e0 l\u2019interpr\u00e9tation finale de la transform\u00e9e, c\u2019est long, c\u2019est compliqu\u00e9, c\u2019est m\u00eame parfois impossible si on ne peut pas se lancer \u00e0 corps perdu dans le sujet\u2026<\/p>\n
Ici, je pars du principe que vous non plus, vous n\u2019\u00eates pas forc\u00e9ment pr\u00eats \u00e0 faire le grand saut\u2026 alors je vais vous pr\u00e9senter une m\u00e9thode li\u00e9e aux transform\u00e9es d\u2019ondelettes discr\u00e8tes, dite Analyse Multi-R\u00e9solution<\/strong> (MRA en anglais) en me concentrant sur le r\u00e9sultat<\/em> plut\u00f4t que sur la m\u00e9thode<\/em>. J\u2019ai tout un tas de supports \u00e9parpill\u00e9s sur la m\u00e9thode<\/em> en question donc il n\u2019est pas exclu que je me fende un jour d\u2019un bon gros billet avec tout plein d\u2019\u00e9quations sur le sujet mais je ne voudrais pas vous faire peur d\u2019embl\u00e9e :-)\u2026<\/small><\/p>\n Notez que dans cette m\u00eame d\u00e9marche (de d\u00e9mocratisation des ondelettes) j\u2019ai \u00e9galement mis en ligne une appli<\/a> permettant de tester tout un tas de m\u00e9thodes d\u2019analyse du signal par les ondelettes (ondelettes discr\u00e8tes ou continues).<\/p>\n Les transform\u00e9es d\u2019ondelettes peuvent servir \u00e0 d\u00e9bruiter<\/strong>, filtrer<\/strong>, compresser<\/strong>, d\u00e9composer<\/strong> le signal et\/ou sa variance en dissociant ses composantes \u00e0 diverses \u00e9chelles<\/strong>.<\/p>\n Pour ma part, je ne m\u2019int\u00e9resse pas ces transform\u00e9es pour la compression du signal mais plut\u00f4t \u00e0 leur usage pour l\u2019analyse de s\u00e9ries<\/strong>. Par ailleurs, je vais ici me concentrer sur l\u2019analyse de s\u00e9ries (temporelles ou spatiales) 1D, excluant de ce fait l\u2019analyse d\u2019images (2D) de mon propos.<\/p>\n J\u2019aimerais commencer en montrant rapidement ce que l\u2019on peut r\u00e9aliser \u00e0 travers les transform\u00e9es d\u2019ondelettes. Consid\u00e9rons le signal suivant, qui correspond \u00e0 une s\u00e9rie de temp\u00e9ratures relev\u00e9es \u00e0 la station m\u00e9t\u00e9o de Lyon Saint-Exup\u00e9ry (temp\u00e9rature moyenne journali\u00e8re, donn\u00e9es disponibles sur ce site<\/a> et que je vous remets ici<\/a>. On devine dans ce signal une forte composante saisonni\u00e8re<\/strong> (\u00e9chelle temporelle \u201cmoyenne\u201d), \u00e0 laquelle se rajoute du \u201cbruit<\/strong>\u201d (\u00e9chelle temporelle fine), et, peut-\u00eatre, une tendance de fond<\/strong> (\u00e9chelle temporelle longue) \u00e0 l\u2019augmentation de la temp\u00e9rature moyenne (on verra si nos analyses le confirment!). De fait, j\u2019ai choisi ces donn\u00e9es car elles se pr\u00eatent tout naturellement \u00e0 une analyse multi-r\u00e9solution!!<\/p>\n La d\u00e9composition en ondelettes discr\u00e8tes permet de r\u00e9aliser une analyse multi-r\u00e9solution (MRA<\/strong>) du signal. Je passe pour l\u2019instant sur le code qui permet de r\u00e9aliser les graphiques suivants\u2026 mais j\u2019y reviendrai!<\/small> C\u2019est \u00e0 dire que l\u2019on va pouvoir d\u00e9composer<\/strong> le signal en diff\u00e9rentes composantes, qui correspondent \u00e0 des \u00e9chelles plus ou moins fines, et que cette d\u00e9composition aura l\u2019int\u00e9r\u00eat d\u2019\u00eatre additive<\/strong>.<\/p>\n Si l\u2019on consid\u00e8re la fa\u00e7on dont le signal est d\u00e9compos\u00e9, \u00e9tape par \u00e9tape, on peut \u00e9crire les \u00e9quations suivantes:<\/p>\n $$ \\begin{aligned} Y &= & & & & & & & &S1 &+ &D1\\ &= & & & & & &S2 &+ &D2 &+ &D1\\ &= & & & &S3 &+ &D3 &+ &D2 &+ &D1\\ &= & & S4 &+ &D4 &+ &D3 &+ &D2 &+ &D1\\ etc. \\end{aligned} $$<\/p>\n A la fin, on arrive toujours \u00e0 une somme de sous-signaux \u00e9gale au signal de d\u00e9part<\/strong>, chaque sous-signal nous informant sur la variation \u00e0 une \u00e9chelle donn\u00e9e, de mani\u00e8re non redondante<\/strong> par rapport aux autres sous-signaux.<\/p>\n Voyez \u00e0 quoi la \u201cpyramide\u201d d\u2019\u00e9quations ci-dessus correspond en terme de repr\u00e9sentation des variations du signal<\/strong> (j\u2019ai repr\u00e9sent\u00e9 les s\u00e9ries verticalement pour une meilleure lisibilit\u00e9):<\/p>\n Il est en fait possible de lisser encore et encore ce signal<\/strong>, jusqu\u2019\u00e0 arriver \u00e0 un smooth qui nous montre la tendance de fond. Ci-dessus je me suis arr\u00eat\u00e9e au niveau 4 car d\u2019un point de vue graphique 4 niveaux de d\u00e9composition \u00e9taient d\u00e9j\u00e0 assez lourds\u2026<\/p>\n Cela \u00e9tant, rien ne nous emp\u00eache de regrouper certains sous-signaux<\/strong>\u2026 Par exemple, on peut essayer de distinguer une composante plut\u00f4t saisonni\u00e8re (S4) au sein des variations journali\u00e8res:<\/p>\n Y<\/em>\u2004=\u2004S<\/em>4\u2005+\u2005(D<\/em>4\u2005+\u2005D<\/em>3\u2005+\u2005D<\/em>2\u2005+\u2005D<\/em>1) Ainsi, en poussant la d\u00e9composition plus loin, et en regroupant un certain nombre de niveaux, on peut essayer d\u2019isoler une tendance de fond:<\/p>\n Y<\/em>\u2004=\u2004S<\/em>10\u2005+\u2005(D<\/em>10\u2005+\u2005…\u2005+\u2005D<\/em>5)\u2005+\u2005(D<\/em>4\u2005+\u2005…\u2005+\u2005D<\/em>1) Bingo! L\u2019examen des graphiques ci-dessus r\u00e9v\u00e8le par exemple:<\/p>\n Pour r\u00e9aliser les transform\u00e9es d\u2019ondelettes discr\u00e8tes sous R j\u2019utilise un package qui s\u2019appelle Voil\u00e0 comment je r\u00e9alise mon analyse multi-r\u00e9solution, \u00e0 l\u2019aide de la fonction Examinons maintenant ce que nous retourne la fonction Cet objet comprend un slot \u201cfilter\u201d qui nous renseigne sur la nature de l\u2019ondelette utilis\u00e9e<\/strong> (car il existe plusieurs ondelettes possibles)! Ici, comme nous n\u2019avons rien pr\u00e9cis\u00e9 dans l\u2019appel \u00e0 la fonction Outre des informations sur la m\u00e9thode de calcul et la nature de l\u2019ondelette (slot \u201cfilter\u201d), l\u2019objet Cela \u00e9tant, quand on travaille dans la logique \u201ctidy\u201d, on a tr\u00e8s envie d\u2019obtenir tous ces sous-signaux comme autant de colonnes d\u2019une table. Pour faciliter cette transformation, j\u2019ai \u00e9crit une petite fonction On est d\u00e8s lors en mesure de faire notre petite tambouille (avec ggplot2, avec dplyr, etc.). Tenez, v\u00e9rifions en deux coups de cuill\u00e8re \u00e0 pot que la d\u00e9composition du signal est bien additive:<\/p>\n En outre, je peux utiliser ce tableau \u00e0 ma guise pour r\u00e9aliser un graphique montrant les composantes de mon choix:<\/p>\n A\u00efe a\u00efe a\u00efe!! Quelle prise de t\u00eate que les transform\u00e9es d\u2019ondelettes quand on n\u2019est pas tr\u00e8s matheux\u2026 Comprendre de A \u00e0 Z comment elles fonctionnent, en d\u00e9roulant tout le raisonnement de la d\u00e9finition d\u2019une ondelette \u00e0 l\u2019interpr\u00e9tation finale de la transform\u00e9e, c\u2019est long, c\u2019est compliqu\u00e9, c\u2019est m\u00eame parfois impossible si on ne peut pas se lancer \u00e0 corps perdu dans le sujet\u2026 Ici, je pars du principe que vous non.. Read More<\/a><\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[5],"tags":[],"class_list":["post-1178","post","type-post","status-publish","format-standard","hentry","category-tous-les-posts"],"_links":{"self":[{"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/posts\/1178","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/comments?post=1178"}],"version-history":[{"count":9,"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/posts\/1178\/revisions"}],"predecessor-version":[{"id":1184,"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/posts\/1178\/revisions\/1184"}],"wp:attachment":[{"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/media?parent=1178"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/categories?post=1178"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/perso.ens-lyon.fr\/lise.vaudor\/wp-json\/wp\/v2\/tags?post=1178"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}Les ondelettes discr\u00e8tes, pour quoi faire?<\/h1>\n
library(tidyverse)\ndata_ex=read.delim(paste0(dat.path,\"TG_STAID000037.txt\"),\n sep=\",\", skip=20, na=-9999)%>% \n mutate(x=lubridate::ymd(DATE),\n y=TG\/10) %>% \n group_by(x) %>% \n summarise(y=mean(y)) %>% \n na.omit()\nggplot(data_ex,aes(x=x,y=y))+\n geom_line(col=\"dark grey\")\n<\/code><\/pre>\n![\"\"](\"..\/..\/lise.vaudor\/Rfigures\/Ondelettes_discretes\/signal_brut-1.png\")
<\/p>\n\n
![\"\"](\"..\/..\/lise.vaudor\/Rfigures\/Ondelettes_discretes\/signal_details_pyramid-1.png\")
<\/p>\n![\"\"](\"..\/..\/lise.vaudor\/Rfigures\/Ondelettes_discretes\/signal_to_S4-1.png\")
<\/p>\n![\"\"](\"..\/..\/lise.vaudor\/Rfigures\/Ondelettes_discretes\/signal_to_S10-1.png\")
<\/p>\n\n
Comment faire, sous R?<\/h1>\n
wavelets<\/code>. Comme \u00e0 mon habitude, je vais aussi utiliser le tidyverse<\/code> pour r\u00e9aliser les quelques traitements suppl\u00e9mentaires dont j\u2019aurai besoin (et patchwork pour organiser mes panneaux de graphiques!)<\/p>\nlibrary(wavelets)\nlibrary(tidyverse)\nlibrary(patchwork)\n<\/code><\/pre>\nmra()<\/code> du package wavelets<\/code>:<\/p>\nmymra=mra(X=data_ex$y, \n n.levels=10,\n boundary=\"reflection\")\n<\/code><\/pre>\n\n
X<\/code> correspond au signal<\/strong> lui-m\u00eame. Pour que la transform\u00e9e d\u2019ondelettes ait un sens, il faut que la s\u00e9rie de donn\u00e9es corresponde \u00e0 des valeurs Y<\/em> r\u00e9guli\u00e8rement espac\u00e9es<\/strong> (ici dans le temps).<\/li>\nn.levels<\/code> correspond au niveau maximal de d\u00e9composition (ce niveau maximal d\u00e9pendra de la nature de l\u2019ondelette et de la longueur de la s\u00e9rie).<\/li>\nboundary<\/code> correspond \u00e0 la m\u00e9thode utilis\u00e9e \u201cau bord du signal\u201d pour calculer les transform\u00e9es d\u2019ondelettes. En effet, pour calculer les coefficients d\u2019ondelettes sur les bords du signal, on est oblig\u00e9s d\u2019allonger artificiellement le signal. Les deux choix possibles sont \u201cperiodic\u201d (on recycle le signal, du d\u00e9but \u00e0 la fin) ou \u201creflection\u201d (on retourne le signal, en miroir). Dans notre cas de s\u00e9ries de temp\u00e9ratures, l\u2019existence d\u2019une \u00e9volution entre le d\u00e9but et la fin du signal (augmentation graduelle de la temp\u00e9rature moyenne) rend la transformation \u201creflection\u201d plus pertinente\u2026<\/li>\n<\/ul>\nmra()<\/code>:<\/p>\nstr(mymra)\n\n## Formal class 'mra' [package \"wavelets\"] with 9 slots\n## ..@ D :List of 10\n## .. ..$ D1 : num [1:32532, 1] 0.418 -1.058 0.292 0.332 -0.51 ...\n## .. ..$ D2 : num [1:32532, 1] -0.538 -0.396 -0.244 0.316 1.119 ...\n## .. ..$ D3 : num [1:32528, 1] 0.1031 -0.207 -0.5185 -0.3569 -0.0928 ...\n## .. ..$ D4 : num [1:32528, 1] 2.249 2 1.697 1.147 0.477 ...\n## .. ..$ D5 : num [1:32512, 1] 0.868 0.899 0.913 0.835 0.716 ...\n## .. ..$ D6 : num [1:32512, 1] 0.861 1.047 1.215 1.3 1.343 ...\n## .. ..$ D7 : num [1:32512, 1] -3.88 -3.89 -3.89 -3.84 -3.76 ...\n## .. ..$ D8 : num [1:32512, 1] -2.47 -2.6 -2.72 -2.8 -2.87 ...\n## .. ..$ D9 : num [1:32256, 1] -2.88 -2.94 -2.99 -3.02 -3.05 ...\n## .. ..$ D10: num [1:31744, 1] 0.62 0.622 0.624 0.625 0.626 ...\n## ..@ S :List of 10\n## .. ..$ S1 : num [1:32532, 1] 13.8 14.3 14.7 15.2 15.7 ...\n## .. ..$ S2 : num [1:32532, 1] 14.3 14.7 15 14.9 14.6 ...\n## .. ..$ S3 : num [1:32528, 1] 14.5 14.4 14.2 13.8 13.2 ...\n## .. ..$ S4 : num [1:32528, 1] 12.2 12.4 12.5 12.6 12.7 ...\n## .. ..$ S5 : num [1:32512, 1] 12.9 13 13.1 13.2 13.3 ...\n## .. ..$ S6 : num [1:32512, 1] 12 12 11.9 11.9 11.9 ...\n## .. ..$ S7 : num [1:32512, 1] 15.9 15.8 15.8 15.7 15.7 ...\n## .. ..$ S8 : num [1:32512, 1] 18.4 18.4 18.5 18.5 18.6 ...\n## .. ..$ S9 : num [1:32256, 1] 11.5 11.5 11.5 11.5 11.5 ...\n## .. ..$ S10: num [1:31744, 1] 11.3 11.3 11.3 11.3 11.3 ...\n## ..@ filter :Formal class 'wt.filter' [package \"wavelets\"] with 7 slots\n## .. .. ..@ L : int 8\n## .. .. ..@ level : int 1\n## .. .. ..@ h : num [1:8] 0.0322 0.0126 -0.0992 -0.2979 0.8037 ...\n## .. .. ..@ g : num [1:8] -0.0758 -0.0296 0.4976 0.8037 0.2979 ...\n## .. .. ..@ wt.class : chr \"Least Asymmetric\"\n## .. .. ..@ wt.name : chr \"la8\"\n## .. .. ..@ transform: chr \"dwt\"\n## ..@ level : int 10\n## ..@ boundary: chr \"reflection\"\n## ..@ series : num [1:32532, 1] 14.2 13.2 15 15.5 15.2 14.4 12.3 11.7 15.2 16.4 ...\n## ..@ class.X : chr \"numeric\"\n## ..@ attr.X : list()\n## ..@ method : chr \"dwt\"\n<\/code><\/pre>\nmra()<\/code>, c\u2019est par d\u00e9faut l\u2019ondelette \u201cLeast Asymetric\u201d de longueur 8 (\u201cla8\u201d) qui a \u00e9t\u00e9 utilis\u00e9e (mais bien s\u00fbr nous aurions pu explicitement demander qu\u2019une autre ondelette soit utilis\u00e9e en utilisant l\u2019argument filter<\/code>):<\/p>\nmymra@filter\n\n## An object of class \"wt.filter\"\n## Slot \"L\":\n## [1] 8\n## \n## Slot \"level\":\n## [1] 1\n## \n## Slot \"h\":\n## [1] 0.03222310 0.01260397 -0.09921954 -0.29785780 0.80373875 -0.49761867 -0.02963553 0.07576571\n## \n## Slot \"g\":\n## [1] -0.07576571 -0.02963553 0.49761867 0.80373875 0.29785780 -0.09921954 -0.01260397 0.03222310\n## \n## Slot \"wt.class\":\n## [1] \"Least Asymmetric\"\n## \n## Slot \"wt.name\":\n## [1] \"la8\"\n## \n## Slot \"transform\":\n## [1] \"dwt\"\n<\/code><\/pre>\nmymra<\/code> comprend l\u2019ensemble des smooths et des details pour les 10 niveaux de d\u00e9composition demand\u00e9s ici. On peut acc\u00e9der \u00e0 chacun de ces sous-signaux de la mani\u00e8re suivante:<\/p>\nmymra@D$D1[1:5]\n\n## [1] 0.4177861 -1.0581337 0.2924646 0.3321684 -0.5095316\n<\/code><\/pre>\nadd_wav_data()<\/code>:<\/p>\nadd_wav_data=function(mydata,mra_obj){\n n=nrow(mydata)\n cut_to_size=function(x,n){\n if(dim(x)[1]>n){\n x=x[1:n]\n }else{\n x=c(x,rep(NA,n-dim(x)[1]))\n }\n return(x)\n }\n if(class(mra_obj)==\"mra\"){\n S=mra_obj@S %>%\n purrr::map(cut_to_size,n) %>% \n bind_cols() %>% \n mutate_all(as.vector)\n D=mra_obj@D %>%\n purrr::map(cut_to_size,n) %>% \n bind_cols() %>%\n mutate_all(as.vector)\n wav_data=bind_cols(S,D)\n }\n mydata=bind_cols(mydata,\n wav_data)\n return(mydata)\n}\n\ndata_ex_wav=add_wav_data(data_ex,mymra) \nstr(data_ex_wav)\n\n## Classes 'tbl_df', 'tbl' and 'data.frame': 16266 obs. of 22 variables:\n## $ x : Date, format: \"1975-04-20\" \"1975-04-21\" \"1975-04-22\" \"1975-04-23\" ...\n## $ y : num 14.2 13.2 15 15.5 15.2 14.4 12.3 11.7 15.2 16.4 ...\n## $ S1 : num 13.8 14.3 14.7 15.2 15.7 ...\n## $ S2 : num 14.3 14.7 15 14.9 14.6 ...\n## $ S3 : num 14.5 14.4 14.2 13.8 13.2 ...\n## $ S4 : num 12.2 12.4 12.5 12.6 12.7 ...\n## $ S5 : num 12.9 13 13.1 13.2 13.3 ...\n## $ S6 : num 12 12 11.9 11.9 11.9 ...\n## $ S7 : num 15.9 15.8 15.8 15.7 15.7 ...\n## $ S8 : num 18.4 18.4 18.5 18.5 18.6 ...\n## $ S9 : num 11.5 11.5 11.5 11.5 11.5 ...\n## $ S10: num 11.3 11.3 11.3 11.3 11.3 ...\n## $ D1 : num 0.418 -1.058 0.292 0.332 -0.51 ...\n## $ D2 : num -0.538 -0.396 -0.244 0.316 1.119 ...\n## $ D3 : num 0.1031 -0.207 -0.5185 -0.3569 -0.0928 ...\n## $ D4 : num 2.249 2 1.697 1.147 0.477 ...\n## $ D5 : num 0.868 0.899 0.913 0.835 0.716 ...\n## $ D6 : num 0.861 1.047 1.215 1.3 1.343 ...\n## $ D7 : num -3.88 -3.89 -3.89 -3.84 -3.76 ...\n## $ D8 : num -2.47 -2.6 -2.72 -2.8 -2.87 ...\n## $ D9 : num -2.88 -2.94 -2.99 -3.02 -3.05 ...\n## $ D10: num 0.62 0.622 0.624 0.625 0.626 ...\n## - attr(*, \"na.action\")= 'omit' Named int 1 2 3 4 5 6 7 8 9 10 ...\n## ..- attr(*, \"names\")= chr \"1\" \"2\" \"3\" \"4\" ...\n<\/code><\/pre>\ndata_ex_wav %>%\n mutate(somme=S10+D10+D9+D8+D7+D6+D5+D4+D3+D2+D1) %>% \n select(y,somme) %>% \n head(n=10)\n\n## # A tibble: 10 x 2\n## y somme\n## <dbl> <dbl>\n## 1 14.2 6.66\n## 2 13.2 4.79\n## 3 15 5.69\n## 4 15.5 5.84\n## 5 15.2 5.30\n## 6 14.4 4.72\n## 7 12.3 2.84\n## 8 11.7 1.28\n## 9 15.2 3.70\n## 10 16.4 2.96\n<\/code><\/pre>\npbase=ggplot(data_ex_wav,aes(x=x,y=y))\npY =pbase+\n geom_line(col=\"dark grey\")+ggtitle(\"Y\")\npDall=pbase +\n geom_line(aes(x=x,y=D4+D3+D2+D1),col=\"lightgreen\")+\n ggtitle(\"D4+...+D1\")\npDseason=pbase+\n geom_line(aes(x=x,y=D10+D9+D8+D7+D6+D5), col=\"slateblue\")+\n ggtitle(\"D10+...+D5\")\npTrend=pbase+\n geom_line(aes(x=x,y=S10),col=\"gold\",size=2)+\n ggtitle(\"S10\")\npY\/pTrend\/pDseason\/pDall\n<\/code><\/pre>\n![\"\"](\"..\/..\/lise.vaudor\/Rfigures\/Ondelettes_discretes\/show_graph-1.png\")
<\/p>\n","protected":false},"excerpt":{"rendered":"