Introduction

Nous étudions les données du taux d’investissement aux Etats-Unis afin de faire de la prévision.

Le taux d’investissement aux États-Unis est généralement mesuré par la formation brute de capital fixe (FBCF). La FBCF mesure la valeur totale des biens durables, tels que les équipements, les bâtiments et les machines, achetés par les entreprises, les gouvernements et les particuliers sur une période donnée.

Les données sur la FBCF sont publiées chaque trimestre par le Bureau of Economic Analysis (BEA) du département du Commerce des États-Unis. Le BEA recueille des données auprès des entreprises, des gouvernements et des particuliers pour calculer la FBCF.

Les données que nous avons récupérées constituent donc une série temporelle. Car une série temporelle se caractérise par un ensemble d’observations mesurées sur des temps distincts à des dates d’observations souvent équidistantes. Dans notre cas nous avons une série annuelle.

L’objectif de ce projet est de prévoir les taux d’investissement des Etats-Unis. Pour faire de la prévision, nous aurons besoin que,le Processus Générateur de Données (PGD), i.e. l’équation qui a générénos données soit stationnaire.

Il existe plusieurs types de PGD comme on peut le voir sur la figure ci-dessous :

Un PGD est stationnaire si :

• si son espérance et son autocovaraince ne dépendent pas du temps
• si sa variance est finie

Un PGD est TS(Trend Stationary) si :

• son espérance dépend du temps
• sa variance ne dépend pas du temps

Un PGD est DS (Difference Stationary) si :

• on espérance dépend du temps
• sa variance dépend du temps

Pour identifier la nature de notre PGD nous réaliserons des tests de racine unitaire. Si nous concluons qu’il n’est pas stationnaire nous devrons le transformer en un processus stationnaire avant de faire la prévision.

I - Visualisation des données

Comme nous pouvons le voir, nous avons les taux d’investissement des Etats-Unis de 1972 jusqu’à 2021.
Nous avons donc 50 observations.

Résumé des données:

library(ggplot2)

# Créer un data frame avec les données de votre série temporelle
df <- data.frame(
  année = c(1972:2021),
  taux_investissement = Y)

ggplot(df, aes(x = année, y = taux_investissement)) +
  geom_bar(stat = "identity", fill = "steelblue4") +
  
  labs(title = "Taux d'investissement annuel de 1972 à 2021 au USA",
       x = "Année",
       y = "Taux d'investissement (%)") +
  theme_minimal()

# Convertir la colonne de la série temporelle en une série temporelle de R
ts_data <- Y
# Calculer les statistiques descriptives
mean_value <- mean(ts_data)
median_value <- median(ts_data)
sd_value <- sd(ts_data)
min_value <- min(ts_data)
max_value <- max(ts_data)
q1 <- quantile(ts_data, probs = 0.25)
q2 <- quantile(ts_data, probs = 0.5)
q3 <- quantile(ts_data, probs = 0.75)
skewness <- moments::skewness(ts_data)
kurtosis <- moments::kurtosis(ts_data)

# Créer un data frame avec les statistiques descriptives
stats_df <- data.frame(
  Statistique = c("Moyenne", "Médiane", "Écart-type", "Minimum", "Maximum", "Q1", "Q2", "Q3", "Skewness", "Kurtosis"),
  Valeur = c(mean_value, median_value, sd_value, min_value, max_value, q1, q2, q3, skewness, kurtosis)
)

# Afficher le tableau
knitr::kable(stats_df, caption = "Statistiques descriptives de la série temporelle du taux d'investissement")

Statistiques descriptives de la série temporelle du taux d’investissement
Statistique	Valeur
Moyenne	21.6138966
Médiane	21.6027209
Écart-type	1.4211751
Minimum	18.3139378
Maximum	24.4266334
Q1	20.6483143
Q2	21.6027209
Q3	22.5329075
Skewness	-0.2595319
Kurtosis	2.4884793

Nos taux d’investissement oscillent entre 18.31% et 24.43% avec une moyenne à 21.61%et une médiane à 21.60%.

Le coefficient de skewness est de -0,26, ce qui indique que la distribution des valeurs de la série est légèrement asymétrique vers la gauche.
Le coefficient d’aplatissement est de 2,49, ce qui indique que la distribution des valeurs de la série est plus aplatie (moins pointue) que la distribution normale.

• Chronogramme

Le chronogramme est une représentation graphique où l’on place notre série en ordonnée et le temps en abscisse. Il nous permettra de voir :
• l’évolution de la série dans le temps
• si la série a une tendance (croissante ou décroissante)
• si la série fluctue autour d’une valeur nulle ou non nulle
• si il y a des effets saisonniers dons notre série

plot(Y,main= "Taux d'investissement annuels des Etats−Unis de 1972 à 2021" ,xlab="années" ,ylab="tx d'investissement en %",type="l")
points(Y, col = "grey38", pch = 19)

box(lwd = 3, col = "black", bty = "l")
abline(h = mean(Y), col = "grey", lty = 3)
text(x = 1981, y = mean(Y), labels = "Valeur moyenne", pos = 1, col = "grey",cex=0.75)

abline(v = c("1975"), col = "red", lty = 3)
text(x = 1973, y = 24, labels = "1975", pos = 1, col = "firebrick4",cex=0.75)


text(x = 1982, y = 22.2, labels = "1982-1983", pos = 1, col = "firebrick4",cex=0.75)
symbols(x = 1982.5, 
        y = 22.5, 
        circles = rep(1.8, 1), 
        add = TRUE,
        inches = FALSE,
        fg = "red",
        lwd = 2)
abline(v = c("1992"), col = "red", lty = 3)
text(x = 1994, y = 20, labels = "1992 ", pos = 1, col = "firebrick4",cex=0.75)


text(x = 2003, y = 21, labels = "2002-2003", pos = 1, col = "firebrick4",cex=0.75)
symbols(x = 2002.5, 
        y = 21.5, 
        circles = rep(2, 1), 
        add = TRUE,
        inches = FALSE,
        fg = "red",
        lwd = 2)
           
abline(v = c("2010"), col = "red", lty = 3)
text(x = 2012, y = 18.6, labels = "2010 ", pos = 1, col = "firebrick4",cex=0.75)

NB : La droite grise horizontale représente la moyenne de nos données.

Il n’y a pas de tendance globale : ni croissante ni décroissante. La série ne semble pas fluctuer autour d’une valeur nulle ni autour d’une valeur non nulle. Il ne semble pas y avoir de saisonnalité, les écarts entre les différents pics ne sont pas égaux. Il ne semble pas y avoir d’hétéroscédasticité.

Quelques faits historiques qui pourraient avoir eu une influence sur nos données :

• 1974 : choc pétrolier et dépôt de bilan de la banque allemande Herstatt
• 1979 : 2ème choc pétrolier et mise en place d’une politique monétaire restrictive aux États-Unis
• 1982 : récession brutale aux Etats-Unis
• 1987 : Krach d’octobre 1987
• 1989 : “Junk Bonds”, Bulle spéculative Japonaise
• 1990 : Invasion du Koweit (pétrol)
• 2001 : Krach boursier de 2001-2002 (“éclatement de la bulle internet”), Attentats du 11 septembre
• 2008 : crise des Subprimes
• 2020 : crise liée au Covid-19

Remarquons que souvent une crise n’a pas d’impact directement sur l’année en cours, elle peut se voir quelques années après.

• (ACF) et (PACF)

Fonction d’autocorrélation (ACF) et d’autocorrélation partielle(PACF) sur la série originale.

Avant de réaliser les tests on regarde l’ACF et le PACF de la série originale afin d’avoir une idée sur l’autocorrélation.

acfY=acf(Y,main="ACF sur la série orginale")

pacfY=pacf(Y,main="PACF sur la série orginale")

D’après l’ACF $𝜌(1)$, $𝜌(2)$,$𝜌(3)$,$𝜌(4)$ sont significatifs .Mais d’après le PACF il n’y a que $𝑎(1)$ et $𝑎(2)$ qui sont significatifs donc il semble y avoir de l’autocorrélation jusqu’à l’ordre 2.

II- Tests de racine unitaire

Les tests de racine unitaire permettent d’identifier la nature de notre Processus Générateur de Données (PGD).Ils vont donc nous permettre d’identifier si notre PGD est stationnaire ou non. Si il n’est pas stationnaire, il a une tendance et il est donc TS(Trend Stationnary) ou DS(Difference Stationnary).

A la vue de notre chronogramme, notre série ne semble pas stationnaire. Les conclusions de nos tests devraient alors nous dire que notre Processus Générateur de Données est soit DS soit TS.

Dans les 4 tests que nous allons réaliser nous prendrons un risque de première espèce : 𝛼 = 0, 05.
Le risque de 1ère espèce représente la probabilité de rejeter à tort$𝐻_0$.

La règle de décision est la suivante : si la p-value est inférieure à 𝛼 on rejette l’hypothèse nulle $(H_0)$. La p-value représente donc la probabilité de rejeter à tort l’hypothese nulle.

Nous réaliserons 4 tests :

• le test de Dickey-Fuller (DF)
• le test de Dickey-Fuller augmenté (ADF)
• le test de Zivot et Andrews (ZA)
• le test de Lee et Strazicich (LS)

Dans un premier temps nous réaliserons le test de Dickey-Fuller et le test de Dickey-Fuller augmenté.

Dans ces deux tests il y a trois spécifications :

• la spécification trend
• la spécification drift
• la spécification none

• Test de Dickey-Fuller (DF)

• None

La spécification drift et trend ne sont pas les bonnes(voir annexe) donc je continue avec None.
On estime :

\[Δ𝑋_𝑡 = (𝜌 − 1)𝑋𝑡−1 + 𝜖_𝑡\]

Et on teste :

\[𝐻0∶ ρ-1=0 \quad VS\quad 𝐻𝑎:| ρ |< 1 \]

summary(ur.df(Y,type="none",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.4005 -0.4742  0.1729  0.4801  1.3538 
## 
## Coefficients:
##          Estimate Std. Error t value Pr(>|t|)
## z.lag.1 -0.001168   0.004695  -0.249    0.805
## 
## Residual standard error: 0.7121 on 48 degrees of freedom
## Multiple R-squared:  0.001288,   Adjusted R-squared:  -0.01952 
## F-statistic: 0.06192 on 1 and 48 DF,  p-value: 0.8046
## 
## 
## Value of test-statistic is: -0.2488 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

On teste alors : \[𝐻_0∶ ρ-1=0 \quad VS\quad 𝐻_𝑎:| ρ |< 1 \]

La règle de décision est : si la statistique t calculé pour $(𝜌 − 1)$ est supérieur à la valeur critique donnée à l’intersection de la ligne tau1 et de la colonne 5%, on accepte $𝐻_0$

plot(ur.df(Y,type= "none",lag=0))

D’après l’ACF et le PACF sur nos résidus il y a de l’autocorrélation dans les aléas $(𝜖_𝑡)$ donc les conclusions du test DF ne sont pas valables.

On réalise donc le test ADF.

• Test de Dickey-Fuller augmenté (ADF)

Le test ADF correspond à un test DF avec des variables explicatives en plus qui représentent la variable retardée jusqu’à l’ordre P.

Pour l’instant on garde la spécification None (comme dans le test DF)

On estimera alors :

\[Δ𝑋_𝑡 = (𝜌 − 1)𝑋_𝑡−1 + 𝛽_0 + 𝛽_1𝑡𝑒𝑛𝑑𝑎𝑛𝑐𝑒_𝑡 + 𝜖𝑡 +\sum_{p=1}^{P} 𝛾𝑝Δ𝑋_{t-p}\] Avant de réaliser le test nous devons trouver la valeur de P qui correspond au nombre de retard que l’on doit ajouter à notre régression pour tenir compte de l’autocorrélation dans les aléas.

On commence par calculer la valeur maximale de P avec la formule de Schwert : $12∗ \left( \frac{T}{100} \right)^{0.25}$ (où l’on prend uniquement la partie entière et où T correspond à la taille de notre échantillon.)

T = length(Y)
#taille=50
pmax<-as.integer(12*(T/100)^(0.25))
#pmax=10

On minimise maintenant le MAIC(p), qui est calculé pour les différentes valeurs de p (de 0 à pmax).

library(CADFtest)
summary(CADFtest(Y,criterion="MAIC",type="none",max.lag.y=pmax))

## Augmented DF test 
##                                             ADF test
## t-test statistic:                          -0.342703
## p-value:                                    0.555133
## Max lag of the diff. dependent variable:    2.000000
## 
## Call:
## dynlm(formula = formula(model), start = obs.1, end = obs.T)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.86565 -0.29418  0.07373  0.33936  0.89112 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## L(y, 1)    -0.001389   0.004053  -0.343   0.5551    
## L(d(y), 1)  0.719266   0.152841   4.706 3.68e-05 ***
## L(d(y), 2) -0.352419   0.153219  -2.300   0.0273 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5411 on 36 degrees of freedom
## Multiple R-squared:  0.3846, Adjusted R-squared:  0.3334 
## F-statistic: 11.12 on 2 and 36 DF,  p-value: 0.000174

On trouve “Max lag of the diff.dependent variable” = 2 on doit ajouter au maximum 2 variables explicatives pour prendre en compte l’autocorrélation.

Comme la valeur que l’on trouve est ≠ 0 on garde P=2 et on ne fait pas le BIC.

Dans un premier temps on va tester les coefficients 𝛾 pour savoir combien de variables explicatives que l’on doit rajouter à notre modèle. La règle de décision est : si la valeur absolue de la statistique t de $𝛾$ est supérieure à 1.64, alors $𝛾$ est significatif. On commence donc par tester $𝛾_2$ :

summary(ur.df(Y,type="none",lag=2))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.91617 -0.35172  0.06866  0.40465  1.12578 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## z.lag.1     -0.0008118  0.0039853  -0.204   0.8395    
## z.diff.lag1  0.6913182  0.1391323   4.969 1.07e-05 ***
## z.diff.lag2 -0.3634306  0.1391050  -2.613   0.0122 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5915 on 44 degrees of freedom
## Multiple R-squared:  0.3603, Adjusted R-squared:  0.3167 
## F-statistic: 8.262 on 3 and 44 DF,  p-value: 0.0001804
## 
## 
## Value of test-statistic is: -0.2037 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

En valeur absolue, la statistique t de $𝛾_2$ = -2.613 est supérieur à 1,64. $𝛾_2$ est donc significatif.

et on teste :

\[𝐻_0∶ ρ-1=0 \quad VS\quad 𝐻_𝑎:| ρ |< 1 \]

La règle de décision est : si la statistique t calculé pour $(𝜌 − 1)$ est supérieur à la valeur critique donnée à l’intersection de la ligne tau1 et de la colonne 5%, on accepte H0.

Dans notre cas, la statistique t calculé vaut -0,2037 et est supérieur à -1.95. Donc on on accepte $𝐻_0$ et on conclut que notre PGD est DS. Donc selon le test ADF il y a une présence de racine unitaire.

• Test de Zivot et Andrews (ZA)

Nous allons réaliser le test de Zivot et Andrews et le test de Lee et Strazicich qui permettent de prendre en compte les changements structurels car nous savons que la dynamique d’un processus économique particulier peut être affectée par des changements structurels qui engendrent une instabilité dans le temps de cette dynamique . Nous allons donc introduire des variables pour modéliser ces changements.

Soit $𝑇_𝐵$ la date du changement structurel,$𝐷𝑈_𝑡$(= $𝛽_0$) la variable à ajouter pour modéliser un changement dans le niveau de la partie déterministe de la série et $𝐷𝑇_𝑡$ (= $𝛽_1$) la variable à ajouter pour modéliser un changement dans la pente de la partie déterministe de la série.

\[ \ DU_t = \left\{ \begin{array}{ll} \ 1 & \mbox{si } \ t ≥T_b+1 \\ \ 0 & \mbox{sinon.} \end{array} \right. \] \[ \ DT_t = \left\{ \begin{array}{ll} \ t-T_b & \mbox{si } \ t ≥T_b+1 \\ \ 0 & \mbox{sinon.} \end{array} \right. \]

Il y a deux spécifications dans le test de Zivot et Andrews :

• la spécification crash (où le niveau de la série est touché) :

\[ 𝑦𝑡 = 𝛽0 + 𝛽1𝑡 + 𝜌𝑦_{𝑡−1} + 𝛾𝐷𝑈_𝑡(𝑇_𝐵) +\sum_{j=1}^{P} 𝛾𝑗Δ𝑦_{y−𝑗} + 𝜖_t \]

• la spécification both (où le niveau de la série et son taux de croissance sont touchés) :

\[ 𝑦𝑡 = 𝛽0 + 𝛽1𝑡 + 𝜌𝑦_{𝑡−1}+ 𝛾_1𝐷𝑈_𝑡(𝑇_𝐵) + 𝛾_2𝐷𝑇_𝑡(𝑇_𝐵) +\sum_{j=1}^{P} 𝛾𝑗Δ𝑦_{y−𝑗} + 𝜖_t \]

avec $𝜖_𝑡\ i.i.d (0, 𝜎2_𝜖)$

• Intercept ou “Crash”

La spécification “both” n’est pas appropriée et “trend” est exclue donc je fais “intercept”.

summary(ur.za(Y, model="intercept",lag=pmax-3))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.23138 -0.15416  0.04002  0.19797  0.73466 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 10.62483    3.03665   3.499 0.001437 ** 
## y.l1         0.50086    0.13060   3.835 0.000577 ***
## trend       -0.03281    0.01115  -2.943 0.006114 ** 
## y.dl1        0.76760    0.15680   4.895  2.9e-05 ***
## y.dl2        0.24125    0.18834   1.281 0.209717    
## y.dl3        0.23152    0.18944   1.222 0.230862    
## y.dl4        0.05954    0.17369   0.343 0.734065    
## y.dl5        0.24145    0.16397   1.473 0.150943    
## y.dl6        0.08707    0.16045   0.543 0.591254    
## y.dl7        0.29979    0.14906   2.011 0.053077 .  
## du           1.13284    0.58541   1.935 0.062143 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4863 on 31 degrees of freedom
##   (8 observations effacées parce que manquantes)
## Multiple R-squared:  0.9122, Adjusted R-squared:  0.8839 
## F-statistic: 32.22 on 10 and 31 DF,  p-value: 1.287e-13
## 
## 
## Teststatistic: -3.8218 
## Critical values: 0.01= -5.34 0.05= -4.8 0.1= -4.58 
## 
## Potential break point at position: 9

Je passe donc au modèle intercept avec lequel j’optient une p-value de 0.062 pour $𝛿_2(dt)$ ,légèrement supérieur à la valeur critique 5% , je suppose que j’accepte le modèle intercept. (sachant que le modèle trend est exclue)

La statistique t calculée = -3.82 est supérieur à -4.8 donc on accepte - 𝐻0.

$𝛽_1$ est significatif car sa p-value = 0.006114 < 0,05.

Nous concluons donc avec le test ZA que notre PGD est DS sans changement structurel, ce qui peut être expliqué par le fait que notre série est exposée à des crises brutales.

En effet, le test de ZA est approprié pour détecter des crises graduelles et non-abruptes.

plot(ur.za(Y,model="intercept",lag=pmax-3))

Pas de date de rupture d’après ZA

• Test de Lee et Strazicich (LS)

Le modèle est le suivant : \[𝑦_𝑡 = 𝛾′𝑍_𝑡 + 𝑒_𝑡\] avec \[𝑒_𝑡 = 𝛽𝑒_{𝑡−1} + 𝜖_𝑡, 𝜖 ∼ 𝒩(0, 𝜎2)\]

et Z la matrice des variables exogènes.

Il y a deux spécifications dans le test de Lee et Strazicich :

• la spécification crash :

\[𝑍_𝑡 = [𝜄, 𝑡𝑒𝑛𝑑𝑎𝑛𝑐𝑒, 𝐷𝑈_{1𝑡}, 𝐷𝑈_{2𝑡}]′\]

• la spécification break : \[𝑍_𝑡 = [𝜄, 𝑡𝑒𝑛𝑑𝑎𝑛𝑐𝑒, 𝐷𝑈_{1𝑡}, 𝐷𝑈_{2𝑡}, 𝐷T_{1𝑡}, 𝐷T_{2𝑡}]′\]

Comme nous avons moins de 50 observations il faut utiliser “bootstrap” dans le test de Lee et Strazicich.

Pour “crash”, nous testons :

\[H_0: y_t=\mu_0+d1B{1t}+d2B{2t}+y{t-1}+v{1t}\]

\[VS\] \[H_a: y_t=\mu_1+\gamma \times trend_t+d1D{1t}+d2D{2t}+v_{2t}\]

Règle de décision : si la valeur critique du 𝜆 estimé (au seuil de 5%) est supérieur à la statistique calculée, alors on rejette $𝐻_0$.

source("C:/Users/Bekim/Desktop/R_USA/LeeStrazicichUnitRootTest.R", encoding = 'UTF-8')
source("C:/Users/Bekim/Desktop/R_USA/LeeStrazicichUnitRootTestParallelization.R", encoding = 'UTF-8')

library(parallel)
library(doSNOW)
myBreaks = 1
myModel = "crash"
myLags <- 5
cl <- makeCluster(max(1, detectCores() - 1))
registerDoSNOW(cl)
myParallel_LS <- ur.ls.bootstrap(y=Y , model = myModel, breaks = myBreaks, lags = myLags, method = "Fixed",pn = 0.1, critval = "bootstrap", print.results = "print")

## [[1]]
## [1] -2.848235
## 
## [1] "First possible structural break at position: 41"
## [1] "The location of the first break - lambda_1: 0.8 , with the number of total observations: 50"
## Critical values - Crash model:
##          1%     5%    10%
## [1,] -4.239 -3.566 -3.211
## [1] "Number of lags used: 5"
## Runtime:
## Time difference of 0.003524748 mins

Le $\lambda$ estimé est 0,8 ce qui donne $T_B=41$ donc la date de rupture correspond à 2013.

On compare donc la statistique = -2.848235 à -3.566. Comme -2.848235 > -3.566, on accepte $𝐻_0$et il y a une présence de racine unitaire.

Donc la conclusion du test LS est que le PGD qui a généré notre série est DS avec une date de rupture

III-PGD de la série différenciée(2)

Test de ru pour la série différenciée à l’odre 1 dans l’annexe

On reprend les tests de racine unitaire pour définir le PGD de la série différenciée à l’ordre 2.

ydiff=diff(Y)
yd2=diff(ydiff)
plot(yd2, main="Chronogramme sur la série différenciée(2)",ylab="Diff(2) Invest USA",xlab="année")
fit2 = lm(yd2 ~ c(1974:2021))
abline(fit2, col="red")

- Tendance constante comme nous le montre la régression en rouge
- Homoscédasticité (pas de forme d’entonoir)

Les hypothèses, conditions d’acceptation et formules de test sont les mêmes que celles utilisées plus haut dans les tests sur la série du taux épargne.

• Test de Dickey-Fuller

• None

summary(ur.df(yd2,type="none",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.56989 -0.43907 -0.07189  0.38787  1.96649 
## 
## Coefficients:
##         Estimate Std. Error t value Pr(>|t|)    
## z.lag.1  -0.9721     0.1462  -6.649 3.07e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7169 on 46 degrees of freedom
## Multiple R-squared:  0.4901, Adjusted R-squared:  0.479 
## F-statistic: 44.21 on 1 and 46 DF,  p-value: 3.065e-08
## 
## 
## Value of test-statistic is: -6.6492 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

La statistique t calculée pour $(ρ−1)$ = -6.6492 < -1,95 , on rejette $H_0$ et le processus qui a généré nos données est stationnaire.

Cependant, les résultats sont à nouveau valables uniquement si les aléas ϵt ne sont pas autocorrélés.

Testons cela avec l’ACF et le PACF des résidus des régressions DF.

plot(ur.df(yd2,type="none",lag=0))

Nous constatons, dans le $PACF$, qu’il y a de l’autocorrélation des aléas puisque $a(4)$ est significatif. La conclusion de notre test DF n’est donc pas valide. Nous allons alors faire le test Dickey-Fuller Augmenté puisque celui-ci tient compte de l’autocorrélation.

• Test de ADF

On garde généralement la même spécification que DF pour ADF.

Dans notre cas, la spécification à garder est none que voici :

Pour déterminer le nombre P de variables explicatives à ajouter, on va commencer par calculer pmax en reprenant la formule de Schwert, décrite plus haut dans le test de Zivot et Andrews. Puis nous allons minimiser le critère d’information de Ng et Perron (2001), le $MAIC(p)$, calculé pour différentes valeurs de p allant de 0 à $pmax$.

T_2= length(yd2)#48
pmax3 = as.integer(12*(T_2/100)^(0.25))
summary(CADFtest(yd2, criterion='MAIC', type='none', max.lag.y=pmax3))

## Augmented DF test 
##                                                 ADF test
## t-test statistic:                          -5.757707e+00
## p-value:                                    2.969035e-07
## Max lag of the diff. dependent variable:    0.000000e+00
## 
## Call:
## dynlm(formula = formula(model), start = obs.1, end = obs.T)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.31368 -0.39115 -0.09644  0.33756  2.01578 
## 
## Coefficients:
##         Estimate Std. Error t value Pr(>|t|)    
## L(y, 1)  -0.9352     0.1624  -5.758 2.97e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6457 on 37 degrees of freedom
## Multiple R-squared:  0.4726, Adjusted R-squared:  0.4583 
## F-statistic:    NA on NA and NA DF,  p-value: NA

On trouve “Max lag of the diff.dependent variable” = 0 on fait donc le BIC.

BIc:

Ce test a pour hypothèses :

$H_0 : \gamma_1 = 0$ VS $H_a : \gamma_1 \ne 0$

Le critère est que si la valeur de la statistique t de $\gamma_1$ soit suppérieure au seuil de 1,64 en valeur absolue, alors on rejette $H_0$

summary(ur.df(yd2,type="none",lag=pmax3,selectlag="BIC"))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.61161 -0.36228 -0.05454  0.34075  1.31396 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## z.lag.1      -1.8926     0.3412  -5.546 3.35e-06 ***
## z.diff.lag1   0.8257     0.2673   3.089  0.00399 ** 
## z.diff.lag2   0.5816     0.2145   2.711  0.01044 *  
## z.diff.lag3   0.3582     0.1413   2.535  0.01601 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5915 on 34 degrees of freedom
## Multiple R-squared:  0.5933, Adjusted R-squared:  0.5454 
## F-statistic:  12.4 on 4 and 34 DF,  p-value: 2.528e-06
## 
## 
## Value of test-statistic is: -5.5463 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

Le critère de minimisation de BIC nous donne 0 variable à ajouter.

La statistique t calculée pour $\gamma_3$ = 2.535 > 1,64 , on rejette $H_0$ et le critère de minimisation de $BIC$ nous donne alors 3 variables à ajouter.

On a donc :

summary(ur.df(yd2,type="none",lag=3))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.63538 -0.38352 -0.01749  0.33877  1.37404 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## z.lag.1      -2.0569     0.3041  -6.765 3.99e-08 ***
## z.diff.lag1   0.8764     0.2355   3.722 0.000607 ***
## z.diff.lag2   0.6877     0.1827   3.763 0.000539 ***
## z.diff.lag3   0.4058     0.1315   3.087 0.003667 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5983 on 40 degrees of freedom
## Multiple R-squared:  0.6452, Adjusted R-squared:  0.6097 
## F-statistic: 18.18 on 4 and 40 DF,  p-value: 1.392e-08
## 
## 
## Value of test-statistic is: -6.7647 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

La statistique t de $𝛾_3$ = 3.087 est supérieur à 1,64 . $𝛾_3$ est donc significatif.

Dans ce cas, la statistique t calculé vaut -6.7647 et est inférieure à -1.95. Donc on on rejette $𝐻_0$ et on conclut que le PGD de notre série différenciée à l’odre 2 est stationnaire.

• Test de ZA

• Both

On commence à nouveau par la spécification both et on applique la procédure de Perron.

summary(ur.za(yd2, model="both",lag=pmax3))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.22890 -0.23785 -0.01839  0.33049  0.98710 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept) -0.64626    0.52154  -1.239  0.22727   
## y.l1        -4.58354    1.54813  -2.961  0.00681 **
## trend        0.03095    0.02339   1.323  0.19828   
## y.dl1        4.22137    1.43056   2.951  0.00697 **
## y.dl2        3.76759    1.33504   2.822  0.00943 **
## y.dl3        3.29584    1.22286   2.695  0.01265 * 
## y.dl4        2.59334    1.05848   2.450  0.02196 * 
## y.dl5        1.97270    0.86033   2.293  0.03091 * 
## y.dl6        1.39246    0.64190   2.169  0.04019 * 
## y.dl7        1.00274    0.45810   2.189  0.03857 * 
## y.dl8        0.70513    0.32241   2.187  0.03872 * 
## y.dl9        0.28123    0.19327   1.455  0.15860   
## du          -1.08514    0.53218  -2.039  0.05260 . 
## dt           0.05470    0.04280   1.278  0.21348   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.611 on 24 degrees of freedom
##   (10 observations effacées parce que manquantes)
## Multiple R-squared:  0.4217, Adjusted R-squared:  0.1085 
## F-statistic: 1.346 on 13 and 24 DF,  p-value: 0.255
## 
## 
## Teststatistic: -3.6066 
## Critical values: 0.01= -5.57 0.05= -5.08 0.1= -4.82 
## 
## Potential break point at position: 32

En valeur absolue, la statistique t de $γ_9$ = 1.455 est inférieure à 1,64. $γ_9$ n’est donc pas significatif. On diminue donc la valeur de lag.

summary(ur.za(yd2, model="both",lag=pmax3-8))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.55294 -0.43657 -0.04653  0.35823  1.79182 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)   2.0831     0.8352   2.494   0.0169 *
## y.l1         -0.2889     0.2010  -1.437   0.1585  
## trend        -0.3470     0.1327  -2.615   0.0125 *
## y.dl1         0.2498     0.1430   1.747   0.0882 .
## du            1.1460     0.5306   2.160   0.0368 *
## dt            0.3429     0.1330   2.577   0.0137 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6879 on 40 degrees of freedom
##   (2 observations effacées parce que manquantes)
## Multiple R-squared:  0.1816, Adjusted R-squared:  0.07933 
## F-statistic: 1.775 on 5 and 40 DF,  p-value: 0.14
## 
## 
## Teststatistic: -6.411 
## Critical values: 0.01= -5.57 0.05= -5.08 0.1= -4.82 
## 
## Potential break point at position: 9

Je refais la même opération plusieurs fois(8fois) pour arriver à la statistique t de $γ_1$ = 1.747est supérieur à 1,64. $γ_1$ est donc significatif.

Nous voyons que $δ_1$(du) et $δ_2$(dt) sont significatifs car la p-value de $δ_1$= 0.0368 < 0,05 et la p-value de $δ_2$= 0.0137 < 0,05.

Le modèle both est donc correcte.

$𝛽_1$ est significatif car sa p-value = 0.0125 < 0,05.

De plus, la statistique t calculée = -6.411 est inférieure à -5,08 donc on rejette - $𝐻_0$.

“Potential break point at position: 9”.

La date de rupture est donc : $T_B$=1982

Nous concluons donc avec le test ZA que notre PGD (sur série diff x2) est TS avec changement structurel (mais possiblement DS avec changement structurel). La date de rupture proposé par le test ZA est $T_B$ = 1982 ce qui correspond à une récession brutale aux Etats-Unis.

plot(ur.za(yd2,model="both",lag=1))

• Test de LS

Nous allons passer au test de Lee et Strazicich qui est plus performant et qui sera celui dont on tirera la conclusion sur le PGD de notre série différenciée à l’odre 2.

Le modèle final que l’on va conserver est le modèle break avec 1 rupture (car on avait δ1 et δ2 significatifs dans notre test de ZA).

source("C:/Users/Bekim/Desktop/R_USA/LeeStrazicichUnitRootTest.R", encoding = 'UTF-8')
source("C:/Users/Bekim/Desktop/R_USA/LeeStrazicichUnitRootTestParallelization.R", encoding = 'UTF-8')

library(parallel)
library(doSNOW)
myBreaks = 1
myModel = "break"
myLags <- 4
cl <- makeCluster(max(1, detectCores() - 1))
registerDoSNOW(cl)
myParallel_LS <- ur.ls.bootstrap(y=yd2 , model = myModel, breaks = myBreaks, lags = myLags, method = "Fixed",pn = 0.1, critval = "bootstrap", print.results = "print")

## [[1]]
## [1] -5.562816
## 
## [1] "First possible structural break at position: 8"
## [1] "The location of the first break - lambda_1: 0.2 , with the number of total observations: 48"
## Critical values - Break model:
##      lambda    1%    5%   10%
## [1,]    0.1 -5.11 -4.50 -4.21
## [2,]    0.2 -5.07 -4.47 -4.20
## [3,]    0.3 -5.15 -4.45 -4.18
## [4,]    0.4 -5.05 -4.50 -4.18
## [5,]    0.5 -5.11 -4.51 -4.17
## [1] "Number of lags used: 4"
## Runtime:
## Time difference of 0.003486701 mins

La valeur de la statistique de test est -5.562816.

Le $\lambda$ estimé est 0,2 ce qui donne $T_B=8$ donc la date de rupture correspond à 1981
La valeur critique au seuil de risque de 5% vaut -4.47.
On a -5.562816 < -4.47.
On rejette donc $H_0$, il y n’a pas de racine unitaire.
La conclusion est que le PGD qui a généré la série différenciée à l’ordre 2 est stationnaire.

IV - Autocorrélation

Nous cherchons à savoir si la série différenciée à l’ordre 2 a de l’autocorrélation.

Fonction d’autocorrélation (ACF) et d’autocorrélation partielle (PACF) sur la série différenciée

• ACF

library(forecast)

Acf(yd2, main="ACF sur la série différenciée à l'ordre 2")

• PACF

library(forecast)

Pacf(yd2, main="ACF sur la série différenciée à l'ordre 2")

D’après l’ACF $ρ(4)$ est significatif et d’après le PACF il n’y a que $a(4)$ qui est significatif donc il semble y avoir de l’autocorrélation jusqu’à l’ordre 4.

• Test de Ljung-Box sur la série différenciée à l’ordre 2

Statistique de Ljung-Box :
\[QK = T \times (T+2) \sum{k=1}^{K} \frac{\hat \rho (k)^2}{T-K}\]

On teste les hypothèses suivantes (à l’ordre k) :

• $𝐻_0$∶ $𝜌_1$ = $𝜌_$2 = $𝜌_3$ = … =$𝜌_𝑘$ = 0 : il n’y a pas autocorrélation (jusqu’à l’ordre k).

versus

• $𝐻_𝑎$: au moins un $𝜌_𝑖$ ≠ 0 : il y a de l’autocorrélation

La règle de décision est : si la p-value est inférieure à 0.05 on rejette $𝐻_0$.

Box.test(yd2, lag=1,type="Ljung-Box")

## 
##  Box-Ljung test
## 
## data:  yd2
## X-squared = 0.039995, df = 1, p-value = 0.8415

D’après le test de Ljung Box il y a de l’autocorrélation. On a rejette $𝐻_0$ à l’ordre 1.

V - Estimation du modèle ARMA(p,q)

• EACF

L’EACF permet de déterminer les valeurs de départ de p et de q pour le modèle ARMA(p,q).

library(TSA)
eacf(yd2)

## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 o o o x o o o o o o o  o  o  o 
## 1 o o o x o o o o o o o  o  o  o 
## 2 x o o x o o o o o o o  o  o  o 
## 3 x o o o o o o o o o o  o  o  o 
## 4 o x o o o o o o o o o  o  o  o 
## 5 x x o o x o o o o o o  o  o  o 
## 6 o o o o o o o o o o o  o  o  o 
## 7 x o x o o o o o o o o  o  o  o

L’EACF nous donne (3,1) comme valeurs de départ de p et q.

Recherche du modèle

Nous estimons :

• ARMA(3,1)

library(forecast)
library(lmtest)

reg<-Arima(Y,order=c(3,2,1), include.constant=TRUE)
coeftest(reg)

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value  Pr(>|z|)    
## ar1  0.610652   0.140297   4.3526 1.346e-05 ***
## ar2 -0.197974   0.161266  -1.2276    0.2196    
## ar3 -0.219477   0.141437  -1.5518    0.1207    
## ma1 -1.000000   0.068979 -14.4972 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

On regarde si les coefficients sont significatifs.

On a le coefficient $ar(3)$ non significatif donc je l’enlève.

• ARMA(2,1)

reg2<-Arima(Y,order=c(2,2,1), include.constant=TRUE)
coeftest(reg2)

## 
## z test of coefficients:
## 
##      Estimate Std. Error  z value  Pr(>|z|)    
## ar1  0.691610   0.135153   5.1173   3.1e-07 ***
## ar2 -0.345629   0.133945  -2.5804  0.009869 ** 
## ma1 -0.999999   0.059531 -16.7979 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

En retirant ar(3) tous nos coefficients sont significatifs.
Avec un BIC de :

BIC(reg2)

## [1] 102.0285

Nous allons essayer de trouver un meilleur modèle en regardant l’eacf.

• ARMA(0,2)

reg3<-Arima(Y,order=c(0,2,2), include.constant=TRUE)
coeftest(reg3)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ma1 -0.41882    0.12494 -3.3520 0.0008022 ***
## ma2 -0.58118    0.11581 -5.0182 5.216e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

BIC(reg3)

## [1] 101.928

Avec le modèle ARMA(0,2) on a un BIC très légèrement inférieur , cependant je privilégie le ARMA(2,1).

VI - Etude des aléas

L’étude des aléas permet de vérifier que les aléas du modèle ARMA que nous avons retenus sont des bruits blancs gaussiens.

Un Bruit Blanc (BB) est une séquence de variables aléatoires iid de moyenne nulle, de variance constante et non autocorrélées.

Nous réaliserons quatre tests :

• le test de Jarque-Bera
• le test de Student
• le test de Ljung-Box
• le test d’Engle
Le test de Jarque-Bera nous permet de savoir les aléas sont gaussiens et les autres si ils sont des bruits blancs.
Si les tests ne sont pas concluant on passera au modèle suivant (i.e. celui qui minimise le BIC parmi les modèles restants).

• Test de Jarque-Bera

— la statistique de test s’écrit:

\[JB= \frac{T-k}{6} (S^{2}+\frac{(K-3)}{4}\]

On cherche à savoir si les aléas sont normalement distribués.
On teste :
• $𝐻_0$∶ les $𝜖$ suivent une Loi Normale

• $𝐻_𝑎$∶ les $𝜖$ ne suivent pas une Loi Normale

La règle de décision est : si la p-value est inférieure à 0.05 on rejette $𝐻_0$.
On effectue le test de Jarque-Bera sur les résidus du modèle ARMA() :

library(tseries)
jarque.bera.test(reg2$res)

## 
##  Jarque Bera Test
## 
## data:  reg2$res
## X-squared = 5.6937, df = 2, p-value = 0.05803

La p-value = 0.05803 est supérieur à 0,05 donc on accepte $H_0$. On conclut dont que les aléas sont normalement distribuées i.e. ils sont gaussiens.

• Test de Student (t.test)

La formule de la statistique t est:

\[\frac{|m(e)|}{σ_e}*\sqrt{T} \]

On cherche à savoir si l’espérance des aléas est nulle.
On teste :
• $𝐻_0$∶ $𝐸(𝜖)$ = 0

• $H_𝑎$∶ $𝐸(𝜖)$ ≠ 0

La règle de décision est : si la p-value est inférieure à 0.05 on rejette $𝐻_0$.

t.test(reg2$res)

## 
##  One Sample t-test
## 
## data:  reg2$res
## t = -0.50274, df = 49, p-value = 0.6174
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.2013642  0.1207740
## sample estimates:
##   mean of x 
## -0.04029509

La p-value = 0.6174 est supérieur à 0,05 donc on accepte $H_0$. On conclut dont que l’espérance des aléas est nulle.

• Test de Ljung-Box sur les résidus centrés-réduits

On cherche à savoir si il y a de l’autocorrélation dans les aléas.
On commence par calculés les résidus centrés-réduits :

Sres = (reg2$res-mean(reg2$res))/sd(reg2$res)

On teste les hypothèses suivantes (à l’ordre k) :
• $𝐻_0$∶ 𝜌1 = 𝜌2 = 𝜌3 = … = 𝜌𝑘 = 0 : il n’y a pas d’autocorrélation (jusqu’à l’ordre k).

Versus

• $𝐻_𝑎$: au moins un $𝜌$𝑖 ≠ 0 : il y a de l’autocorrélation
La règle de décision est : si la p-value est inférieure à 0.05 on rejette $𝐻_0$ .

Box.test(Sres, lag=12, type="Ljung-Box")

## 
##  Box-Ljung test
## 
## data:  Sres
## X-squared = 8.6371, df = 12, p-value = 0.7336

La $p-value$ est de 0.7336 ce qui est supérieur à 0,05 donc nous ne rejetons pas $H_0$. Il n’y a pas d’autocorrélation dans les aléas jusqu’a l’ordre 12.

• Test d’Engle sur les résidus centrés - réduits

On cherche à savoir si il existe des clusters de volatilité dans les données.

Le test ARCH(p) consiste à estimer par MCO l’équation suivante :

\[e^{2}_t=a+c_1e^{2}_{t−1}+c_2e^{2}_{t−2}+⋯+c_pe^{2}_{t−p}+erreur_t\]

$e$ les résidus de la régression et $a$ une constante.

Statistique de test:

\[LM=T×R^{2}\] On teste :
• $𝐻_0$∶ $a_1$ = $a_2=$ = $a_3$ = … = $𝛼_𝑝$ = 0 : il n’y a pas d’effets ARCH et pas de clusters de volatilité

Versus

• $𝐻_𝑎$: au moins un $𝛼_𝑖$ ≠ 0 : il y a des effets ARCH et des clusters de volatilité

Si $H_0$ est vrai :

-σ2t=α0
-non existence de cluster de volatilité
-pas d’effet ARCH
-homoscédasticité conditionnelle

Si $H_a$ est vrai :

-σ2t=α0+α1e2t−1
-existence de cluster de volatilité
-présence d’effet ARCH
-hétéroscédasticité conditionnelle
-Sous $H_0$, $LM∼χ2(p)$

La règle de décision est : si la p-value est inférieure à 0.05 on rejette $𝐻_0$.

On effectue le test sur les résidus standardisés.

library(FinTS)
ArchTest(Sres, lag=12)

## 
##  ARCH LM-test; Null hypothesis: no ARCH effects
## 
## data:  Sres
## Chi-squared = 4.7552, df = 12, p-value = 0.9657

La p−valueest de 0.9657 ce qui est supérieur à 0,05 donc on accepte $H_0$ . Il n’y a donc pas d’effet ARCH ni de cluster de volatilité.

On peut alors en conclure que les aléas sont des bruits blancs.

VII - Prévision

Le modèle ARMA(2,1) est celui qui minimise le critère d’information BIC (il vaut 102.0285) c’est donc celui que l’on va utiliser pour la prévision.

BIC(reg2)

## [1] 102.0285

Prévision :

prev = forecast(reg2, h=4, level=0.95)
prev

##      Point Forecast    Lo 95    Hi 95
## 2022       20.97294 19.79903 22.14686
## 2023       20.86341 18.53645 23.19037
## 2024       20.85444 17.64948 24.05941
## 2025       20.87919 17.06277 24.69560

De 2022 à 2024 les valeurs de la prévision passent de 20.97%à 20.85% et de 2024 à 2025 la prévision augmente légérement à 20.87%.

plot(prev, lwd=1, main="Prévision sur quatre ans du taux d'investissement américain", ylab="Taux d'investissement USA", xlab="Années")

Puisque les aléas sont normalement distribués, on peut également définir un intervalle de confiance à 95% pour les prévisions. On a donc, à 5% de chances de se tromper, un taux d’épargne compris :

entre 19.79% et 22.14% en 2022.

entre 18.53% et 23.19% en 2023.

entre 17.64% et 24.05% en 2024.

entre 17.06% et 24.69% en 2025.

L’intervalle de confiance s’élargit avec le temps pour prendre en compte l’ incertitude accrue.

On remarque avec le temps, que il y a d’incertitude dans la prédiction, car il peut y avoir des facteurs imprévus qui influencent la série.

Annexe

Tests de DF non concluants

•Test sur la série originel

On avait commencé par la spécification trend :

• Trend

On estime :

\[Δ𝑋_𝑡 = (𝜌 − 1)𝑋𝑡−1 + 𝛽0 + 𝛽1𝑡𝑒𝑛𝑑𝑎𝑛𝑐𝑒_𝑡 + 𝜖_𝑡\]

Et on teste :

\[𝐻0∶ 𝛽1 = 0 \hspace{0.25cm} et \hspace{0.25cm} ρ-1=0 \quad Vs\quad 𝐻𝑎∶ 𝛽1 ≠ 0 \hspace{0.25cm}et \hspace{0.25cm}| ρ |< 1 \]

La règle de décision est : si la p-value de $𝛽_1$ < $𝛼$ on rejette$H_0$

summary(ur.df(Y,type="trend",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression trend 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.3656 -0.3670  0.1051  0.4255  1.2598 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  3.913309   1.886278   2.075   0.0436 *
## z.lag.1     -0.170842   0.081815  -2.088   0.0423 *
## tt          -0.009327   0.008214  -1.136   0.2620  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6955 on 46 degrees of freedom
## Multiple R-squared:  0.08665,    Adjusted R-squared:  0.04694 
## F-statistic: 2.182 on 2 and 46 DF,  p-value: 0.1243
## 
## 
## Value of test-statistic is: -2.0881 1.4613 2.1821 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -4.15 -3.50 -3.18
## phi2  7.02  5.13  4.31
## phi3  9.31  6.73  5.61

La p-value de $𝛽_1$ = 0.2620 est supérieur à $𝛼$ = 0,05 donc on accepte $H_0$ et $𝛽_1$ n’est pas significatif.

• Drift

La spécification trend n’est pas la bonne je continue avec Drift.
On estime :

\[Δ𝑋_𝑡 = (𝜌 − 1)𝑋𝑡−1 + 𝛽0 + 𝜖_𝑡\]

Et on teste :

\[𝐻_0∶ 𝛽_0 = 0 \hspace{0.25cm} et \hspace{0.25cm} ρ-1=0 \quad Vs\quad 𝐻_𝑎∶ 𝛽_0 ≠ 0 \hspace{0.25cm}et \hspace{0.25cm}| ρ |< 1 \]

summary(ur.df(Y,type="drift",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression drift 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.45892 -0.39576  0.06371  0.37661  1.39704 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)  2.63937    1.52105   1.735   0.0893 .
## z.lag.1     -0.12271    0.07019  -1.748   0.0870 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6977 on 47 degrees of freedom
## Multiple R-squared:  0.06105,    Adjusted R-squared:  0.04107 
## F-statistic: 3.056 on 1 and 47 DF,  p-value: 0.08697
## 
## 
## Value of test-statistic is: -1.7481 1.5378 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.58 -2.93 -2.60
## phi1  7.06  4.86  3.94

La p-value de $𝛽_0$ = 0.0893 est supérieur à $𝛼$ = 0,05 donc on accepte $H_0$ et $𝛽_0$ n’est pas significatif.

Test de ZA non concluant

• Both

On commencera par le modèle both car il est plus complet.
On testera alors :
• $𝐻_0∶ 𝜌 = 1$, processus DS sans changement structurel

Versus

• $𝐻_𝑎∶ |𝜌| < 1$, processus TS avec un unique changement structurel qui survient à la date $𝑇_𝐵$.

Règle de décision : si la statistique t calculé qui dépend de la date de rupture est supérieure à la valeur critique à 5%, on est enclin à accepter $𝐻_0$.

On prend une valeur de pmax pour le lag puis on diminue le nombre de γ de tel sorte que le dernier ait une statistique t en valeur absolue > 1.64 .

summary(ur.za(Y, model="both",lag=pmax))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.15024 -0.17848  0.06542  0.19871  0.51427 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 25.89381    5.01529   5.163 2.75e-05 ***
## y.l1        -0.06586    0.20592  -0.320 0.751858    
## trend       -0.10757    0.02828  -3.804 0.000864 ***
## y.dl1        1.23695    0.18326   6.750 5.56e-07 ***
## y.dl2        0.64514    0.22261   2.898 0.007896 ** 
## y.dl3        0.65935    0.21114   3.123 0.004626 ** 
## y.dl4        0.44917    0.19601   2.292 0.031010 *  
## y.dl5        0.49467    0.17415   2.841 0.009035 ** 
## y.dl6        0.43539    0.17221   2.528 0.018454 *  
## y.dl7        0.38785    0.15832   2.450 0.021966 *  
## y.dl8        0.24413    0.15893   1.536 0.137588    
## y.dl9        0.06901    0.14760   0.468 0.644334    
## y.dl10       0.21968    0.14834   1.481 0.151648    
## du           1.13772    0.41631   2.733 0.011597 *  
## dt          -0.01442    0.02459  -0.586 0.563134    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4218 on 24 degrees of freedom
##   (11 observations effacées parce que manquantes)
## Multiple R-squared:  0.9425, Adjusted R-squared:  0.9089 
## F-statistic: 28.08 on 14 and 24 DF,  p-value: 1.765e-11
## 
## 
## Teststatistic: -5.176 
## Critical values: 0.01= -5.57 0.05= -5.08 0.1= -4.82 
## 
## Potential break point at position: 31

En valeur absolue, la statistique t de $𝛾_{10}$ = 1,481 est inférieure à 1,64. $𝛾_{10}$ n’est donc pas significatif. On diminue donc la valeur de lag.

summary(ur.za(Y, model="both",lag=pmax-3))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.23138 -0.15416  0.04002  0.19797  0.73466 
## 
## Coefficients: (1 not defined because of singularities)
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 10.62483    3.03665   3.499 0.001437 ** 
## y.l1         0.50086    0.13060   3.835 0.000577 ***
## trend       -0.03281    0.01115  -2.943 0.006114 ** 
## y.dl1        0.76760    0.15680   4.895  2.9e-05 ***
## y.dl2        0.24125    0.18834   1.281 0.209717    
## y.dl3        0.23152    0.18944   1.222 0.230862    
## y.dl4        0.05954    0.17369   0.343 0.734065    
## y.dl5        0.24145    0.16397   1.473 0.150943    
## y.dl6        0.08707    0.16045   0.543 0.591254    
## y.dl7        0.29979    0.14906   2.011 0.053077 .  
## du           1.13284    0.58541   1.935 0.062143 .  
## dt                NA         NA      NA       NA    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4863 on 31 degrees of freedom
##   (8 observations effacées parce que manquantes)
## Multiple R-squared:  0.9122, Adjusted R-squared:  0.8839 
## F-statistic: 32.22 on 10 and 31 DF,  p-value: 1.287e-13
## 
## 
## Teststatistic: -3.8218 
## Critical values: 0.01= -5.57 0.05= -5.08 0.1= -4.82 
## 
## Potential break point at position: 9

Je refais la même opération plusieurs fois(3fois) pour arriver à la statistique t de $𝛾_7$ = 2.011 est supérieur à 1,64. $𝛾_7$est donc significatif.

Dans notre sortie, la valeur de “dt” est “NA”, ce qui signifie que le test n’a pas détecté de rupture structurelle dans la série temporelle. Cela peut être interprété comme une indication que la série est stationnaire ou qu’elle possède une racine unitaire sans rupture structurelle.

Nous voyons que $𝛿_1(du)$ et $𝛿_2(dt)$ ne sont pas significatif car la p-value de 𝛿1 = 0,062 > 0,05 et la p-value de $𝛿_2$=NA.
Le modèle both n’est donc pas approprié.

LS non concluant

• Break

On commence par un lag égal à 5 et la méthode break.

• $𝐻_0 ∶ 𝑦_𝑡 = 𝜇_0 + 𝑑_1𝐵_{1𝑡} + 𝑑_2𝐵_{2𝑡} + 𝑑_3𝐷_{1𝑡} + 𝑑_4𝐷_{2𝑡} + 𝑦_{𝑡−1} + 𝑣_{1𝑡}$

• $𝐻_𝑎 ∶ 𝑦_𝑡 = 𝜇_1 + 𝛾 × 𝑡𝑟𝑒𝑛𝑑_𝑡 +𝑑_1D_{1𝑡} + 𝑑_2D_{2𝑡} + 𝑑_3𝐷T_{1𝑡} + 𝑑_4𝐷T_{2𝑡} + 𝑣_{2t}$

On teste alors :

• $𝐻_0$: présence de racine unitaire dans le PGD

• $H_𝑎$: pas de présence de racine unitaire dans le PGD

Règle de décision : si la valeur critique du 𝜆 estimé (au seuil de 5%) est supérieur à la statistique calculée, alors on rejette $𝐻_0$.

TESTS-Série diffénrenciée à l’odre 1

Les tests de racine unitaire nous on mené à la conclusion d’un PGD DS. Nous allons donc différencier notre série une première fois.

•Différenciation

ydiff=diff(Y)
plot(ydiff, main="Chronogramme sur la série différenciée",ylab="Diff Invest USA",xlab="année")
fit1 = lm(ydiff ~ c(1973:2021))
abline(fit1, col="green")

- Tendance constante comme nous le montre la régression en vert
- Homoscédasticité

Les hypothèses, conditions d’acceptation et formules de test sont les même que celles utilisées plus haut dans les tests sur la série du taux épargne.

• Test de DF

• Trend

summary(ur.df(ydiff,type="trend",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression trend 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.87313 -0.34202  0.09263  0.43851  1.14056 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0298777  0.1862583  -0.160 0.873276    
## z.lag.1     -0.5006169  0.1287174  -3.889 0.000329 ***
## tt           0.0005769  0.0066176   0.087 0.930918    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6352 on 45 degrees of freedom
## Multiple R-squared:  0.2517, Adjusted R-squared:  0.2184 
## F-statistic: 7.568 on 2 and 45 DF,  p-value: 0.001468
## 
## 
## Value of test-statistic is: -3.8893 5.0494 7.5682 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -4.15 -3.50 -3.18
## phi2  7.02  5.13  4.31
## phi3  9.31  6.73  5.61

La p-value de $𝛽_1$ = 0.93 est supérieur à $𝛼$ = 0,05 donc on accepte $H_0$ et $𝛽_1$ n’est pas significatif.

• Drift

summary(ur.df(ydiff,type="drift",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression drift 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.86653 -0.34891  0.08805  0.43750  1.13277 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.01574    0.09069  -0.174 0.862947    
## z.lag.1     -0.50065    0.12732  -3.932 0.000282 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6283 on 46 degrees of freedom
## Multiple R-squared:  0.2516, Adjusted R-squared:  0.2353 
## F-statistic: 15.46 on 1 and 46 DF,  p-value: 0.0002815
## 
## 
## Value of test-statistic is: -3.9322 7.7372 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.58 -2.93 -2.60
## phi1  7.06  4.86  3.94

La p-value de $𝛽_0$ = 0.862947 est supérieur à 𝛼 = 0,05 donc on accepte $H_0$ et $𝛽_0$ n’est pas significatif.

• None

On commence à nouveau par la spécification trend. Après l’avoir testée et rejetée, on teste la spécification drift qui elle aussi est rejetée (cf. annexe). On passe alors à la spécification none.

On estime : $\Delta X_t=(\rho-1)X_{t-1}+\epsilon_t$

Et on teste :

$H_0:\rho-1=0$ VS $H_a:|\rho|<1$

Si la statistique t calculée pour $(\rho-1)$ est inférieure à la valeur critique donnée à l’intersection de la ligne tau1 et de la colonne 5%, alors on rejette $H_0$.

Dans notre cas :

summary(ur.df(ydiff,type="none",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.88189 -0.36493  0.07198  0.42182  1.11706 
## 
## Coefficients:
##         Estimate Std. Error t value Pr(>|t|)    
## z.lag.1  -0.5003     0.1260  -3.971 0.000244 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6217 on 47 degrees of freedom
## Multiple R-squared:  0.2512, Adjusted R-squared:  0.2353 
## F-statistic: 15.77 on 1 and 47 DF,  p-value: 0.000244
## 
## 
## Value of test-statistic is: -3.9711 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

La statistique t calculée pour $(\rho-1)$ = -3.971 < -1,95 , on rejette $H_0$ et le processus qui a généré nos données est stationnaire.

Cependant, les résultats sont à nouveau valables uniquement si les aléas $\epsilon_t$ ne sont pas autocorrélés.

Testons cela avec l’ACF et le PACF des résidus des régressions DF.

plot(ur.df(ydiff,type="none",lag=0))

Nous constatons, dans le PACF, qu’il y a de l’autocorrélation des aléas puisque $a(4)$ et$a(13)$ sont significatif. La conclusion de notre test DF n’est donc pas valide. Nous allons alors faire le test Dickey-Fuller Augmenté puisque celui-ci tient compte de l’autocorrélation.

• Test de ADF

Le test de Dickey-Fuller Augmenté (ADF) est un test de Dickey-Fuller (DF) avec des variables explicatives en plus qui sont la variable dépendante retardée jusqu’à l’ordre P.
P étant le nombre de retards que nous ajoutons dans les régressions pour tenir compte de l’autocorrélation des aléas.
On garde généralement la même spécification que DF pour ADF.

Dans notre cas, la spécification à garder est none que voici :

Pour déterminer le nombre P de variables explicatives à ajouter, on va commencer par calculer $pmax$ en reprenant la formule de Schwert, décrite plus haut dans le test de Zivot et Andrews. Puis nous allons minimiser le critère d’information de Ng et Perron (2001), le $MAIC(p)$, calculé pour différentes valeurs de p allant de 0 à $pmax$.

T_dep= length(ydiff)

pmax2 = as.integer(12*(T_dep/100)^(0.25))
summary(CADFtest(ydiff, criterion='MAIC', type='none', max.lag.y=pmax2))

## Augmented DF test 
##                                                ADF test
## t-test statistic:                          -3.288869517
## p-value:                                    0.001647139
## Max lag of the diff. dependent variable:    0.000000000
## 
## Call:
## dynlm(formula = formula(model), start = obs.1, end = obs.T)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.83001 -0.33223  0.04011  0.40366  1.12206 
## 
## Coefficients:
##         Estimate Std. Error t value Pr(>|t|)   
## L(y, 1)  -0.4526     0.1376  -3.289  0.00165 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5692 on 37 degrees of freedom
## Multiple R-squared:  0.2262, Adjusted R-squared:  0.2053 
## F-statistic:    NA on NA and NA DF,  p-value: NA

On trouve “Max lag of the diff.dependent variable” = 0 on fait donc le $BIC$.

• Méthode 1:

summary(CADFtest(ydiff,criterion="BIC",type="none",max.lag.y=pmax2))

## Augmented DF test 
##                                             ADF test
## t-test statistic:                          -4.183801
## p-value:                                    0.000100
## Max lag of the diff. dependent variable:    1.000000
## 
## Call:
## dynlm(formula = formula(model), start = obs.1, end = obs.T)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.87015 -0.34936  0.02842  0.30840  0.92455 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## L(y, 1)     -0.6127     0.1465  -4.184   0.0001 ***
## L(d(y), 1)   0.3599     0.1524   2.362   0.0237 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.537 on 36 degrees of freedom
## Multiple R-squared:   0.33,  Adjusted R-squared:  0.2928 
## F-statistic:    NA on NA and NA DF,  p-value: NA

• Méthode 2:

summary(ur.df(ydiff,type="none",lag=pmax,selectlag="BIC"))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.87015 -0.34936  0.02842  0.30840  0.92455 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1     -0.6127     0.1465  -4.184 0.000176 ***
## z.diff.lag   0.3599     0.1524   2.362 0.023707 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.537 on 36 degrees of freedom
## Multiple R-squared:   0.33,  Adjusted R-squared:  0.2928 
## F-statistic: 8.867 on 2 and 36 DF,  p-value: 0.0007395
## 
## 
## Value of test-statistic is: -4.1838 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61

Dans la méthode 1 on trouve le lag dans “Max lag of the diff. dependent variable” et dans la méthode 2 on voit la statistique t calculée pour $\gamma_1$ = 2.362 > 1,64 , on rejette $H_0$ et le critère de minimisation de BIC nous donne alors 1 variable à ajouter.

Comme la valeur que l’on trouve est ≠ 0 on prend P=1.(et dernière valeur significatif )

summary(ur.df(ydiff,type="none",lag=1))

Dans notre cas, la statistique t calculé vaut -4.9031 et est inférieure à -1.95. Donc on on rejette 𝐻0 et on conclut que le PGD de notre série différenciée est stationnaire.

On va maintenant passer aux tests prenant en compte des changements structurels, en commençant par Zivot et Andrews.

•Test de ZA

On commence à nouveau par la spécification both et on applique la procédure de Perron.

summary(ur.za(ydiff, model="both",lag=pmax2))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.24707 -0.28615  0.02902  0.24724  0.82136 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)  
## (Intercept) -0.351351   0.454427  -0.773    0.447  
## y.l1        -0.575062   0.591847  -0.972    0.341  
## trend        0.012958   0.018459   0.702    0.490  
## y.dl1        1.073460   0.533392   2.013    0.056 .
## y.dl2        0.781625   0.476593   1.640    0.115  
## y.dl3        0.660763   0.459597   1.438    0.164  
## y.dl4        0.382851   0.422280   0.907    0.374  
## y.dl5        0.382812   0.371714   1.030    0.314  
## y.dl6        0.278136   0.333184   0.835    0.412  
## y.dl7        0.252551   0.280516   0.900    0.377  
## y.dl8        0.188034   0.231244   0.813    0.424  
## y.dl9       -0.003125   0.203874  -0.015    0.988  
## y.dl10       0.047242   0.185050   0.255    0.801  
## du          -0.848585   0.528820  -1.605    0.122  
## dt           0.064864   0.050562   1.283    0.212  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5857 on 23 degrees of freedom
##   (11 observations effacées parce que manquantes)
## Multiple R-squared:  0.5379, Adjusted R-squared:  0.2567 
## F-statistic: 1.913 on 14 and 23 DF,  p-value: 0.08114
## 
## 
## Teststatistic: -2.6613 
## Critical values: 0.01= -5.57 0.05= -5.08 0.1= -4.82 
## 
## Potential break point at position: 35

En valeur absolue, la statistique t de $γ_10$ = 0.255 est inférieure à 1,64.$γ_10$ n’est donc pas significatif. On diminue donc la valeur de lag.

summary(ur.za(ydiff, model="both",lag=pmax2-8))

## 
## ################################ 
## # Zivot-Andrews Unit Root Test # 
## ################################ 
## 
## 
## Call:
## lm(formula = testmat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.97555 -0.20640  0.06932  0.26356  1.00278 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  1.08020    0.36879   2.929  0.00565 **
## y.l1         0.05496    0.15991   0.344  0.73294   
## trend       -0.09469    0.02990  -3.166  0.00299 **
## y.dl1        0.42835    0.13376   3.202  0.00271 **
## y.dl2        0.28984    0.13680   2.119  0.04053 * 
## du           0.78688    0.32721   2.405  0.02103 * 
## dt           0.09129    0.03216   2.838  0.00716 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5337 on 39 degrees of freedom
##   (3 observations effacées parce que manquantes)
## Multiple R-squared:  0.512,  Adjusted R-squared:  0.4369 
## F-statistic: 6.819 on 6 and 39 DF,  p-value: 5.325e-05
## 
## 
## Teststatistic: -5.9098 
## Critical values: 0.01= -5.57 0.05= -5.08 0.1= -4.82 
## 
## Potential break point at position: 19

Je refais la même opération plusieurs fois(8fois) pour arriver à la statistique t de $γ_2$ = 2.119 est supérieur à 1,64. $γ_2$ est donc significatif.

La $p-value$ de $\delta_1$ et $\delta_2$ sont significatif car leur p-value(0.02103 et 0.00716) sont < 0.05 .

De plus, On a pour $\beta_1$, une p-value 0.00299 < 0,05 : $\beta_1$ est significatif.

La statistique de test est -5.9098 < -5.08 la valeur critique à 5% . On rejette donc $H_0$ :$|\rho| < 1$.
En ce qui concerne la date de rupture, on a “Potential break point at position: 19”. La date de rupture est donc : $T_B=1991$

Nous concluons donc avec le test ZA que notre PGD (sur série diff) est TS avec changement structurel (mais possiblement DS avec changement structurel). La date de rupture proposée par le test ZA est $𝑇_𝐵$ = 1991 ce qui correspondrait à l’invasion du Koweit.

plot(ur.za(ydiff,model="both",lag=pmax2-8))

• Test de LS

Nous allons passer au test de Lee et Strazicich qui est plus performant et qui sera celui dont on tirera la conclusion sur le PGD de notre série diférencié à l’odre 1.

Le modèle final que l’on va conserver est le modèle break avec 1 rupture (car on avait δ1 et δ2 significatifs dans ZA).
Je n’ai pas réalisé le test LS à deux dates de ruptures car cela affecte la puissance de calcul.

source("C:/Users/Bekim/Desktop/R_USA/LeeStrazicichUnitRootTest.R", encoding = 'UTF-8')
source("C:/Users/Bekim/Desktop/R_USA/LeeStrazicichUnitRootTestParallelization.R", encoding = 'UTF-8')

library(parallel)
library(doSNOW)
myBreaks = 1
myModel = "break"
myLags <- 4
cl <- makeCluster(max(1, detectCores() - 1))
registerDoSNOW(cl)
myParallel_LS <- ur.ls.bootstrap(y=ydiff , model = myModel, breaks = myBreaks, lags = myLags, method = "Fixed",pn = 0.1, critval = "bootstrap", print.results = "print")

## [[1]]
## [1] -3.556002
## 
## [1] "First possible structural break at position: 40"
## [1] "The location of the first break - lambda_1: 0.8 , with the number of total observations: 49"
## Critical values - Break model:
##      lambda    1%    5%   10%
## [1,]    0.1 -5.11 -4.50 -4.21
## [2,]    0.2 -5.07 -4.47 -4.20
## [3,]    0.3 -5.15 -4.45 -4.18
## [4,]    0.4 -5.05 -4.50 -4.18
## [5,]    0.5 -5.11 -4.51 -4.17
## [1] "Number of lags used: 4"
## Runtime:
## Time difference of 0.003393765 mins

Le $𝜆$ estimé est 0.8. On compare donc la statistique =-3.556002 à -4,51. Comme -3.556002 > -4,51, on accepte $𝐻_0$ et il y a une présence de racine unitaire. Donc la conclusion du test LS est que le PGD qui a généré notre série est DS avec une date de rupture dans la constante et dans la tendance en 2012.

Pour conclure cette partie sur les tests de racine unitaire, nous retenons le fait que le PGD qui a généré la série différenciée à l’ordre 1 est DS avec une date de rupture.

TESTS-Série diffénrenciée à l’odre 2

Test de DF

• Trend

summary(ur.df(yd2,type="trend",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression trend 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.6024 -0.4432 -0.1089  0.3976  1.9794 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.042035   0.217321   0.193    0.848    
## z.lag.1     -0.971023   0.149518  -6.494 6.32e-08 ***
## tt          -0.001487   0.007884  -0.189    0.851    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7327 on 44 degrees of freedom
## Multiple R-squared:  0.4905, Adjusted R-squared:  0.4674 
## F-statistic: 21.18 on 2 and 44 DF,  p-value: 3.604e-07
## 
## 
## Value of test-statistic is: -6.4943 14.1221 21.1813 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -4.15 -3.50 -3.18
## phi2  7.02  5.13  4.31
## phi3  9.31  6.73  5.61

La p-value de 𝛽1 = 0.851 est supérieur à 𝛼 = 0,05 donc on accepte $H_0$ et $𝛽_1$ n’est pas significatif.(on passe à drift)

• Drift

summary(ur.df(yd2,type="drift",lag=0))

## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression drift 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.57624 -0.44543 -0.07825  0.38151  1.96014 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.006357   0.105722   0.060    0.952    
## z.lag.1     -0.972076   0.147804  -6.577 4.33e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7248 on 45 degrees of freedom
## Multiple R-squared:  0.4901, Adjusted R-squared:  0.4788 
## F-statistic: 43.25 on 1 and 45 DF,  p-value: 4.334e-08
## 
## 
## Value of test-statistic is: -6.5768 21.629 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.58 -2.93 -2.60
## phi1  7.06  4.86  3.94

La p-value de $𝛽_0$ = 0.952 est supérieur à 𝛼 = 0,05 donc on accepte $H_0$ et $𝛽_0$ n’est pas significatif.(on passe à none)

Initialisation de la série

Convertir la série en taux grâce au pib (voir code)

#Prendre les valeurs de 1972 a 2021 
PIBusa= slice(PIBusa, -c(1:2))

#diviser valeur d'inv par pib * 100 pour l'avoir en taux
Y=ts((Iusa$Value/PIBusa$Value)*100,freq=1,start=1972,end=2021)

Visualisation de la série

Y # vérifie les valeurs

## Time Series:
## Start = 1972 
## End = 2021 
## Frequency = 1 
##  [1]  18.10227  18.79561  17.97266  16.71187  16.82914  17.79188  19.53423
##  [8]  20.01125  20.10827  20.60408  18.64693  18.77015  20.69743  21.16697
## [15]  20.63735  19.84664  19.73593  20.17562  19.46639  18.06257  17.74055
## [22]  18.03587  18.41757  18.58343  18.99204  19.22174  19.67576  20.74027
## [29]  21.71004  20.81805  19.26890  18.99797  19.56768  20.66081  21.54968
## [36]  22.31710  20.84163  17.45763  16.95622  17.35371  18.07093  18.24608
## [43]  19.05729  19.11441  18.54922  18.80977  20.39429  19.23850 290.52091
## [50] 293.16679

Sources

Liste crises monétaires utilisé pour le chronogramme.
Bureau of Economic Analysisl pour comprendre le taux d’investissement.
Cours Econométrie des séries temporelles M1

Taux d’investissement des USA

Albin Trshani

2023-04-16

Introduction

I - Visualisation des données

Résumé des données:

• Chronogramme

• (ACF) et (PACF)

II- Tests de racine unitaire

• Test de Dickey-Fuller (DF)

• None

• Test de Dickey-Fuller augmenté (ADF)

• Test de Zivot et Andrews (ZA)

• Intercept ou “Crash”

• Test de Lee et Strazicich (LS)

III-PGD de la série différenciée(2)

• Test de Dickey-Fuller

• None

• Test de ADF

• Test de ZA

• Both

• Test de LS

IV - Autocorrélation

V - Estimation du modèle ARMA(p,q)

Recherche du modèle

VI - Etude des aléas

• Test de Jarque-Bera

• Test de Student (t.test)

• Test de Ljung-Box sur les résidus centrés-réduits

• Test d’Engle sur les résidus centrés - réduits

VII - Prévision

Annexe

Tests de DF non concluants

• Trend

• Drift

Test de ZA non concluant

• Both

LS non concluant

• Break

TESTS-Série diffénrenciée à l’odre 1

• Test de DF

• Trend

• Drift

• None

• Test de ADF

•Test de ZA

• Test de LS

TESTS-Série diffénrenciée à l’odre 2

Test de DF

Initialisation de la série

Sources