4.3 Manipulation et description

Tout comme pour les statistiques plus classiques, nous commençons toujours notre exploration de séries chronologiques par une partie purement descriptive, soit numérique, soit graphique. De plus, nous savons que nous ne devons pas espérer obtenir d’emblée des données qui se présentent de manière parfaite. Un petit peu de manipulation de données est parfois nécessaire… et nous allons découvrir quelques fonctions utiles dans ce contexte pour les séries chronologiques.

Comme nous l’avons déjà vu, le graphe de base (obtenu par plot() ou en utilisant geom_line()) est une ligne brisée qui relie les différents points de la série, avec le temps sur l’axe des abscisses, et l’axe des ordonnées représentant les observations de la variable quantitative mesurée au fil du temps.

4.3.1 Statistiques glissantes

Les descripteurs statistiques tels que la moyenne, la médiane, l’écart type, etc. restent utiles pour les séries temporelles. Par contre, la valeur de ces descripteurs est susceptible de varier au cours du temps, et cette variation donne des indications utiles sur les mécanismes biologiques sous-jacents (par exemple, des fluctuations, des changements brutaux de régimes, etc.). Il nous faut donc calculer ces descripteurs pour des intervalles de temps précis le long de l’axe temporel de notre série. Les statistiques glissantes (sliding statistics en anglais) correspondent précisément à l’analyse de blocs successifs de données suivant un axe spatio-temporel. Dans {pastecs}, elles se calculent à l’aide de la fonction stat.slide().

Exemple

Nous reprenons ici l’exemple du manuel de {pastecs}. On peut calculer des statistiques glissantes pour la série ClausocalanusA de marbio par groupe de 10 stations, les imprimer et faire un graphique de la “moyenne glissante” à l’aide des instructions suivantes :

marbio_ts <- ts(read("marbio", package = "pastecs"))
statsl <- stat.slide(1:68, marbio_ts[, "ClausocalanusA"],
  xmin = 0, n = 7, deltat = 10)
statsl
# # A data.frame: 10 x 7
#    ` `      `[0,10[` `[10,20[` `[20,30[` `[30,40[` `[40,50[` `[50,60[` `[60,70[`
#    <rownam>    <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>
#  1 xmin         0          10        20        30        40        50        60 
#  2 xmax        10          20        30        40        50        60        70 
#  3 nbr.val      9          10        10        10        10        10         9 
#  4 nbr.null     2           0         0         0         0         0         0 
#  5 nbr.na       0           0         0         0         0         0         0 
#  6 min          0         160       158        96       136        56        52 
#  7 max         12         832      1980      1204      2484       824       656 
#  8 median       4         344       732       752       260       340       120 
#  9 mean         4.78      351.      810.      682       494.      365.      219.
# 10 std.dev      4.79      197.      620.      350.      711.      235.      202.
plot(statsl, stat = "mean", leg = TRUE, lpos = c(55, 2500),
  xlab = "Station", ylab = "ClausocalanusA")

Toujours pour la même série, on peut calculer toutes les statistiques sur des intervalles irréguliers, puis représenter l’étendue (minimum, maximum) et la médiane pour chaque intervalle comme suit (ces intervalles correspondent aux différentes masses d’eaux croisées le long du transect, voir légende du graphique) :

statsl2 <- stat.slide(1:68, marbio_ts[, "ClausocalanusA"],
  xcut = c(0, 17, 25, 30, 41, 46, 70),
  basic = TRUE, desc = TRUE, norm = TRUE, pen = TRUE, p = 0.95)
statsl2
# # A data.frame: 29 x 6
#    ` `        `[0,17[` `[17,25[` `[25,30[` `[30,41[` `[41,46[` `[46,70[`
#    <rownames>    <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>
#  1 xmin              0        17        25        30        41        46
#  2 xmax             17        25        30        41        46        70
#  3 nbr.val          16         8         5        11         5        23
#  4 nbr.null          2         0         0         0         0         0
#  5 nbr.na            0         0         0         0         0         0
#  6 min               0       158      1040        96       136        52
#  7 max             832       424      1980      1204       264      2484
#  8 range           832       266       940      1108       128      2432
#  9 sum            2785      2151      6716      7060      1088      9236
# 10 median           12       212      1204       744       256       308
# # … with 19 more rows
plot(statsl2, stat = "median", xlab = "Stations", ylab = "Effectifs",
  main = "Clausocalanus A")     # Médiane
lines(statsl2, stat = "min")    # Minimum
lines(statsl2, stat = "max")    # Maximum
lines(c(17, 17), c(-50,2600), col = 4, lty = 2) # Séparations des masses d'eaux
lines(c(25, 25), c(-50,2600), col = 4, lty = 2)
lines(c(30, 30), c(-50,2600), col = 4, lty = 2)
lines(c(41, 41), c(-50,2600), col = 4, lty = 2)
lines(c(46, 46), c(-50,2600), col = 4, lty = 2)
text(c(8.4, 21, 27.5, 35, 43.5, 57.2), 2300,
  labels = c("Zone périphérique", "D1", "C", "Front", "D2", "Zone centrale")) # Labels
legend(0, 1900, c("série", "médian", "étendue"), col = 1:3, lty = 1)

Enfin, on peut extraire différentes statistiques, les valeurs de la série initiale (y), les valeurs de temps de la série (x), et le vecteur temps de coupure des périodes (xcut) à partir de l’objet stat.slide renvoyé :

statsl2$stat[c("mean", "median", "min", "max"), ]
# # A data.frame: 4 x 6
#   ` `        `[0,17[` `[17,25[` `[25,30[` `[30,41[` `[41,46[` `[46,70[`
#   <rownames>    <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>
# 1 mean           174.      269.     1343.      642.      218.      402.
# 2 median          12       212      1204       744       256       308 
# 3 min              0       158      1040        96       136        52 
# 4 max            832       424      1980      1204       264      2484
statsl2$y
#  [1]    0    4    2    8    4    0    1   12   12  196  160  406  448  328  372
# [16]  832  360  200  206  218  158  193  424  392 1980 1348 1144 1204 1040 1204
# [31] 1080  424  744   96  888  608  776  240  760  240  256  264  136  264  168
# [46]  516  136  480 2484  824  376  352  160  640  280  488  144  328   56  112
# [61]  408  308   68   52   84  164  656  120
statsl2$x
#  [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
# [26] 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
# [51] 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68
statsl2$xcut
# [1]  0 17 25 30 41 46 70
À vous de jouer !

Effectuez maintenant les exercices du tutoriel C04La_ts_intro (Introduction aux séries chronologiques).

BioDataScience3::run("C04La_ts_intro")

4.3.2 Manipulations de ts

Les jeux de données bricolés artificiellement sont parfois bien utiles pour explorer de nouvelles techniques. En effet, puisque nous les construisons nous-mêmes, nous savons dès le départ de quoi ils sont faits. Et donc, nous pouvons vérifier que le résultat obtenu est bien conforme. Par exemple, nous pouvons construire une série de données mensuelle commençant en avril 1998 (start = c(1998, 4)), contenant 100 valeurs et comportant une composante sinusoïdale d’une période annuelle ainsi qu’une composante aléatoire ayant une distribution normale de moyenne nulle et d’écart type 0.5 :

tser <- ts(sin((1:100) / 6 * pi) + rnorm(100, sd = 0.5), start = c(1998, 4), 
  frequency = 12)
tser
#              Jan         Feb         Mar         Apr         May         Jun
# 1998                                      0.43284443  0.55796734  1.70265599
# 1999 -1.24312825 -0.49738897  1.19660403  0.05225123  0.57422522  2.21443653
# 2000 -1.27975774 -0.39864209  0.19718918 -0.25121656  0.79224043  0.67165187
# 2001 -1.33085705 -0.31227960 -0.34076874  1.26087821  0.90076680  0.72463865
# 2002 -1.49102561 -1.23752020  0.36432449  0.16888829  0.49321406  0.41112583
# 2003 -0.33160770  0.03132431  0.75403988  0.54097729  0.91892553  1.01564832
# 2004 -0.66010261 -0.55657839  0.90346085  0.88632106  0.28931583  0.52236703
# 2005 -0.18028826 -1.67848525  0.75703973  0.37294084  1.04893710  1.20347449
# 2006 -0.72877221 -0.71386908 -0.33470161  0.97725885  0.35056006  2.09328093
#              Jul         Aug         Sep         Oct         Nov         Dec
# 1998  0.71754359  1.65962493 -0.34002388 -0.93349797 -0.08686532 -1.06923582
# 1999  1.30384411  1.35349511 -0.14072954 -0.58534641 -1.44827029 -1.04973955
# 2000  1.27543290  0.27594795 -0.10605802 -0.42613394 -1.26309700  0.04021914
# 2001  0.43546515  0.99636669 -0.31451471 -0.36733500 -1.35331820 -1.11501552
# 2002  0.81847016  0.36884872  0.51316559 -0.83152133 -0.71954553 -0.82510491
# 2003  1.63474275  1.08403393 -0.04882395 -0.01570617 -0.54439103 -1.40656027
# 2004  1.66861579  0.30669569  0.53616651  0.03691720 -0.86264360 -1.26470062
# 2005  0.86734660 -0.04116614 -0.72024479 -0.40289493 -1.21954480 -1.17567484
# 2006  0.89088635

Créons maintenant une série multiple en prenant notre série de départ, ainsi qu’une version décalée de cinq mois à l’aide de la fonction stats::lag() et représentons cela graphiquement :

mtser <- ts.intersect(tser, stats::lag(tser, 5))
plot(mtser, xlab = "Temps", main = "Série tser et série décalée de 5 mois")

4.3.2.1 Manipulation du temps

Les fonctions tsp(), start(), end() et frequency() permettent de récupérer les différents paramètres de temps de la série :

tsp(tser)
# [1] 1998.25 2006.50   12.00
start(tser)
# [1] 1998    4
end(tser)
# [1] 2006    7
frequency(tser)
# [1] 12

Nous avons déjà vu la fonction time() qui permet de reconstituer le vecteur temps pour la série. La fonction cycle() indique l’appartenance de chaque donnée de la série à un cycle. Elle permet, par exemple, de séparer les données par mois si la base de temps est l’année à l’aide de la fonction split(). Ensuite, il est possible de traiter ou de représenter séparément les statistiques mois par mois pour en faire, par exemple, des boites de dispersions parallèles :

tser.cycle <- cycle(tser)
tser.cycle
#      Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
# 1998               4   5   6   7   8   9  10  11  12
# 1999   1   2   3   4   5   6   7   8   9  10  11  12
# 2000   1   2   3   4   5   6   7   8   9  10  11  12
# 2001   1   2   3   4   5   6   7   8   9  10  11  12
# 2002   1   2   3   4   5   6   7   8   9  10  11  12
# 2003   1   2   3   4   5   6   7   8   9  10  11  12
# 2004   1   2   3   4   5   6   7   8   9  10  11  12
# 2005   1   2   3   4   5   6   7   8   9  10  11  12
# 2006   1   2   3   4   5   6   7
boxplot(split(tser, tser.cycle), names = month.abb, col = "cornsilk")

Une autre manipulation courante de l’axe temporel consiste à agréger les données en des pas de temps plus longs. aggregate() permet de réduire le nombre d’observations en diminuant la fréquence. Par exemple pour transformer la série tser qui a un pas de temps d’un mois en une série de pas de temps d’un trimestre, en calculant des moyennes trimestrielles, nous ferons (puisqu’il y a quatre trimestres dans une année, nous indiquons 4 ici) :

aggregate(tser, 4, mean)
#             Qtr1        Qtr2        Qtr3        Qtr4
# 1998              0.89782259  0.67904821 -0.69653303
# 1999 -0.18130440  0.94697099  0.83886989 -1.02778542
# 2000 -0.49373688  0.40422525  0.48177428 -0.54967060
# 2001 -0.66130180  0.96209455  0.37243904 -0.94522291
# 2002 -0.78807377  0.35774273  0.56682815 -0.79205726
# 2003  0.15125216  0.82518371  0.88998425 -0.65555249
# 2004 -0.10440672  0.56600131  0.83715933 -0.69680901
# 2005 -0.36724459  0.87511748  0.03531189 -0.93270486
# 2006 -0.59244763  1.14036661

Notez encore une fois que R assume que l’unité de temps est l’année et il crée des intitulés particuliers pour des séries de fréquence égale à 4 (Qtr1 -> Qtr4) ou à 12 (intitulé des mois en abrégé).

Si nous ne souhaitons pas utiliser la série sur toute sa longueur, nous pouvons employer window() pour extraire une sous-série contenue dans une fenêtre de temps spécifique sans modifier la fréquence des observations. Par exemple, pour extraire les années 1999 à 2001 complètes de tser, on utilisera :

window(tser, start = c(1999, 1), end = c(2001, 12))
#              Jan         Feb         Mar         Apr         May         Jun
# 1999 -1.24312825 -0.49738897  1.19660403  0.05225123  0.57422522  2.21443653
# 2000 -1.27975774 -0.39864209  0.19718918 -0.25121656  0.79224043  0.67165187
# 2001 -1.33085705 -0.31227960 -0.34076874  1.26087821  0.90076680  0.72463865
#              Jul         Aug         Sep         Oct         Nov         Dec
# 1999  1.30384411  1.35349511 -0.14072954 -0.58534641 -1.44827029 -1.04973955
# 2000  1.27543290  0.27594795 -0.10605802 -0.42613394 -1.26309700  0.04021914
# 2001  0.43546515  0.99636669 -0.31451471 -0.36733500 -1.35331820 -1.11501552

La fonction stats::lag() permet de décaler la série en arrière dans le temps (ou en avant si on fournit une valeur de décalage négative) ; nous l’avons déjà utilisée plus haut. La fonction diff() calcule la différence entre une série et elle-même décalée dans le temps de k observations. Les fonction ts.intersect() et ts.union() permettent de réunir deux ou plusieurs séries au sein d’une même matrice multivariée, avec une échelle de temps commune. Les fréquences respectives des séries à assembler doivent être identiques. Là où ts.intersect() retient uniquement l’intervalle de temps commun à toutes les séries, ts.union() conservera l’ensemble de l’intervalle temporel, toutes séries confondues.

Exercez-vous !

Manipulez et décrivez les différentes séries d’exemples. Notez par exemple que EEG est constituée de mesures effectuées à un rythme de 256 par seconde, mais que l’axe du temps reprend simplement les observations les unes après les autres (ou si vous préférez, une “unité” de temps de 1/256 sec et une fréquence de 1). Comment feriez-vous pour transformer cette série avec une unité de temps de la seconde ? Ou de la milliseconde ?