INTRODUCTION

Une série temporelle, {yt, t = 1, 2, ··· , T}, est une suite d’observations d’un phénomène, indicées par une date, un temps. Le moment auquel l’observation est faite est généralement une information importante sur le phénomène observé.

Dans le cadre de notre étude statistique,nous serons amené à modéliser le trafic voyageur de la SNCF en 2e classe de 1963 à 1969. Plus précisement, il est question de réaliser une analyse descriptive de cette série temporelle dans un premier temps puis de prédire par la technique statistique de Holt Winters le nombre de voyageurs pour les 12 mois de l’année 1970.

Codes, packages et fonctions R

PARTIE I- ANALYSE DESCRIPTIVE DE SERIE TEMPORELLE

1- Construction et présentation de la série temporelle

# Construction de la série de données
a<-c(1750,1560,1820,2090,1910,2410,3140,2850,2090,1850,1630,2420,1710,1600,1800,2120,2100,2460,3200,2960,2190,1870,1770,2270,1670,1640,1770,2190,2020,2610,3190,2860,2140,1870,1760,2360,1810,1640,1860,1990,2110,2500,3030,2900,2160,1940,2750,2330,1850,1590,1880,2210,2110,2480,2880,2670,2100,1920,1670,2520,1834,1792,1860,2138,2115,2485,2581,2639,2038,1936,1784,2391,1798,1850,1981,2085,2120,2491,2834,2725,1932,2058,1856,2553)
# Affichage de la serie de données par mois de 1963 à 1969
serie_t <- ts(a, frequency = 12, start = c(1963, 1))
serie_t # Le graphique nous présente une série additive. 
##       Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
## 1963 1750 1560 1820 2090 1910 2410 3140 2850 2090 1850 1630 2420
## 1964 1710 1600 1800 2120 2100 2460 3200 2960 2190 1870 1770 2270
## 1965 1670 1640 1770 2190 2020 2610 3190 2860 2140 1870 1760 2360
## 1966 1810 1640 1860 1990 2110 2500 3030 2900 2160 1940 2750 2330
## 1967 1850 1590 1880 2210 2110 2480 2880 2670 2100 1920 1670 2520
## 1968 1834 1792 1860 2138 2115 2485 2581 2639 2038 1936 1784 2391
## 1969 1798 1850 1981 2085 2120 2491 2834 2725 1932 2058 1856 2553
# Structure de la serie temporelle
str(serie_t)
##  Time-Series [1:84] from 1963 to 1970: 1750 1560 1820 2090 1910 2410 3140 2850 2090 1850 ...
# Structure de la serie temporelle
attributes(serie_t)
## $tsp
## [1] 1963.000 1969.917   12.000
## 
## $class
## [1] "ts"
# Résumé statistique de la serie temporelle
summary(serie_t)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1560    1850    2095    2172    2481    3200

2- Visualisation graphique de la série temporelle

# Visualisation graphique en courbe de la serie temporelle
plot.ts(serie_t,xlab="Années",ylab="Passagers",las=1, main="Tendance de la série temporelle")

# Visualisation graphique en histogramme de la serie temporelle
hist(serie_t,xlab="Passagers",ylab="Fréquence",las=1, main="Fréquence de la série temporelle") #Les valeurs entre 1500 et 2000 sont les plus representatives dans la série de données. 

# Visualisation graphique sur une boîte à moustâche de la serie temporelle
boxplot(serie_t,ylab="Passagers",las=1, main="Boîte à moustâche de la série temporelle") #Les valeurs entre 1500 et 2000 sont les plus representatives dans la série de données. 

3- Indicateurs d’analyses de la série temporelle

A- Indicateurs de tendance centrale

# Moyenne
mean(serie_t) #Les valeurs entre 1500 et 2000 sont les plus representatives dans la série de données. 
## [1] 2171.738

B- Indicateurs de dispersion

# Variance
var(serie_t) #Les valeurs entre 1500 et 2000 sont les plus representatives dans la série de données. 
## [1] 178575.5
# Ecart-Type
sd(serie_t) #Les valeurs entre 1500 et 2000 sont les plus representatives dans la série de données. 
## [1] 422.5819

C-Indicateurs de forme

library(moments)
skewness(serie_t)
## [1] 0.7182608
kurtosis(serie_t) 
## [1] 2.585744

C-Indicateurs de dépendance

# Autocorellation
acf(serie_t,lag.max=12,plot = TRUE)

acf(serie_t,lag.max=12,plot = FALSE)
## 
## Autocorrelations of series 'serie_t', by lag
## 
## 0.0000 0.0833 0.1667 0.2500 0.3333 0.4167 0.5000 0.5833 0.6667 0.7500 0.8333 
##  1.000  0.495  0.031 -0.199 -0.178 -0.351 -0.458 -0.349 -0.129 -0.144  0.038 
## 0.9167 1.0000 
##  0.435  0.762
# Autocorellation partielle
pacf(serie_t,lag.max=12,plot = TRUE)

pacf(serie_t,lag.max=12,plot = FALSE)
## 
## Partial autocorrelations of series 'serie_t', by lag
## 
## 0.0833 0.1667 0.2500 0.3333 0.4167 0.5000 0.5833 0.6667 0.7500 0.8333 0.9167 
##  0.495 -0.285 -0.112  0.006 -0.416 -0.235 -0.125 -0.210 -0.483 -0.036  0.245 
## 1.0000 
##  0.366

PARTIE II- PREDICTIONS DU NOMBRE DE PASSAGERS POUR LES 12 PROCHAINS MOIS

1- Choix de la technique Holt-Winters pour la prédiction

Nous allons utiliser la technique de Holt-Winters. Cette technique permet de mettre à jour les prédictions en t+1 sur la base de moyennes pondérées des valeurs passées. Ce lissage exponentiel triple, décompose la variable d’intérêt en deux composants, une tendance (T) et une saisonnalité (S), auxquels on peut ajouter un terme d’erreur (E). Il existe deux grandes versions de ce modèle : une version additive qui est indiquée lorsque la variance est stable dans le temps, et une version multiplicative plus adéquate lorsque la variance croît/décroît dans le temps. Le lissage expo triple Hotwinter comporte les éléments suivants : (parametre gamma, parametre de la saisonalité c’est-à-dire les moments de pic et de flop et le phénonmene tendanciel).

La fonction Holt-Winters de notre série nous donne les meilleurs parametres alpha, beta et 12 saisonalités. Nous en déduisons que la série est saisonniere.

HoltWinters(serie_t)
## Holt-Winters exponential smoothing with trend and additive seasonal component.
## 
## Call:
## HoltWinters(x = serie_t)
## 
## Smoothing parameters:
##  alpha: 0.01771231
##  beta : 0.3099923
##  gamma: 0.3346422
## 
## Coefficients:
##            [,1]
## a   2161.165644
## b     -1.650786
## s1  -420.947398
## s2  -485.242941
## s3  -327.845171
## s4   -99.715999
## s5  -109.920570
## s6   276.199906
## s7   661.565776
## s8   539.981819
## s9  -159.893503
## s10 -239.975485
## s11 -328.243309
## s12  254.439558

Si les données prédites sont proches de la réalité alors nous pourrons utiliser cela pour faire nos prédictions. Mais pour cela il nous a fallu trouver un système de pondération qui me permet de coller aux valeurs réelles. Plus α , 𝜷, 𝜸 sont petits, plus on tient compte du passé lointain pour la prévision. Plus α , 𝜷, 𝜸 sont grands, plus on tient compte du passé récent pour la prévision.

# Utilisons les meilleures valeurs de alpha, beta et gamma produite par l'algorithme Holt-Winters
xlisse <- HoltWinters(serie_t,alpha=0.01771231,beta=0.3099923,gamma=0.3346422,seasonal="add")
xlisse
## Holt-Winters exponential smoothing with trend and additive seasonal component.
## 
## Call:
## HoltWinters(x = serie_t, alpha = 0.01771231, beta = 0.3099923,     gamma = 0.3346422, seasonal = "add")
## 
## Smoothing parameters:
##  alpha: 0.01771231
##  beta : 0.3099923
##  gamma: 0.3346422
## 
## Coefficients:
##            [,1]
## a   2161.165663
## b     -1.650788
## s1  -420.947417
## s2  -485.242966
## s3  -327.845193
## s4   -99.716016
## s5  -109.920588
## s6   276.199888
## s7   661.565771
## s8   539.981807
## s9  -159.893516
## s10 -239.975508
## s11 -328.243329
## s12  254.439535

HoltWinters nous a fait les prévisions pour la même période couverte par nos données d’origine.

xlisse$fitted
##              xhat    level      trend     season
## Jan 1964 1678.069 2111.439  5.6920163 -439.06250
## Feb 1964 1567.418 2117.697  5.8673397 -556.14583
## Mar 1964 1765.292 2124.141  6.0462343 -364.89583
## Apr 1964 2087.143 2130.802  6.2368059  -49.89583
## May 1964 2067.476 2137.621  6.4172110  -76.56250
## Jun 1964 2435.064 2144.614  6.5957904  283.85417
## Jul 1964 3174.322 2151.652  6.7327036 1015.93750
## Aug 1964 2891.651 2158.839  6.8736926  725.93750
## Sep 1964 2139.277 2166.924  7.2489770  -34.89583
## Oct 1964 1907.286 2175.071  7.5274813 -275.31250
## Nov 1964 1684.782 2181.938  7.3227550 -504.47917
## Dec 1964 2474.082 2190.770  7.7906607  275.52083
## Jan 1965 1773.050 2194.946  6.6701125 -428.56627
## Feb 1965 1660.460 2199.791  6.1042972 -545.43581
## Mar 1965 1858.038 2205.533  5.9919599 -353.48674
## Apr 1965 2176.379 2209.966  5.5085699  -39.09537
## May 1965 2155.428 2215.715  5.5833593  -65.87134
## Jun 1965 2515.791 2218.900  4.8397703  292.05086
## Jul 1965 3255.144 2225.409  5.3570431 1024.37820
## Aug 1965 2983.016 2229.612  4.9993595  748.40495
## Sep 1965 2218.534 2232.432  4.3239177  -18.22240
## Oct 1965 1951.689 2235.365  3.8927142 -287.56901
## Nov 1965 1764.788 2237.811  3.4441871 -476.46668
## Dec 1965 2453.024 2241.170  3.4178953  208.43607
## Jan 1966 1783.407 2242.941  2.9071289 -462.44039
## Feb 1966 1697.211 2246.319  3.0531411 -552.16119
## Mar 1966 1868.671 2248.358  2.7390160 -382.42623
## Apr 1966 2219.017 2250.944  2.6914049  -34.61790
## May 1966 2140.624 2249.579  1.4339436 -110.38838
## Jun 1966 2574.755 2250.470  1.2657947  323.01885
## Jul 1966 3254.232 2250.412  0.8553387 1002.96446
## Aug 1966 2954.888 2247.296 -0.3758469  707.96774
## Sep 1966 2201.233 2245.948 -0.6772175  -44.03760
## Oct 1966 1929.215 2244.540 -0.9036143 -314.42134
## Nov 1966 1764.942 2243.828 -0.8443986 -478.04071
## Dec 1966 2442.853 2260.431  4.5642366  177.85761
## Jan 1967 1813.242 2262.996  3.9445986 -453.69896
## Feb 1967 1700.771 2267.592  4.1464259 -570.96717
## Mar 1967 1888.038 2269.776  3.5382171 -385.27660
## Apr 1967 2166.767 2273.172  3.4940836 -109.89934
## May 1967 2160.708 2277.432  3.7314628 -120.45509
## Jun 1967 2582.164 2280.265  3.4530395  298.44575
## Jul 1967 3214.057 2281.909  2.8920893  929.25609
## Aug 1967 2969.867 2278.884  1.0578896  689.92535
## Sep 1967 2216.450 2274.630 -0.5885849  -57.59148
## Oct 1967 1959.875 2271.979 -1.2279768 -310.87623
## Nov 1967 2114.360 2270.045 -1.4469181 -154.23762
## Dec 1967 2397.602 2260.727 -3.8867592  140.76125
## Jan 1968 1814.178 2259.009 -3.2147107 -441.61600
## Feb 1968 1645.660 2256.145 -3.1058741 -607.37930
## Mar 1968 1865.410 2255.631 -2.3023673 -387.91878
## Apr 1968 2155.213 2253.233 -2.3320721  -95.68797
## May 1968 2111.046 2250.596 -2.4265829 -137.12368
## Jun 1968 2510.698 2248.239 -2.4048715  264.86290
## Jul 1968 3062.280 2245.379 -2.5459683  819.44659
## Aug 1968 2820.475 2234.309 -5.1885231  591.35454
## Sep 1968 2123.851 2225.906 -6.1849438  -95.87047
## Oct 1968 1887.560 2218.200 -6.6563222 -323.98376
## Nov 1968 1905.706 2212.402 -6.3903559 -300.30554
## Dec 1968 2377.793 2203.856 -7.0586059  180.99531
## Jan 1969 1754.945 2197.031 -6.9860894 -435.10020
## Feb 1969 1624.783 2190.808 -6.7496889 -559.27511
## Mar 1969 1792.837 2188.047 -5.5130952 -389.69714
## Apr 1969 2080.041 2185.867 -4.4799530 -101.34612
## May 1969 2041.198 2181.475 -4.4527246 -135.82387
## Jun 1969 2430.814 2178.418 -4.0200499  256.41574
## Jul 1969 2833.017 2175.464 -3.6895860  661.24265
## Aug 1969 2699.809 2171.792 -3.6841888  531.70102
## Sep 1969 2040.917 2168.554 -3.5458710 -124.09085
## Oct 1969 1850.874 2163.079 -4.1438998 -308.06095
## Nov 1969 1819.285 2162.604 -3.0066367 -340.31220
## Dec 1969 2342.779 2160.247 -2.8050447  185.33672
# Faisons la représentation graphique de notre lissage exponentiel
plot(xlisse)

2- Prédiction du nombre de passagers de Janvier à Décembre de l’année suivante

Nous allons prédire le nombre de voyageurs sans intervalle de confiance et avec intervalle de confiance. Ceci permettra de mieux nous ajuster en terme de prise de décisions. L’algorithme nous a permis de calculer les valeurs optimales des parametres.

A- Prédiction du nombre de passagers sans intervalle de confiance

xlisse <- HoltWinters(serie_t,alpha=0.01771231,beta=0.3099923,gamma=0.3346422,seasonal="add")
serie_t.prev<-predict(xlisse,n.ahead=12,prediction.interval=FALSE) 
serie_t.prev
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 1970 1738.567 1672.621 1828.368 2054.846 2042.991 2427.461 2811.176 2687.941
##           Sep      Oct      Nov      Dec
## 1970 1986.415 1904.682 1814.764 2395.796

B- Prédiction du nombre de passagers avec intervalle de confiance

xlisse <- HoltWinters(serie_t,alpha=0.01771231,beta=0.3099923,gamma=0.3346422,seasonal="add")
serie_t.prev<-predict(xlisse,n.ahead=12,prediction.interval=TRUE) 
serie_t.prev
##               fit      upr      lwr
## Jan 1970 1738.567 2083.281 1393.854
## Feb 1970 1672.621 2017.427 1327.815
## Mar 1970 1828.368 2173.316 1483.420
## Apr 1970 2054.846 2399.996 1709.697
## May 1970 2042.991 2388.411 1697.571
## Jun 1970 2427.461 2773.231 2081.690
## Jul 1970 2811.176 3157.387 2464.965
## Aug 1970 2687.941 3034.693 2341.189
## Sep 1970 1986.415 2333.817 1639.013
## Oct 1970 1904.682 2252.854 1556.510
## Nov 1970 1814.764 2163.834 1465.693
## Dec 1970 2395.796 2745.903 2045.688

C- Representation graphique de la prévision

# Représentation graphique de la prévision
par(mfrow=c(1,2), mar=c(3,3,3,3))
ts.plot(serie_t, serie_t.prev, col=c("black","red","green"))
library("forecast")
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
plot(forecast(xlisse))

par(mfrow=c(1,1) , mar=c(0,0,0,0))

3- Validation du modèle de lissage par l’analyse des résidus

Il s’agit de vérifier si les résidus du modèle choisi suivent un processus bruit blanc Dans le cas contraire, le modèle n’est pas bon. Il va falloir choisir un autre ou reprendre les étapes précédentes. La satisfaction de la seule condition de bruit blanc suffit pour valider le modèle. Mais si les résidus suivent en plus un processus gaussien (loi normal), le modèle serait encore plus intéressant.

Nous allons nous interesser à l’analyse des erreurs de prévision dans l’échantillon de données.

# faire le lissage exponentielle 
xlisse <- HoltWinters(serie_t)
# récupérer les résidus 
residus=residuals(xlisse)
residus
##               Jan          Feb          Mar          Apr          May
## 1964   31.9310897   32.5815098   34.7081814   32.8566131   32.5241019
## 1965 -103.0501318  -20.4596420  -88.0382830   13.6211473 -135.4275057
## 1966   26.5927513  -57.2106051   -8.6712478 -229.0173989  -30.6244264
## 1967   36.7581691 -110.7711115   -8.0378878   43.2331339  -50.7083350
## 1968   19.8220788  146.3401552   -5.4100165  -17.2129344    3.9542472
## 1969   43.0548925  225.2168951  188.1629093    4.9590424   78.8016792
##               Jun          Jul          Aug          Sep          Oct
## 1964   24.9355674   25.6778635   68.3493562   50.7230867  -37.2861517
## 1965   94.2092571  -65.1437868 -123.0160555  -78.5337220  -81.6888035
## 1966  -74.7550262 -224.2318904  -54.8876538  -41.2329033   10.7847831
## 1967 -102.1640616 -334.0569227 -299.8671485 -116.4503981  -39.8750646
## 1968  -25.6975024 -481.2800577 -181.4748967  -85.8506130   48.4396175
## 1969   60.1863588    0.9830088   25.1914038 -108.9170565  207.1261049
##               Nov          Dec
## 1964   85.2181836 -204.0818890
## 1965   -4.7884300  -93.0242552
## 1966  985.0575202 -112.8526770
## 1967 -444.3605065  122.3980770
## 1968 -121.7062373   13.2072267
## 1969   36.7153074  210.2211246
#Graphe des résidus : Vérification de la blancheur des résidus 
plot(residus)

# Autocorrelation  
acf(residus, lag.max=12, na.action = na.pass)

Box.test(residus, lag=20, type="Ljung-Box")
## 
##  Box-Ljung test
## 
## data:  residus
## X-squared = 20.89, df = 20, p-value = 0.4036
hist(residus)

shapiro.test(residus)
## 
##  Shapiro-Wilk normality test
## 
## data:  residus
## W = 0.77549, p-value = 3.97e-09
# La moyenne des erreurs de prédiction
mean(residus, na.rm=TRUE)
## [1] -18.5739

CONCLUSION