D2
Gaby
2025-05-18
#Commande pour effacer tout
#rm(list=ls())
setwd("~/USHERBROOKE/ESCUELA/ÉTÉ/MQG813 - Statistiques décisionnelles avancées/D2")
#Verifie si une librairie est installee, l'installe au besoin et la charge
list.of.packages <- c("modelr", "rio", "tidyverse", "DescTools", "pastecs",
"psych", "pander", "gmodels", "vcd", "fastDummies",
"questionr", "PerformanceAnalytics", "Hmisc", "FSA",
"car", "robustHD", "hexView", "ggplot2", "dplyr",
"forecast", "latticeExtra", "numbers", "zoo", "lubridate",
"stringr", "fastDummies", "descr","tinytex","readr")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)) {install.packages(new.packages)}
lapply(list.of.packages, require, character.only = TRUE)## Loading required package: modelr## Loading required package: rio## Loading required package: tidyverse## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.0.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
## Loading required package: DescTools
##
## Loading required package: pastecs
##
##
## Attaching package: 'pastecs'
##
##
## The following objects are masked from 'package:dplyr':
##
## first, last
##
##
## The following object is masked from 'package:tidyr':
##
## extract
##
##
## Loading required package: psych
##
##
## Attaching package: 'psych'
##
##
## The following objects are masked from 'package:DescTools':
##
## AUC, ICC, SD
##
##
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
##
##
## Loading required package: pander
##
## Loading required package: gmodels
##
## Registered S3 method overwritten by 'gdata':
## method from
## reorder.factor DescTools
##
## Loading required package: vcd
##
## Loading required package: grid
##
## Loading required package: fastDummies
##
## Loading required package: questionr
##
##
## Attaching package: 'questionr'
##
##
## The following object is masked from 'package:psych':
##
## describe
##
##
## Loading required package: PerformanceAnalytics
##
## Loading required package: xts
##
## Loading required package: zoo
##
##
## Attaching package: 'zoo'
##
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
##
##
## ######################### Warning from 'xts' package ##########################
## # #
## # The dplyr lag() function breaks how base R's lag() function is supposed to #
## # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or #
## # source() into this session won't work correctly. #
## # #
## # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
## # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop #
## # dplyr from breaking base R's lag() function. #
## # #
## # Code in packages is not affected. It's protected by R's namespace mechanism #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning. #
## # #
## ###############################################################################
##
##
## Attaching package: 'xts'
##
##
## The following objects are masked from 'package:pastecs':
##
## first, last
##
##
## The following objects are masked from 'package:dplyr':
##
## first, last
##
##
##
## Attaching package: 'PerformanceAnalytics'
##
##
## The following object is masked from 'package:vcd':
##
## Kappa
##
##
## The following object is masked from 'package:graphics':
##
## legend
##
##
## Loading required package: Hmisc
##
##
## Attaching package: 'Hmisc'
##
##
## The following objects are masked from 'package:questionr':
##
## describe, wtd.mean, wtd.table, wtd.var
##
##
## The following object is masked from 'package:psych':
##
## describe
##
##
## The following objects are masked from 'package:DescTools':
##
## %nin%, Label, Mean, Quantile
##
##
## The following objects are masked from 'package:dplyr':
##
## src, summarize
##
##
## The following objects are masked from 'package:base':
##
## format.pval, units
##
##
## Loading required package: FSA
##
## ## FSA v0.9.6. See citation('FSA') if used in publication.
## ## Run fishR() for related website and fishR('IFAR') for related book.
##
##
## Attaching package: 'FSA'
##
##
## The following object is masked from 'package:psych':
##
## headtail
##
##
## Loading required package: car
##
## Loading required package: carData
##
## Registered S3 methods overwritten by 'car':
## method from
## hist.boot FSA
## confint.boot FSA
##
##
## Attaching package: 'car'
##
##
## The following object is masked from 'package:FSA':
##
## bootCase
##
##
## The following object is masked from 'package:psych':
##
## logit
##
##
## The following object is masked from 'package:DescTools':
##
## Recode
##
##
## The following object is masked from 'package:dplyr':
##
## recode
##
##
## The following object is masked from 'package:purrr':
##
## some
##
##
## Loading required package: robustHD
##
## Loading required package: perry
##
## Loading required package: parallel
##
##
## Attaching package: 'perry'
##
##
## The following object is masked from 'package:modelr':
##
## mape
##
##
## Loading required package: robustbase
##
## Loading required package: hexView
##
## Loading required package: forecast
##
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
##
##
## Attaching package: 'forecast'
##
##
## The following object is masked from 'package:DescTools':
##
## BoxCox
##
##
## Loading required package: latticeExtra
##
## Loading required package: lattice
##
##
## Attaching package: 'latticeExtra'
##
##
## The following object is masked from 'package:vcd':
##
## rootogram
##
##
## The following object is masked from 'package:ggplot2':
##
## layer
##
##
## Loading required package: numbers
##
##
## Attaching package: 'numbers'
##
##
## The following object is masked from 'package:PerformanceAnalytics':
##
## Omega
##
##
## The following object is masked from 'package:psych':
##
## omega
##
##
## The following objects are masked from 'package:DescTools':
##
## GCD, LCM, Primes
##
##
## Loading required package: descr
##
##
## Attaching package: 'descr'
##
##
## The following object is masked from 'package:questionr':
##
## freq
##
##
## The following object is masked from 'package:gmodels':
##
## CrossTable
##
##
## Loading required package: tinytex## [[1]]
## [1] TRUE
##
## [[2]]
## [1] TRUE
##
## [[3]]
## [1] TRUE
##
## [[4]]
## [1] TRUE
##
## [[5]]
## [1] TRUE
##
## [[6]]
## [1] TRUE
##
## [[7]]
## [1] TRUE
##
## [[8]]
## [1] TRUE
##
## [[9]]
## [1] TRUE
##
## [[10]]
## [1] TRUE
##
## [[11]]
## [1] TRUE
##
## [[12]]
## [1] TRUE
##
## [[13]]
## [1] TRUE
##
## [[14]]
## [1] TRUE
##
## [[15]]
## [1] TRUE
##
## [[16]]
## [1] TRUE
##
## [[17]]
## [1] TRUE
##
## [[18]]
## [1] TRUE
##
## [[19]]
## [1] TRUE
##
## [[20]]
## [1] TRUE
##
## [[21]]
## [1] TRUE
##
## [[22]]
## [1] TRUE
##
## [[23]]
## [1] TRUE
##
## [[24]]
## [1] TRUE
##
## [[25]]
## [1] TRUE
##
## [[26]]
## [1] TRUE
##
## [[27]]
## [1] TRUE
##
## [[28]]
## [1] TRUE
##
## [[29]]
## [1] TRUEinstall_formats()
#Chargement des librairies et fonctions necessaires
source("MQG813_librairie_temporel.R")## Loading required package: nlme
##
## Attaching package: 'nlme'
##
## The following object is masked from 'package:forecast':
##
## getResponse
##
## The following object is masked from 'package:dplyr':
##
## collapse
##
## Loading required package: knitr
## Loading required package: ggpubr
##
## Attaching package: 'ggpubr'
##
## The following object is masked from 'package:forecast':
##
## gghistogram
##
## Loading required package: vars
## Loading required package: MASS
##
## Attaching package: 'MASS'
##
## The following object is masked from 'package:dplyr':
##
## select
##
## Loading required package: strucchange
## Loading required package: sandwich
##
## Attaching package: 'strucchange'
##
## The following object is masked from 'package:stringr':
##
## boundary
##
## Loading required package: urca
## Loading required package: lmtest
##
## Attaching package: 'vars'
##
## The following object is masked from 'package:DescTools':
##
## PhiMISE EN SITUATION Le jeu de donn´ees Demande_´Electricité_2019.csv
contient l’historique dela demande quotidienne d’´electricité (en
Gigawattheures) au Québec entre avril 2019 et juillet 2019 . Le jeu de
données Température_2019-2022.csv contient quant
a lui les donn´ees quotidiennes de températurea Montréal,
qui serviront de référence pour approximer la température au Québec. On
souhaite d´evelopper un mod`ele de r´egression lin´eaire afin d’´etablir
des pr´evisions sur la demande quotidienne d’électricité. Les donn´ees
du mois de juillet serviront de plage d’essai.
On utilisera un seuil de signification de α = 0, 05.
library(tidyverse)
library(lubridate)
library(rio)
dem<- import("Demande_Électricité_2019.csv")
temp<-import("temperature_2019-2022.csv")
#ajout les variable t, mois et annee
#Création de la variable mois et de la variable année
dem <- dem%>% mutate(t = row_number(),
mois = as.factor(month(DATE, label=T, abbr=F)),
année = as.factor(year(DATE)),jour_semaine = wday(DATE, label = TRUE, abbr = FALSE, week_start = 1))QUESTION 1 (3 points) Pr´esentez et commentez le graphe s´equentiel et les graphes saisonniers des demandes quotidiennes entre avril 2019 et juin 2019.
Graphe sequentiel
commentaire: Le graphique séquentiel montre une tendance générale à la baisse de la demande entre avril et juin 2019. Cette tendance est non linéaire, comme le suggère la courbe de lissage LOESS (en rouge), qui permet de mieux visualiser les variations progressives dans le temps. On observe des fluctuations(patterns), ce qui laisse supposer une saisonnalité hebdomadaire ou des effets spécifiques à certains jours (on supponse que pour les weekend, on voit des valeurs plus bas que pour le rest de la semaine)
dem<-dem%>% mutate(DATE=as.Date(DATE))
dem<-dem%>% mutate(t=as.numeric(t))
#dem_03_06<-dem %>% filter(DATE<"2019-07-01")
app<-c(1:91)
essai<-c(92:122)
#GRAPHIQUE
ggplot(data = dem[app,] , aes(x=DATE, y=DEMANDE)) +
geom_line() +
geom_point() +
geom_smooth(method="loess", se=FALSE, color="red") +
labs(title="Demande quitidiennes d'électricité entre avril et juin 2019",
x="Date",
y="Demande (GWh)") +
theme_minimal() +
scale_x_date(date_labels = "%b %Y", date_breaks = "1 month") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))## `geom_smooth()` using formula = 'y ~ x'
Graphe saisonnier
commentaire: Le graphique saisonnie est divisé par jour de la semaine: La demande est plus élevée les lundis et mardis, probablement en lien avec la reprise des activités économiques. Elle diminue généralement en fin de semaine, en particulier les samedis et dimanches, ce qui reflète un comportement de travail typique, dans les weekends(Moins de consommation industrielle) . La ligne rouge horizontale, qui représente la moyenne hebdomadaire, permet de visualiser les écarts entre jours ouvrables et week-ends.
#moyenne mensuel
dem$moy<-NA
dem[app,]<-dem[app,]%>% group_by(jour_semaine) %>%
mutate(moy=mean(DEMANDE))%>% ungroup()
#Graphique sequentiel saisonnier
ggplot(data = dem[app,], aes (x=DATE, y=DEMANDE)) +
geom_line() +
facet_wrap(~ jour_semaine, ncol=7) +
geom_point() +
geom_hline(aes(yintercept = moy,group=mois), color = "red") + labs(title="Demande quotidienne d'électricité entre avril et juin 2019",
x="Date",
y="Demande (GWh)") +
theme_bw()
QUESTION 2 (3 points) Premier
modele : On ´etudiera d’abord le modéle dans lequel seule la variable incrémentale t servira de variable explicative. (a) Pr´esentez un r´esum´e de ce premier modele
`a l’aide de la fonction summary().
Le modele explique 72.4% de la variation de la demande Le modele me dit :Ma demande debase d’ electricité est 526.6 et on enleve 1.823GWh de la demande a chaque fois que t augmente (donc a chaque periode qu’on augmente). Aussi on voit que la variable t est tres significatif (p<0.0001) .
################
### Modele 1 ###
################
#Premier modele de regression lineaire
modele1 <- lm(DEMANDE ~ t,
data = dem[app,])
#Resume des statistiques du modele (2 possibilites)
pander(summary(modele1))| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 526.6 | 6.316 | 83.37 | 2.89e-86 |
| t | -1.823 | 0.1192 | -15.29 | 1.243e-26 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 91 | 29.88 | 0.7243 | 0.7212 |
#Intervalles de confiance des coefficients de la regression
#pander(confint(modele1, level = 0.95))- Pr´esentez le corrélogramme des r´esidus du
mod
ele. Quelle variable seriez-vous tent´e d’ajouter au modele ?
On choisit d’ajouter le décalage (lag) 1 en premier, car c’est celui qui présente la plus forte autocorrélation dans le corrélogramme. Certes, on observe aussi des autocorrélations notables aux décalages 2 et 7, et même des pics plus faibles aux décalages 3, 5, 6 et 8. Cependant, commencer par lag(1) est une approche logique, en espérant que cela contribuera à réduire également les autres autocorrélations résiduelles.
#Pour commencer on doit calculer les previsions,residus et residus standardises
dem<-ajouter_previsions(dem,modele1,plage_connue=app,
intervalles=95,residus=T,
noms_personnalises="1")
#Corrélogramme des residus
correlogramme(dem[app,"RES_1"], lag.max = 25)
## lag ac pac
## 1 1 0.68395311 0.683953113
## 2 2 0.32877570 -0.261206376
## 3 3 0.23931686 0.270955614
## 4 4 0.22641285 -0.042568153
## 5 5 0.16598712 0.005455624
## 6 6 0.23901777 0.297140067
## 7 7 0.34183649 0.028511088
## 8 8 0.24054671 -0.155980244
## 9 9 0.09571062 0.044133016
## 10 10 0.03895944 -0.106360376
## 11 11 0.03472555 0.033201802
## 12 12 0.03266384 0.034809746
## 13 13 0.11210865 0.085482866
## 14 14 0.21785472 0.094166997
## 15 15 0.19396155 -0.046615279
## 16 16 0.10632253 0.021735122
## 17 17 -0.01014018 -0.188503893
## 18 18 -0.03592494 0.112576330
## 19 19 0.01260555 0.027874420
## 20 20 0.06975496 -0.049644665
## 21 21 0.08456308 0.017926875
## 22 22 0.07090315 -0.027639333
## 23 23 0.05790869 0.046312587
## 24 24 -0.02902867 -0.089376910
## 25 25 -0.12377066 -0.093911573#HISTOGRAME des residus
#ggplot(dem[app,], aes(x=RES_1)) +
# geom_histogram(bins=10)
#Test de normalite des residus (Shapiro-Wilk)
#pander(shapiro.test(dem[["RES_1"]][app]))
#Multicolinearite (plus d'une variable explicative)
#pander(vif(modele1)) # fonctionne pas on a petite nombre de variables
#Indices de performance
#indices_performance(prev$DEMANDE, prev$PREV_1)QUESTION 3 (10 points) Deuxieme modele : On tentera
d’ajouter des variables pertinentes au premier
modele en basant nos choix sur le corrélogramme des résidus. (a) Ajoutez la variable explicative proposéea
la question 2 (lag1) (b). Pr´esentez le corrélogramme des nouveaux
résidus et ajoutez une autre variable explicative si besoin est. Répétez
ces étapes jusqu’`a ce que le corrélogramme ne présente plus aucune
autocorr´elation significative.
#Deuxieme modele de regression lineaire, on ajoute lag1
model2<-lm(DEMANDE~ t+lag(DEMANDE,1),data=dem[app,])
#Ajouter les prévisions, les residus et les residus standardises
dem <- ajouter_previsions(dem, model2,residus=T,noms_personnalises="2")
#Correlogramme des residus
correlogramme(dem[app,"RES_2"], lag.max = 25)## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).
## lag ac pac
## 1 1 0.193957020 0.193957020
## 2 2 -0.288371592 -0.338733857
## 3 3 -0.137458064 0.003661452
## 4 4 0.077459349 0.017203234
## 5 5 -0.100458844 -0.199751714
## 6 6 0.056596544 0.193791964
## 7 7 0.243717177 0.132589248
## 8 8 0.043382646 -0.044224043
## 9 9 -0.112905101 0.071032356
## 10 10 -0.061845431 -0.059063862
## 11 11 0.006948458 -0.010347223
## 12 12 -0.021218968 0.010346087
## 13 13 -0.012178846 -0.070392614
## 14 14 0.140158310 0.143873050
## 15 15 0.112255240 0.037258271
## 16 16 0.092521409 0.167553432
## 17 17 -0.145398407 -0.153343071
## 18 18 -0.096793947 0.031232421
## 19 19 -0.030114080 -0.043489356
## 20 20 0.148958328 0.111137092
## 21 21 0.102452417 0.031131228
## 22 22 0.013261010 -0.040091627
## 23 23 0.053855900 0.145968097
## 24 24 0.001148652 -0.014574442
## 25 25 -0.149980891 -0.1011841563EME ESSAIS ON FAIT LAG 1 ET 2
#Troiseme modele de regression lineaire, on ajoute lag1+lag2
model3<-lm(DEMANDE~ t+lag(DEMANDE,1) +lag(DEMANDE,2),data=dem[app,])
#Ajouter les prévisions, les residus et les residus standardises
dem <- ajouter_previsions(dem, model3,residus=T,noms_personnalises="3")
#Graphe sequentiel des residus standardises
#ggplot(data=dem[app,], aes(x=t, y = ZRES_3)) +
# geom_point() +
#geom_line() +
#geom_hline(yintercept = 0)
#Correlogramme des residus
correlogramme(dem[app,"RES_3"], lag.max = 25)## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).
## lag ac pac
## 1 1 0.081932694 0.08193269
## 2 2 -0.174439692 -0.18237695
## 3 3 0.023891729 0.05804043
## 4 4 0.228026618 0.19634943
## 5 5 -0.081604658 -0.11905394
## 6 6 0.053579531 0.15690292
## 7 7 0.242121867 0.19671210
## 8 8 0.027914298 -0.04144687
## 9 9 -0.078371499 0.04039571
## 10 10 -0.008674383 -0.06222334
## 11 11 0.060567929 -0.02890939
## 12 12 0.011971475 0.04307279
## 13 13 -0.028865022 -0.07668207
## 14 14 0.127621448 0.12690013
## 15 15 0.051399126 0.01520248
## 16 16 0.125542180 0.18408027
## 17 17 -0.128138305 -0.13943231
## 18 18 -0.022633797 -0.01985218
## 19 19 -0.034447347 -0.07296971
## 20 20 0.146382258 0.09076859
## 21 21 0.067360025 0.06807600
## 22 22 -0.028589935 -0.07262947
## 23 23 0.054018006 0.10007483
## 24 24 0.024960463 0.00453503
## 25 25 -0.104424982 -0.10347838##maintenant on a autocorrelation en lag 2 et 7.. on ajoute lag 2 pour tenter regler le probleme
#HISTOGRAME des residus
#ggplot(dem[app,], aes(x=RES_3)) +
# geom_histogram(bins=10)
#Test de normalite des residus (Shapiro-Wilk)
#pander(shapiro.test(dem[["RES_3"]][app]))
#Multicolinearite (plus d'une variable explicative)
#pander(vif(model3)) # fonctionne pas on a petite nombre de variables
#Indices de performance
#indices_performance(dem[app,"DEMANDE"], dem[app,"PREV_3"])4EME ESSAIS ON FAIT LAG 1 +lag2 + lag 7
#Quatrieme modele de regression lineaire, on ajoute lag1+lag2 + lag 7
model4<-lm(DEMANDE~ t+lag(DEMANDE,1) +lag(DEMANDE,2)+ lag(DEMANDE,7),data=dem[app,])
#Ajouter les prévisions, les residus et les residus standardises
dem <- ajouter_previsions(dem, model4,residus=T,noms_personnalises="4")
#Correlogramme des residus
correlogramme(dem[app,"RES_4"], lag.max = 25)## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).
## Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).
## lag ac pac
## 1 1 0.1362909928 0.136290993
## 2 2 -0.0027341767 -0.021712730
## 3 3 -0.0164242366 -0.013338183
## 4 4 0.0479274367 0.052991440
## 5 5 -0.1210325759 -0.138353315
## 6 6 -0.0003449274 0.038732428
## 7 7 -0.0623497607 -0.071675769
## 8 8 -0.1529026853 -0.148211030
## 9 9 -0.0672525708 -0.009688560
## 10 10 -0.1140429092 -0.140469031
## 11 11 -0.0801535140 -0.046843760
## 12 12 0.0579455408 0.076534185
## 13 13 -0.0144013294 -0.087890397
## 14 14 0.0884792505 0.123866125
## 15 15 0.0020188373 -0.070987340
## 16 16 0.1090598830 0.079791395
## 17 17 0.0286623905 0.022115128
## 18 18 0.1026723605 0.021997849
## 19 19 -0.0194184021 -0.009932553
## 20 20 0.1375351078 0.136417661
## 21 21 0.1114124039 0.089497592
## 22 22 -0.0129245858 -0.001768205
## 23 23 0.0622269149 0.127248583
## 24 24 0.0080001866 -0.002829505
## 25 25 -0.1067215232 -0.054410110##maintenant on a autocorrelation en lag 2 et 7.. on ajoute lag 2 pour tenter regler le problemeCe nouveau mod`ele peut-il ˆetre consid´er´e comme autoregressif ? OUI, On peut le considerer autoregressif car il y a des variables explicatives qui sont des valeurs passées de la variable cible
Pr´esentez un r´esum´e de ce nouveau mod
elea l’aide de la fonction summary(). Que se passe-t-il de particulier avec la variable t ?
varaiable t:la variable t, qui était significative dans le modèle 1, ne l’est plus dans le modèle autorégressif. En effet, sa valeur p est de 0.322 (> 0.05), ce qui indique qu’elle n’apporte plus d’information significative une fois que les valeurs passées de la demande (lags) ont été intégrées dans le modèle.
R2: La demande est explique par le modele 86.12%
Cela signifie que la tendance temporelle captée par t est en réalité absorbée par les effets de mémoire de la série temporelle, notamment par lag(DEMANDE, 1), qui est très significatif (p < 0.0001).
pander(summary(model4))| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 109.4 | 50.49 | 2.166 | 0.0333 |
| t | -0.2931 | 0.1926 | -1.522 | 0.132 |
| lag(DEMANDE, 1) | 0.8192 | 0.104 | 7.879 | 1.503e-11 |
| lag(DEMANDE, 2) | -0.2735 | 0.1014 | -2.697 | 0.008558 |
| lag(DEMANDE, 7) | 0.2283 | 0.07476 | 3.054 | 0.003078 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 84 | 18.64 | 0.8612 | 0.8542 |
#Intervalles de confiance des coefficients de la regression
#pander(confint(model4, level = 0.95))Analysez les hypotheses de validité de ce modele. On peut décomposer ces hypothèses de la façon suivante :
➟ Taille de l’échantillon : 84 observations pour 4 variables explicatives, on respect le minimun (40 varaibles) ➟ LHomoscedasticite:Pas respecté.. on voit heteroscedasticité La relation à modéliser n’est pas linéaire :On peut observer les donnes sont pas autour de 0 (ily a une type de “u”)on peut visualiser un cone , qui montre des variation sur la variance ➟Valeur aberrantes:on a 1 valeur aberrante qui est sousestimé ➟ Normalite : Les résidus sont distribués normalement :L’ histograme des residues montre une distribution symetrique et centree, le test de Shapiro-Wilk confirme la normalité : avec une p-valeur non significative (p:0.03298 <0.05), on rejette l’hypothèse de normalité.
➟ LIndependance de residues (= aucun autocorrelation significatif), deja confirme avec lajout des variables lag dephasee ➟ Multicolinearite: le vif est > 5 dans toutes les varaibles exxcept lag7 (4.6)– Dans lensemble le modele foctionne mais on peut pas se fier de chaque variable separement
Cela remet en question la validité du modèle, mais c’est à l’utilisateur, selon le contexte, de juger s’il reste pertinent ou non.
#Graphe sequentiel des residus standardises
ggplot(data=dem[app,], aes(x=t, y = ZRES_4)) +
geom_point() +
geom_line() +
geom_hline(yintercept = 0)## Warning: Removed 7 rows containing missing values or values outside the scale range
## (`geom_point()`).## Warning: Removed 7 rows containing missing values or values outside the scale range
## (`geom_line()`).
#Correlogramme des residus
correlogramme(dem[app,"RES_4"], lag.max = 25)## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).
## lag ac pac
## 1 1 0.1362909928 0.136290993
## 2 2 -0.0027341767 -0.021712730
## 3 3 -0.0164242366 -0.013338183
## 4 4 0.0479274367 0.052991440
## 5 5 -0.1210325759 -0.138353315
## 6 6 -0.0003449274 0.038732428
## 7 7 -0.0623497607 -0.071675769
## 8 8 -0.1529026853 -0.148211030
## 9 9 -0.0672525708 -0.009688560
## 10 10 -0.1140429092 -0.140469031
## 11 11 -0.0801535140 -0.046843760
## 12 12 0.0579455408 0.076534185
## 13 13 -0.0144013294 -0.087890397
## 14 14 0.0884792505 0.123866125
## 15 15 0.0020188373 -0.070987340
## 16 16 0.1090598830 0.079791395
## 17 17 0.0286623905 0.022115128
## 18 18 0.1026723605 0.021997849
## 19 19 -0.0194184021 -0.009932553
## 20 20 0.1375351078 0.136417661
## 21 21 0.1114124039 0.089497592
## 22 22 -0.0129245858 -0.001768205
## 23 23 0.0622269149 0.127248583
## 24 24 0.0080001866 -0.002829505
## 25 25 -0.1067215232 -0.054410110##maintenant on a autocorrelation en lag 2 et 7.. on ajoute lag 2 pour tenter regler le probleme
#HISTOGRAME des residus
ggplot(dem[app,], aes(x=RES_4)) +
geom_histogram(bins=10) ## Warning: Removed 7 rows containing non-finite outside the scale range
## (`stat_bin()`).
#Test de normalite des residus (Shapiro-Wilk)
pander(shapiro.test(dem[["RES_4"]][app]))| Test statistic | P value |
|---|---|
| 0.9677 | 0.03298 * |
#Multicolinearite (plus d'une variable explicative)
pander(vif(model4)) # fonctionne pas on a petite nombre de variables| t | lag(DEMANDE, 1) | lag(DEMANDE, 2) | lag(DEMANDE, 7) |
|---|---|---|---|
| 5.272 | 6.139 | 5.991 | 4.346 |
#Indices de performance
indices_performance(dem[app,"DEMANDE"], dem[app,"PREV_4"])## N MAD MAPE
## Valeur de l'indice 84 13.32427 0.029881#Multicolinearite (plus d'une variable explicative)
pander(vif(model4))| t | lag(DEMANDE, 1) | lag(DEMANDE, 2) | lag(DEMANDE, 7) |
|---|---|---|---|
| 5.272 | 6.139 | 5.991 | 4.346 |
QUESTION 4 (8 points) Troisieme modele : On tentera
d’am´eliorer la qualit´e de notre
modele en ajoutant la temp´erature observ´ee parmi les variables explicatives. (a) Laquelle des variables de temp´erature est la plus significative lorsqu’on les ajoute individuellement au modele
final de la question 3 ? Parmi les 3 variables, TEMP_MAX est celle qui
ameliore le plus le modele. Elle est la plus
significative(p:7.282e-08<0.001), son effet est négatif et confirmé
par l’intervalle de confia nce, et elle permet d’obtenir le meilleur R²
ajusté (0.899).
#on doit fusioner les deux tables
temp$DATE <- as.Date(temp$DATE)
dem_temp <- dem %>%
left_join(temp, by = "DATE") %>%
mutate(t = row_number())
app<-c(1:91)
essai<-c(92:122)
# modele avec temp_moy
modele_temp_moy <- lm(DEMANDE ~ t + lag(DEMANDE, 1) + lag(DEMANDE, 2) + lag(DEMANDE, 7) + TEMP_MOY, data = dem_temp[app,])
summary(modele_temp_moy)##
## Call:
## lm(formula = DEMANDE ~ t + lag(DEMANDE, 1) + lag(DEMANDE, 2) +
## lag(DEMANDE, 7) + TEMP_MOY, data = dem_temp[app, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -36.559 -10.441 -0.501 8.624 44.011
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 109.36421 44.40543 2.463 0.0160 *
## t 0.49296 0.23296 2.116 0.0375 *
## lag(DEMANDE, 1) 0.62925 0.09928 6.338 1.39e-08 ***
## lag(DEMANDE, 2) -0.14966 0.09268 -1.615 0.1104
## lag(DEMANDE, 7) 0.30378 0.06752 4.499 2.35e-05 ***
## TEMP_MOY -3.46836 0.70587 -4.914 4.83e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.39 on 78 degrees of freedom
## (7 observations deleted due to missingness)
## Multiple R-squared: 0.894, Adjusted R-squared: 0.8872
## F-statistic: 131.6 on 5 and 78 DF, p-value: < 2.2e-16pander(confint(modele_temp_moy, level = 0.95))| 2.5 % | 97.5 % | |
|---|---|---|
| (Intercept) | 20.96 | 197.8 |
| t | 0.02917 | 0.9568 |
| lag(DEMANDE, 1) | 0.4316 | 0.8269 |
| lag(DEMANDE, 2) | -0.3342 | 0.03485 |
| lag(DEMANDE, 7) | 0.1694 | 0.4382 |
| TEMP_MOY | -4.874 | -2.063 |
MODELE AVEC TEMP MAX
# modele avec temp_max
modele_temp_max <- lm(DEMANDE ~ t + lag(DEMANDE, 1) + lag(DEMANDE, 2) + lag(DEMANDE, 7) + TEMP_MAX, data = dem_temp[app,])
pander(summary(modele_temp_max))| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 109.3 | 42.15 | 2.593 | 0.01134 |
| t | 0.4365 | 0.2022 | 2.159 | 0.03395 |
| lag(DEMANDE, 1) | 0.6003 | 0.09428 | 6.367 | 1.229e-08 |
| lag(DEMANDE, 2) | -0.1256 | 0.08823 | -1.424 | 0.1585 |
| lag(DEMANDE, 7) | 0.3254 | 0.06451 | 5.045 | 2.89e-06 |
| TEMP_MAX | -2.783 | 0.4679 | -5.947 | 7.282e-08 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 84 | 15.56 | 0.9045 | 0.8984 |
pander(confint(modele_temp_max, level = 0.95))| 2.5 % | 97.5 % | |
|---|---|---|
| (Intercept) | 25.4 | 193.2 |
| t | 0.03393 | 0.8391 |
| lag(DEMANDE, 1) | 0.4126 | 0.788 |
| lag(DEMANDE, 2) | -0.3013 | 0.05004 |
| lag(DEMANDE, 7) | 0.197 | 0.4539 |
| TEMP_MAX | -3.714 | -1.851 |
MODELE AVEC TEMP_MIN
# modele avec temp_min
modele_temp_min <- lm(DEMANDE ~ t + lag(DEMANDE, 1) + lag(DEMANDE, 2) + lag(DEMANDE, 7) + TEMP_MIN, data = dem_temp[app,])
summary(modele_temp_min)##
## Call:
## lm(formula = DEMANDE ~ t + lag(DEMANDE, 1) + lag(DEMANDE, 2) +
## lag(DEMANDE, 7) + TEMP_MIN, data = dem_temp[app, ])
##
## Residuals:
## Min 1Q Median 3Q Max
## -43.923 -8.472 -1.405 10.991 62.614
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 109.025622 49.885478 2.186 0.03185 *
## t 0.002985 0.256995 0.012 0.99076
## lag(DEMANDE, 1) 0.771483 0.106443 7.248 2.65e-10 ***
## lag(DEMANDE, 2) -0.244720 0.101591 -2.409 0.01836 *
## lag(DEMANDE, 7) 0.241717 0.074275 3.254 0.00168 **
## TEMP_MIN -1.542494 0.900211 -1.713 0.09060 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 18.41 on 78 degrees of freedom
## (7 observations deleted due to missingness)
## Multiple R-squared: 0.8662, Adjusted R-squared: 0.8576
## F-statistic: 101 on 5 and 78 DF, p-value: < 2.2e-16pander(confint(modele_temp_min, level = 0.95))| 2.5 % | 97.5 % | |
|---|---|---|
| (Intercept) | 9.711 | 208.3 |
| t | -0.5087 | 0.5146 |
| lag(DEMANDE, 1) | 0.5596 | 0.9834 |
| lag(DEMANDE, 2) | -0.447 | -0.04247 |
| lag(DEMANDE, 7) | 0.09385 | 0.3896 |
| TEMP_MIN | -3.335 | 0.2497 |
- S’il est impossible d’utiliser la temp´erature observ´ee la journ´ee
mˆeme
a laquelle on souhaite ´etablir une pr´evision, pourrait-on se r´ef´erera la temp´erature observ´ee pr´ec´edemment ? oui, *(utilise temp lag )A l’aide du corr´elogramme crois´e entre la variable rep´er´ee en (a) et les r´esidus du modele final de la question 3, seriez-vous tent´e d’ajouter une variable d´ephas´ee de la temp´erature observ´ee dans votre mod`ele de r´egression lin´eaire ? Si oui, quel ordre de d´ephasage serait sugg´er´e par le corr´elogramme ?
Le corrélogramme croisé entre TEMP_MAX et RES_4 présente une corrélation significative à un décalage de 8 jours. Cela suggère que la température maximale observée 8 jours avant pourrait aider à mieux expliquer la demande d’électricité actuelle (mes residues). On pourrait donc tester l’ajout de la variable lag(TEMP_MAX, 8) dans un futur modèle.
#prev4 <- prev4 %>%
# left_join(temp, by = "DATE")
#app <- prev4$DATE < "2019-07-01"
#Correlations croisees
correlogramme("TEMP_MAX","RES_4", data = dem_temp[app, ])## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_hline()`).
## lag cc1 cc2
## 0 0 -0.2941123187 -0.294112319
## 1 1 0.0392900819 -0.255885237
## 2 2 0.1496562379 -0.089587805
## 3 3 0.1046773608 0.015519939
## 4 4 0.0741533817 0.018413278
## 5 5 0.0707307339 0.009947674
## 6 6 0.0568055149 -0.003769565
## 7 7 0.1742491630 0.063238080
## 8 8 0.2601820839 0.074923915
## 9 9 0.1305451227 0.059353662
## 10 10 0.1255983289 0.088290894
## 11 11 0.1129265679 0.042265015
## 12 12 0.0428530907 -0.016685180
## 13 13 0.0513239417 -0.052371079
## 14 14 0.0227418137 -0.067467500
## 15 15 0.0392305232 -0.135553685
## 16 16 -0.0553076461 -0.113531292
## 17 17 -0.0084143259 -0.151998692
## 18 18 0.0205771611 -0.149042804
## 19 19 -0.0002972628 -0.129850526
## 20 20 -0.0230928000 -0.099048385
## 21 21 0.0788206063 -0.089801527
## 22 22 0.0830445339 -0.043121847
## 23 23 0.0718742861 -0.117789980
## 24 24 0.1433515545 -0.061394766
## 25 25 0.2223393916 -0.001779372# on veut expliquer les residues aver le temp max dephasse - On d´ecide de conserver le mod
ele repéréa la sous-question (a). (DONC AVEC TEMP_MAX)
Pr´esentez le graphe s´equentiel des pr´evisions statiques sur la plage d’essai(test). Faites de mˆeme avec les pr´evisions dynamiques. Ces deux graphes devront aussi contenir les valeurs r´eelles et les intervalles de pr´evision.
Graphe sequentiel avec previsions et intervalles de prevision statiques
#Ajout des previsisons
prev_essai <- ajouter_previsions(dem_temp, modele_temp_max,
intervalles=0.95, residus=T,noms_personnalises="4")
#prevision de toutes les donnes
#Graphe sequentiel avec previsions et intervalles de prevision statiques
#sur la ESSAI
ggplot(data = prev_essai[essai,], aes(x = DATE, y = DEMANDE)) +
geom_point() +
geom_line() +
geom_line(aes(y = PREV_4), color = "red") +
geom_line(aes(y = IP_95_INF_4), color = "blue", linetype = "dashed") +
geom_line(aes(y = IP_95_SUP_4), color = "blue", linetype = "dashed")
#rouge previsions et poiunt noir sont les occupancy reels... comment le modele se comporte avec des donnes inconnuGraphe dynamique avec previsions et intervalles de prevision dynamique
#Ajout des previsisons
prev_dynamique <- ajouter_previsions(dem_temp, modele_temp_max, plage_connue=app, #
intervalles=0.95, residus=T,noms_personnalises="5")
#Graphe sequentiel avec previsions et intervalles de prevision dynamique
#sur la ESSAI
ggplot(data = prev_dynamique[essai,], aes(x = DATE, y = DEMANDE)) +
geom_point() +
geom_line() +
geom_line(aes(y = PREV_5), color = "red") +
geom_line(aes(y = IP_95_INF_5), color = "blue", linetype = "dashed") +
geom_line(aes(y = IP_95_SUP_5), color = "blue", linetype = "dashed")
(d) Selon vous, pourquoi les pr´evisions dynamiques semblent-elles
donner de moins bons r´esultats que les pr´evisions statiques ?
Les prévisions statiques s’appuient sur les valeurs réelles observées pour générer la prédiction suivante. Cela limite la propagation des erreurs.
En revanche, les prévisions dynamiques utilisent les valeurs prédites des périodes précédentes comme entrées pour les périodes suivantes. Cela fait effet cumulatif des erreurs, ce qui rend les prévisions dynamiques souvent moins fiables sur la durée
QUESTION 5 (6 points) Quatrieme modele : Un
collegue affirme qu’il a trouv´e un meilleur modele, avec
un r2aj plus ´elev´e que le vˆotre. Dans son
modele, il propose plutˆot d’utiliser les variables explicatives suivantes : le jour de la semaine, la temp´erature observ´ee rep´er´eea
la question 4 (a) et la demande d’´electricit´e de la journ´ee
pr´ec´edente. (a) Pr´esentez un r´esum´e du
modele propos´e par votre collegue `a l’aide de la fonction
summary(). Vous devrez s´electionner un jour de r´ef´erence appropri´e.
Jai pris le samedi comme jour de reference
#jour_semaine as factor
dem_temp$jour_semaine <- factor(dem_temp$jour_semaine,ordered=F)
#chisisr 1 jour de la semaine comme reference
dem_temp$jour_semaine <- relevel(dem_temp$jour_semaine, ref = "Saturday") # avec monday toute -
#Premier modele de regression lineaire
modele_college <- lm(DEMANDE ~ jour_semaine + TEMP_MAX+ lag(DEMANDE, 1),
data = dem_temp[app,])
#Resume des statistiques du modele
pander(summary(modele_college))| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 161.2 | 26.8 | 6.017 | 4.908e-08 |
| jour_semaineMonday | 35.07 | 6.451 | 5.436 | 5.609e-07 |
| jour_semaineTuesday | 23.28 | 6.191 | 3.76 | 0.0003199 |
| jour_semaineWednesday | 23.55 | 6.15 | 3.829 | 0.0002525 |
| jour_semaineThursday | 23.18 | 6.146 | 3.772 | 0.0003072 |
| jour_semaineFriday | 19.84 | 6.143 | 3.23 | 0.001792 |
| jour_semaineSunday | 9.816 | 6.189 | 1.586 | 0.1166 |
| TEMP_MAX | -2.314 | 0.3651 | -6.338 | 1.236e-08 |
| lag(DEMANDE, 1) | 0.676 | 0.0474 | 14.26 | 8.372e-24 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 90 | 15.66 | 0.926 | 0.9187 |
#Ajout des previsisons
prev_collegue <- ajouter_previsions(dem_temp, modele_college, #
intervalles=0.95, residus=T,noms_personnalises="7")
#Indices de performance
#MAPE(prev_collegue$PREV_7, prev_collegue$DEMANDE,app=app )
indices_performance(prev_collegue$DEMANDE, prev_collegue$PREV_7,app=app)## N MAD MAPE
## Apprentissage 90 12.08765 0.027074
## Essai 31 22.43960 0.053458# mon code
# modele avec temp_max
modele_temp_max <- lm(DEMANDE ~ t + lag(DEMANDE, 1) + lag(DEMANDE, 2) + lag(DEMANDE, 7) + TEMP_MAX, data = dem_temp[app,])
pander(summary(modele_temp_max))| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 109.3 | 42.15 | 2.593 | 0.01134 |
| t | 0.4365 | 0.2022 | 2.159 | 0.03395 |
| lag(DEMANDE, 1) | 0.6003 | 0.09428 | 6.367 | 1.229e-08 |
| lag(DEMANDE, 2) | -0.1256 | 0.08823 | -1.424 | 0.1585 |
| lag(DEMANDE, 7) | 0.3254 | 0.06451 | 5.045 | 2.89e-06 |
| TEMP_MAX | -2.783 | 0.4679 | -5.947 | 7.282e-08 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 84 | 15.56 | 0.9045 | 0.8984 |
pander(confint(modele_temp_max, level = 0.95))| 2.5 % | 97.5 % | |
|---|---|---|
| (Intercept) | 25.4 | 193.2 |
| t | 0.03393 | 0.8391 |
| lag(DEMANDE, 1) | 0.4126 | 0.788 |
| lag(DEMANDE, 2) | -0.3013 | 0.05004 |
| lag(DEMANDE, 7) | 0.197 | 0.4539 |
| TEMP_MAX | -3.714 | -1.851 |
#Ajouter les prévisions, les residus et les residus standardises
prev_temp <- ajouter_previsions(dem_temp, modele_temp_max,intervalles=95,residus=T,noms_personnalises="2")
#indices
indices_performance(prev_temp$DEMANDE, prev_temp$PREV_2,app=app)## N MAD MAPE
## Apprentissage 84 11.97307 0.027461
## Essai 31 15.32504 0.036606- Pourquoi serait-il pr´ef´erable de baser le choix du mod`ele sur les r´esultats de la plage d’essai plutˆot que sur les r´esultats de la plage d’apprentissage ? Il vaut mieux choisir le modèle en regardant les résultats sur la plage d’essai, car cela montre s’il est capable de bien prévoir de nouvelles données,sur les donnes inconnu. Cela aide aussi à éviter le surapprentissage.