Le(s) papier(s) d’origine

Abordons une autre anomalie exotique documentée comme affectant les rendements boursiers. Celle documentée par Kamstra et al. (2000) que l’on nomme la Daylight Saving Anomaly. Il s’agit de la tendance qu’aurait les marchés à baisser sur la première séance de cotation suivant les changements d’heures. Ces changements d’heures sont des mesures réglementaires imposant sur un pays l’avance ou le recule des horloges d’une heure afin de permettre des économies d’énergie. L’idée est de limiter les coûts d’éclairage en permettant une meilleure adéquation entre l’activité humaine (la journée de travail ordinaire) et la présence de la lumière naturelle (le jour). Le mécanisme sous-jacent de l’anomalie serait que suite à ces changements d’heures les agents intervenant sur le marché auraient eu plus fréquemment un sommeil perturbé du fait du changement de rythme imposé. Cela les amènerait à être de mauvaise humeur, augmenterait leur anxiété, ce qui se concrétiserait dans une plus grande aversion au risque (temporaire) le temps de s’habituer au décalage (une journée). Ainsi, les cours de la première journée de cotation suite aux changements d’heures entraînerait une pression vendeuse conduisant à une tendance à la baisse des prix.

Encore une fois, comme la précédente (les phases de la lune), cette anomalie exotique fait l’objet de contestations. Je vous renvoie comme précédemment à Kim and Shamsuddin (2023) pour en avoir un aperçu. J’ajouterais ici Pinegar (2002) et Müller et al. (2009) qui traitent/critiquent plus spécifiquement de cette seule anomalie (ce que nous ferons également).

Chargeons les packages et les données

De notre coté, nous reprendrons le contexte de l’article d’origine en travaillant sur l’indice Standard & Poor’s 500 sur une période assez longue : entre 1967 et fin 2022. Sur cette période, les règles de changements d’heure (heure d’été/heure d’hivers) ont évolué. La première réglementation allant dans ce sens pour les Etats Unis, le Daylight Saving Act, entre en vigueur début 1967 (notre point de départ).

Mais avant d’aller plus loin, chargeons le tidyverse. Comme d’habitude, nous utiliserons la syntaxe qui en est issue pour la suite. Chargeons également gt qui nous servira à la mise en forme de certaines tables

library(tidyverse)
library(gt)

Ceux-ci étant fait. Voyons pour les données. Afin d’obtenir les cours de l’indice directement sous la forme d’une data frame, utilisons la fonction tq_get() du package tidyquant. Précisons le ticker de l’indice, ^GSPC, ainsi que les dates de début et de fin de notre série temporelle. Limitons les informations à la date et cours ajustés (adjusted) et renommons ce dernier SP500 pour plus de clarté.

dat<-tidyquant::tq_get("^GSPC",from="1966-12-30",to="2022-12-30") |> 
  select(date,adjusted) |> 
  rename(SP500=adjusted)

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

head(dat)

## # A tibble: 6 × 2
##   date       SP500
##   <date>     <dbl>
## 1 1966-12-30  80.3
## 2 1967-01-03  80.4
## 3 1967-01-04  80.6
## 4 1967-01-05  81.6
## 5 1967-01-06  82.2
## 6 1967-01-09  82.8

Traitement des données

Maintenant que nous avons les données, mettons les en forme afin de pouvoir effectuer nos tests.

Commençons par calculer les rendements journaliers continus de l’indice. Après tout, notre hypothèse de recherche porte sur eux. Ces rendements seraient en moyenne plus petits le premier jour de cotation après un weekend de changement d’heure.

dat<-dat |> 
  mutate(ret=log(SP500/lag(SP500))) 
head(dat)

## # A tibble: 6 × 3
##   date       SP500       ret
##   <date>     <dbl>     <dbl>
## 1 1966-12-30  80.3 NA       
## 2 1967-01-03  80.4  0.000622
## 3 1967-01-04  80.6  0.00211 
## 4 1967-01-05  81.6  0.0130  
## 5 1967-01-06  82.2  0.00708 
## 6 1967-01-09  82.8  0.00764

Passons maintenant à la construction de notre variable teste. Il s’agit d’une binaire (dummy) marquant le fait qu’un jour de cotation donné est (1) ou non (0) le premier jour de cotation qui suit un changement d’heure (qui intervient le weekend précédent).

Pour ce faire, commençons par mettre à plat l’évolution sur la période d’étude du Daylight Saving Act qui entre en vigueur en 1967. Celui-ci s’applique encore de nos jours. En fait, il n’a été suspendu qu’une fois pour une année en 1974. Il consiste en la programmation d’un changement d’heure (avance d’une heure, recule d’une heure) à deux moments de l’année à l’automne et au printemps pour limiter les coûts d’éclairage. Par contre, les dates des changements d’heure ont au fils du temps changé. Ainsi, jusqu’en 2006, le changement d’heure d’automne était réalisé le dernier dimanche d’octobre. Depuis, 2007 ce changement s’opère le premier dimanche de novembre. Concernant le changement d’heure du printemps, on a eu deux modifications. Jusqu’en 1986, le changement s’opèrait le dernier dimanche d’avril. Entre 1987 et 2006, il s’est opéré le premier dimanche de ce même mois. Depuis 2007, il s’opère le second dimanche de mars.

Vous voyez les choses ne sont pas si simples et la variable à coder va nécessiter un peu de travail…

Commençons par transformer notre variable de date de cotation en un série d’autres variables que nous pourrons, par la suite, mobiliser plus facilement pour l’établir. Sortons le jours de la semaine (lundi, mardi, mercredi…), le jours du mois (1,2,3,…), le nombre de jours dans le mois considéré (28, 29, 30, 31), le mois et l’année. Créons une variables pour chacune de ces informations. Elles servirons par la suite.

A cette fin, nous utiliserons une série de fonctions issues du même package lubridate qui est inclue dans le tidyverse. Pour le jours de la semaine, utilisons la fonction wday(). Demandons lui de commencer la semaine le lundi et de ne pas afficher le texte du jours mais son numéro (1 pour lundi; 2 pour mardi …). Cela est plus pratique pour la suite. Pour le jour du mois, utilisons la fonction day(). Pour le mois, la fonction month(). Pour le nombre de jours dans le mois la fonction days_in_month(). Enfin, pour l’année, la fonction year().

dat<-dat |> 
  mutate(
    jours_semaine=wday(date,week_start = 1,label = FALSE, abbr = FALSE),
    jours_du_mois=day(date),
    mois=month(date),
    nb_j_mois=days_in_month(date),
    annee=year(date))
head(dat)

## # A tibble: 6 × 8
##   date       SP500       ret jours_semaine jours_du_mois  mois nb_j_mois annee
##   <date>     <dbl>     <dbl>         <dbl>         <int> <dbl>     <int> <dbl>
## 1 1966-12-30  80.3 NA                    5            30    12        31  1966
## 2 1967-01-03  80.4  0.000622             2             3     1        31  1967
## 3 1967-01-04  80.6  0.00211              3             4     1        31  1967
## 4 1967-01-05  81.6  0.0130               4             5     1        31  1967
## 5 1967-01-06  82.2  0.00708              5             6     1        31  1967
## 6 1967-01-09  82.8  0.00764              1             9     1        31  1967

Travaillons maintenant à marquer les jours de cotation venant potentiellement juste après un weekend de changement d’heure. Autrement-dit, ceux arrivant immédiatement après le premier dimanche du mois, ceux arrivant immédiatement après le second dimanche du mois et ceux arrivant immédiatement après le dernier dimanche du mois. Créons une variable binaire pour chacun de ces cas (c’est ici, à mon avis, la tâche le plus compliquée dans notre réplication).

Pour ce faire, établissons une série de conditions cumulatives. Utilisons ainsi un empilement de fonction ifelse(). La première condition sera la même pour les trois variables. Il s’agit de marquer d’un 0 les jours de cotation qui ne sont pas les premiers d’une semaine. Ce faisant, nous utilisons la relation suivante: si le jours de la semaine (marqué par un chiffre 1 pour le lundi, 2 pour le mardi, …) et plus grand que celui qui le précède c’est que l’on n’est pas sur le premier jour d’une semaine de cotation. Il y a ici 2922 observations dans ce cas. Pour le reste, les choses divergent légèrement d’une variable à l’autre.

Afin d’identifier le jours suivant le premier dimanche, on commence par marquer d’un 0 les jours dont la différence entre le numéro (jours du mois) et le numéro du jour de la semaine est plus grand ou égale à 8. Autrement-dit, on exclut les jours qui ne seraient pas dans la premier semaine de cotation du mois plus un jour (le dimanche qui précède). Les jours non exclus, ceux non marqués jusque là par un 0, sont alors 777. Sur ces observations, on exclut les cas où il n’y a pas eu de dimanches avant ces jours dans le mois. C’est le cas quand le numéro du jours du mois et plus grand que numéro du jours de la semaine (si ce n’était pas le cas cela marquerait le fait que le mois à commencer un lundi donc que le premier jour de cotation après le premier dimanche du mois arrivera la semaine qui suit).

L’identification du jours suivant le second dimanche du mois suit schématiquement la même procédure. Deux ajustements sont néanmoins à réaliser. Pour le second critère d’exclusion, on marque d’un 0 les dates qui ne seraient pas dans les deux premières semaines du mois. Celles pour qui la différence entre le jours du mois et le jours de la semaine est supérieure ou égale à 15. Pour le dernier critère, il s’agit d’exclure les dates correspondant à la première semaine, celles intervenant après le premier dimanche mais avant le second. Pour cela, nous nous assurons que la différence entre le numéro du jours du mois et le numéro du jours de la semaine soit supérieur à 7.

L’identification du premier jours de cotation précédent le dernier dimanche du mois est un peu plus compliquée. Les deux étapes suivant l’identification du premier jour de cotation de la semaine sont différentes de celles mobilisées précédemment. La première consiste à marquer les jours qui sont dans la dernière semaine du mois. Cela est réalisée à partir de la condition suivante: si la différence entre le numéro du jour du mois et du numéro du jours de la semaine plus 7 (une semaine qui fini un dimanche jour 7) est supérieure au nombre de jours total dans le mois, on est bien sur des jours de la dernière semaine. On marque ainsi les cas où le premier jour ouvré arrivant après le dernier dimanche du mois et dans le même mois que le dimanche. Ils sont alors 567. Reste donc à pointer les cas où le premier jours de cotation après le dernier dimanche du mois se retrouve le mois suivant. C’est ce que fait la dernière étape. Elle pointe les jours de la première semaine du mois, ceux pour qui le jour du mois est plus petit ou égale au jour de la semaine.

dat_<-dat |> 
  mutate(
    j_ap_prem_dim=ifelse(jours_semaine>lag(jours_semaine),0,
                         ifelse(jours_du_mois-jours_semaine>=8,0,
                                ifelse(jours_du_mois-jours_semaine>0,1,0))),
    j_ap_sec_dim=ifelse(jours_semaine>lag(jours_semaine),0,
                         ifelse(jours_du_mois-jours_semaine>=15,0,
                                ifelse(jours_du_mois-jours_semaine>7,1,0))),
    j_ap_der_dim=ifelse(jours_semaine>lag(jours_semaine),0,
                        ifelse(jours_du_mois-jours_semaine+7>nb_j_mois,1,
                              ifelse(jours_du_mois-jours_semaine<=0,1,0))),
    lundi=jours_semaine==1
  )
head(dat_)

## # A tibble: 6 × 12
##   date       SP500       ret jours_semaine jours_du_mois  mois nb_j_mois annee
##   <date>     <dbl>     <dbl>         <dbl>         <int> <dbl>     <int> <dbl>
## 1 1966-12-30  80.3 NA                    5            30    12        31  1966
## 2 1967-01-03  80.4  0.000622             2             3     1        31  1967
## 3 1967-01-04  80.6  0.00211              3             4     1        31  1967
## 4 1967-01-05  81.6  0.0130               4             5     1        31  1967
## 5 1967-01-06  82.2  0.00708              5             6     1        31  1967
## 6 1967-01-09  82.8  0.00764              1             9     1        31  1967
## # ℹ 4 more variables: j_ap_prem_dim <dbl>, j_ap_sec_dim <dbl>,
## #   j_ap_der_dim <dbl>, lundi <lgl>

Avant d’aller plus loin vérifions que l’on a bien pour chaque variable un jour marqué par mois. Pour ce faire, commençons par compter le nombre de mois dans la base. Excluons du décompte la première observation (le 30 décembre 1966) qui n’a été mobilisée que pour le calcul du premier rendement.

dat$date[-1] |> format("%m-%Y") |> unique() |> length()

## [1] 672

La base contient 672 mois et donc chacune de nos variables “jours après” devrait compter ce même nombre de marquages (de 1). Voyons ce qu’il en est en sommant une à une les valeurs de ces variables.

dat_ |> summarise(sum_ap_prem_dim=sum(j_ap_prem_dim,na.rm=TRUE),
                  sum_ap_sec_dim=sum(j_ap_sec_dim,na.rm=TRUE),
                  sum_ap_der_dim=sum(j_ap_der_dim,na.rm=TRUE))

## # A tibble: 1 × 3
##   sum_ap_prem_dim sum_ap_sec_dim sum_ap_der_dim
##             <dbl>          <dbl>          <dbl>
## 1             672            672            672

C’est bien le cas. Il semble qu’elles soient correctement codées. Elles vont servir à établir notre variable teste.

Il s’agit alors de coder notre variable teste qui marque le premier jour de cotation suite à un weekend de changement d’heure. Nommons la day_light_sav. Il s’agit d’une variable binaire. La difficulté est ici associée à l’évolution de la réglementation dans le temps. Afin de faciliter le codage, utilisons la fonction case_when() du package dplyr qui est inclus dans le tidyverse. Elle permet d’articuler très facilement de multiples conditions via une syntaxe simplifiée (notamment par rapport à ifelse()). Celle-ci consiste à poser en série les conditions de codage via des “phrases” posant dans un premier temps la condition logique à tester et d’indiquer suite au symbole ~ la valeur qui sera inscrite si cette condition est vraie. Ainsi, par exemple, dans notre cas, la première condition posée est le fait que l’année soit 1974, ce qui est indiqué si la variable annee est égale à 1974, dans ce cas la valeur sera 0. 1974 est une année où il n’y a pas de changement d’heure. Pour les autres couples conditions/valeur, nous les organisons en fonction du fait qu’elles concernent l’heure d’hivers ou l’heure d’été. Suite à l’énoncé de l’ensemble de ces phrases de codage, il demeure un ensemble de cas ne répondant à aucune des conditions énoncées. Nous utilisons l’option .default pour leur attribuer la valeur 0.

dat_<-dat_ |> 
  mutate(
    day_light_sav=case_when(
      # pas de changement d'heure
      annee==1974~0,
      # passages à l'heure "d'hivers"
      annee>=2007&mois==11&j_ap_prem_dim==1~1,
      annee<=2006&mois==10&j_ap_der_dim==1~1,
      # passages à l'heure "d'été"
      annee>=2007&mois==3&j_ap_sec_dim==1~1,
      annee>=1987&annee<=2006&mois==4&j_ap_prem_dim==1~1,
      annee<=1986&mois==4&j_ap_der_dim==1~1,
     .default = 0))

Voyons combien de jours sont ici marqués comme des premiers jours de cotation suite à un changement d’heure.

dat_$day_light_sav |> sum()

## [1] 109

Sur nos 14 096 observations (jours de cotation), on en a que 109 qui sont issues des premiers jours de cotation suite à un changement d’heure.

Analyses et tests statistiques

Maintenant que nos données sont prêtes, procédons aux tests de notre hypothèse (des rendements sont en moyenne plus petits le premier jour de cotation précédant un changement d’heure). Pour ce faire, procédons directement à une série de régressions. Commençons simplement par le modèle associant uniquement notre variable day_light_sav et le modèle incluant en plus la variable de contrôle lundi. Présentons les dans la même table.

reg1<-lm(ret~day_light_sav ,data=dat_)
reg2<-lm(ret~day_light_sav+lundi ,data=dat_)

modelsummary::modelsummary(list(reg1,reg2),
                           stars=c('*' = 0.1, '**' = 0.05,'***'=0.01),
                           gof_map = c("nobs","r.squared","F"))

Notons que lorsque la variable teste est seule dans la régression, le coefficient associé apparaît bien négatif et statistiquement significatif mais au seuil de 10% (ce qui est limité). Par ailleurs, le R2 de la régression est très petit ce qui signifie que notre variable n’explique qu’une part négligeable de la variance des rendements quotidiens. L’introduction de la variable de contrôle effet lundi détruit cette faible significativité. Elle est elle-même négative et statistiquement significative au seuil de 1%.

Notez que le test proposé par la régression 1 (la régression simple) et simplement l’équivalent d’un test de Student simple faisant l’hypothèse d’égalité des variances sur les deux sous-échantillons comparés.

t.test(ret~day_light_sav,data=dat_,var.equal=TRUE)$p.value

## [1] 0.07842915

C’est ici une hypothèse forte. Les rendements boursiers quotidiens sont bien souvent caractérisés par un certain degré d’autocorrélation et d’hétéoscédasticité. Autrement-dit, les rendements d’un jours ne sont généralement pas totalement indépendant de ceux de la veille et la variance de ces rendements peut être différente d’une période à l’autre. En période de crise, elle est plus importante qu’en période normale. Sur notre période d’étude, les marchés ont connu plusieurs crises. Les tests de Student pratiqués dans les régressions précédentes ne sont donc pas fiables. Il faut ajuster la matrice de corrélations des erreurs pour corriger ces défauts. Pour ce faire, mobilisons la correction de Newey-West.

modelsummary::modelsummary(list(reg1,reg2),
                           vcov="NeweyWest",
                           stars=c('*' = 0.1, '**' = 0.05,'***'=0.01),
                           gof_map = c("nobs","r.squared","F"))

Le correctif détruit la significativité (qui était déjà faible) de la variable teste dans le premier modèle et n’affecte pas le diagnostique vis-à-vis de la variable de contrôle lundi introduite dans le modèle 2. Tout ceci ne plaide pas en faveur de l’existence d’un effet changement d’heure.

Néanmoins, on pourrait penser que la diffusion de l’étude de Kamstra, Kramer et Levi (avant sa publication) à partir de 1998 ait conduit à la réalisation d’arbitrages qui ont détruit l’anomalie (plombant de fait la significativité sur l’ensemble de la période). Voyons cela en découpant l’échantillon.

Commençons par présenter les résultats brutes sans ajustement des tests. Intégrons en surplus une spécification présentant une variable d’interaction entre notre variable teste day_light_sav et une variable indicatrice marquant le fait que les observations soit ou non celles d’avant 1998 (1) ou celles d’après (0). Si elle était statistiquement signifcative, cela montrerait une différence d’effet du changement d’heure avant et après 1998.

dat_<-dat_ |> 
  mutate(avant_1998=annee<=1997)

reg3<-lm(ret~day_light_sav ,data=dat_|> filter(annee<=1997))
reg4<-lm(ret~day_light_sav+lundi ,data=dat_|> filter(annee<=1997))
reg5<-lm(ret~day_light_sav ,data=dat_|> filter(annee>=1998))
reg6<-lm(ret~day_light_sav+lundi ,data=dat_|> filter(annee>=1998))
reg7<-lm(ret~day_light_sav*avant_1998+lundi ,data=dat_)

modelsummary::modelsummary(list(reg3,reg4,reg5,reg6,reg7),
                           fmt=5,
                           stars=c('*' = 0.1, '**' = 0.05,'***'=0.01),
                           gof_map = c("nobs","r.squared","F"),
                           output = "gt") |> 
              tab_spanner(label = "Avant 1998", columns = c("(1)", "(2)")) |> 
              tab_spanner(label = "Après 1998", columns = c("(3)", "(4)")) |> 
              tab_spanner(label = "Échantillon complet", columns = "(5)") 

Sur la période avant 1998, dans les deux régressions (avec et sans la variable de contrôle), la variable teste est statistiquement significative, au seuil de 1% quand elle est seule et au seuil de 5% avec la variable de contrôle (lundi). Cela semble plaider en la faveur de la présence de l’anomalie avant la publication de l’article. Sur la période qui suit (1998 et après), ces coefficients ne sont pas statistiquement significatifs (ici aucun coefficient ne l’est). Plus surprenant, c’est également le cas de la variable d’interaction dans la dernière régression. Elle est bien négative comme le laissait penser la régression 2 mais non significativement différente de 0 à un niveau acceptable. Cela pose question. Si l’effet de notre variable teste est significatif sur la première période et pas sur la seconde, la variable d’interaction devrait être significative.

Procédons aux mêmes tests mais cette fois en prenant en compte les éventuels problèmes associés à l’autocorrélation et à l’hétéroscédasticité via le correctif de Newey-West.

modelsummary::modelsummary(list(reg3,reg4,reg5,reg6,reg7),
                           fmt=5,
                           stars=c('*' = 0.1, '**' = 0.05,'***'=0.01),
                           vcov= list(
                             sandwich::NeweyWest(reg3, prewhite = FALSE, lag = 5),
                             sandwich::NeweyWest(reg4, prewhite = FALSE, lag = 5),
                             sandwich::NeweyWest(reg5, prewhite = FALSE, lag = 5),
                             sandwich::NeweyWest(reg6, prewhite = FALSE, lag = 5),
                             sandwich::NeweyWest(reg7, prewhite = FALSE, lag = 5)),
                           gof_map = c("nobs","r.squared","F"),
                           output = "gt") |> 
  tab_spanner(label = "Avant 1998", columns = c("(1)", "(2)")) |> 
  tab_spanner(label = "Après 1998", columns = c("(3)", "(4)")) |> 
  tab_spanner(label = "Échantillon complet", columns = "(5)") 

Dans cette configuration, seule la variable de contrôle (lundi) apparaît statistiquement significative mais uniquement sur la période pré-1998 et sur l’ensemble de la période d’étude. Elle ne l’est pas sur la période post-1998. Quoi qu’il en soit notre variable teste ne l’ai plus du tout dans aucun des modèles. L’hypothèse de recherche, l’existence d’un effet changement d’heure, ne peut donc ici pas être validée.

Attention, cela ne veut pas forcément dire qu’elle doit être rejetée autrement-dit que nous avons une preuve qu’elles soient fausses. Il se peut simplement que notre configuration de test manque de puissance et donc in fine ne permette pas de trancher.

Examinons cela. Commençons par une configuration simple qui correspondrait à un test de Student pratiqué sur les rendements des premiers jours de cotation post changements d’heure et ceux des autres jours. Si on considère l’ensemble des observations, il s’agit de calculer la puissance de ce test pratiqué sur un groupe de 109 observations (les post changements) contre 13 987 autres observations/jours de cotation (14096-109). Concernant la taille de l’effet reprenons celles mis en évidence dans les régressions simple soit -0.002. Pour l’écart type des observations, calculons la simplement. Pour faciliter la procédure, mobilisons la fonction pwr.t2n.test() du package pwr qui permet le calcul de la puissance de nombreux types de test (la fonction power.t.test() n’est pas ici adaptée compte tenu de la différence de taille des groupes comparés). Retenons un seuil de confiance de 5%.

pwr::pwr.t2n.test(
  n1 = 109 , 
  n2 = 13987 , 
  d = -0.002 / sd(dat_$ret, na.rm = TRUE), # d de Cohen = delta / sd
  sig.level = 0.05,
  alternative = "two.sided")

## 
##      t test power calculation 
## 
##              n1 = 109
##              n2 = 13987
##               d = 0.1863643
##       sig.level = 0.05
##           power = 0.4912992
##     alternative = two.sided

La puissance du test est très basse : 49,12 %. On est très loin du seuil minimal de 80 % généralement admis en économétrie et en statistique. Mais qu’est-ce que cela signifie concrètement ? La puissance d’un test est égale à \(1 - \beta\), où \(\beta\) représente la probabilité d’une erreur de type II (ou risque de deuxième espèce). Ici, nous n’avons que 49,12 % de chances de ne pas commettre une erreur de type II, c’est-à-dire de ne pas faire de faux négatif (manquer la mise en évidence du phénomène alors qu’il existe réellement). Autrement dit, si l’anomalie existe bel et bien sur les marchés financiers, notre test a plus d’une chance sur deux (environ 51 %) de passer à côté sans la détecter.

Voyons ce qu’il en est sur la période avant la diffusion de la recherche (avant 1998). Comptons le nombre de changement d’heure et le nombre de jours de cotation hors premier jours post changement d’heure. On a 60 évènements et 7806 jours de cotation en dehors de ceux-ci. Encore une fois reprenons l’écart de rendements journaliers moyens à partir de la régression dédiées (-0,00302). Et calculons l’écart type des rendements sur la sous période.

pwr::pwr.t2n.test(
  n1 = 60 , 
  n2 = 7806 , 
  d = -0.00320 / sd(dat_$ret[which(dat_$annee<=1997)], na.rm = TRUE), 
  sig.level = 0.05,
  alternative = "two.sided")

## 
##      t test power calculation 
## 
##              n1 = 60
##              n2 = 7806
##               d = 0.3529968
##       sig.level = 0.05
##           power = 0.7774356
##     alternative = two.sided

Même si la situation est moins dramatique que précédemment, la puissance du test, ici 77,74%, demeure en dessus de 80%. Et donc, même sur cette période, notre test ne permet pas de trancher.

Ces diagnostiques doivent nous alarmer sur la capacité de ce type d’études à répondre à la question de recherche. Néanmoins, ici ils ne concernent que de les tests de la différence de moyennes. Voyons ce qu’il en est pour les régressions. Pour cela, nous allons mobiliser la fonction pwr.f2.test() du même package (pwr). Traitons de la période pré 1998.

r2_full <-summary(reg4)$r.squared 
r2_red  <-summary(lm(ret~lundi ,data=dat_|> filter(annee<=1997)))$r.squared 
  
f2 <- (r2_full - r2_red) / (1 - r2_full)

pwr::pwr.f2.test(u = 1, # 1 seule variable testée (day_light_sav)
                 v = 7805 - 2 - 1,# N - k - 1 : N est la taille de l'échantillon 
                                  # et k est le nombre de variables du modèle
                 f2 = f2, sig.level = 0.05)

## 
##      Multiple regression power calculation 
## 
##               u = 1
##               v = 7802
##              f2 = 0.0004968297
##       sig.level = 0.05
##           power = 0.503581

Sur la période (reg4, N = 7805), l’introduction de la variable test fait passer le \(R^2\) de 0,25 % (avec le seul contrôle lundi) à 0,30 %, soit un apport marginal de variance expliquée de seulement \(\Delta R^2 = 0,05\ \%\) (\(f^2 \approx 0,0005\)). En calculant la puissance du test \(F\) associé à cette amélioration marginale du modèle (pwr.f2.test avec u = 1, la variable test, et v = 7802 pour les degrés de liberté), on obtient une puissance d’environ 50,3 %. Le constat est sans appel. Malgré un échantillon en apparence très vaste (7 805 jours de cotation), la faible fréquence de l’évènement testé fait chuter la puissance du modèle multivarié au niveau d’un simple pile ou face.

Résultat de tous cela, on ne peut pas dire grand chose sur l’hypothèse d’un effet changement d’heure. “Absence de preuve n’est pas preuve d’absence.” Les données ne permettent tout simplement pas de répondre à la question. Le puissance des tests possibles est trop faible (moins de 80%). Un résultat non significatif après corrections rigoureuses signifie simplement que vos données ne permettent pas de répondre à la question. À l’inverse, un résultat hautement significatif obtenu avec une faible puissance statistique est le symptôme d’une anomalie méthodologique, pas d’une découverte économique.

Mais somme nous complètement dépourvu face à ce type de configuration?

Pour éviter de travailler avec un groupe de contrôle géant et bruyant (13 987 jours ordinaires), la solution classique consiste à construire un groupe de contrôle ciblé par appariement / propensity score matching ou coarsened exact matching. On pourrait ici par exemple sélectionner pour chaque lundi de changement d’heure (\(N_1 =\) 109) seulement 2 ou 3 lundis “témoins” aux caractéristiques comparables (même mois, même niveau de volatilité, même contexte macroéconomique). On passe ainsi d’une comparaison de 109 cas traités contre cas 13 987 non traités à une comparaison équilibrée 109 contre 327 cas de contrôle non traités. On élimine le bruit de fond des 13 600 autres jours qui n’apportent que de la variance résiduelle et réduisent le ratio signal/bruit.