R Notebook : Régression Linéaire
Introduction
Les modèles de régression linéaire par les Moindres Carrés Ordinaires (MCO) sont des modèles d’estimation de relations entre une variable endogène \(Y_t\) (ou variable expliquée), et une ou plusieurs variables dites exogènes \(X_{n,t}\) (ou variables explicatives). C’est l’un des (voire, le) modèles les plus simples en économétrie, mais qui peut être très efficace s’il est bien mené et que les données s’y prêtent.
Notamment, une des hypothèses qui est posée par le modèle de regression linéaire par les MCO est l’exogénéité des variables explicatives. Autrement dit, on teste l’impact de \(X_t\) sur \(Y_t\) sous l’hypothèse que \(Y_t\) n’influence pas \(X_t\). Cela implique l’absence de causalité inverse et/ou de simultanéité entre \(X_t\) et \(Y_t\). En général, si on pense qu’il y a cette causalité inverse ou ce biais de simultanéité, on préfère opter pour d’autres types de modèles (modèles VAR notamment) même si des spécifications alternatives de modèle peuvent parfois régler le problème (modèles ARDL).
Une autre hypothèse est qu’il n’y a pas de biais de variable omise : le modèle explique suffisamment \(Y_t\) pour que les résidus, qui représente la dynamique de \(Y_t\) si elle n’était impactée par rien d’autre que par ses propres chocs, soient non-autocorrélés et suivent une loi Normale centrée réduite de variance constante. Autrement dit, la dynamique des résidus représentent la dynamique de \(Y_t\) a son état stationnaire.
Enfin, le biais de variable omise pourrait inciter à inclure la totalité des variables disponibles dans le monde. Cependant, cela peut créer un biais important : le bad control. Par exemple, si on veut estimer l’impact de \(X_{1,t}\) et de \(X_{2,t}\) sur \(Y_t\) mais qu’en même temps, \(X_1\) impacte \(X_2\), alors l’estimation du coefficient associé à \(X_2\) va être biaisé par ce lien. Il faut donc respecter le principe de parsimonie : en faire suffisamment, mais ne pas en faire trop.
Modèle de Régression Linéaire Univarié
\[\begin{equation} Y_t = a_0 + \beta X_t + \varepsilon_t \end{equation}\]Modèle de Régression Linéaire Multivarié
\[\begin{equation} Y_t = a_0 + \beta_1 X_{1,t}+ \beta_2 X_{2,t} + \beta_3 X_{3,t} + \varepsilon_t \end{equation}\]Modèle de Régression Non-Linéaire Univarié
\[\begin{equation} Y_t = a_0 + \beta_1 X_{1,t}+ \beta_2 X_{2,t} + \beta_3 X_{2,t}^2 + \varepsilon_t \end{equation}\]Courbe de Phillips
On va étudier ici la validité de la courbe de Phillips pour la France. La courbe de Phillipps illustre la relation inverse entre taux de chômage et taux d’inflation démontrée empiriquement par l’économiste néo-zélandais William Phillips en 1958. Le concept a ensuite été repris par Samuelson et Solow dans le développement des modèles néo-keynésiens.
Afin de justifier cette relation sur le plan théorique, l’explication la plus populaire se fonde sur les pouvoirs de négociation des salariés. Lorsque le taux de chômage est haut, les salariés sont “facilement” remplaçables dans les entreprises. Ils ne peuvent donc pas négocier de hausse de salaires car les entreprises les remplaçeraient par d’autres au chômage. Si les salaires n’augmentent pas, les prix pratiqués par les entreprises n’ont pas de raison d’augmenter (parce que les coûts n’augmentent pas d’une part, et parce que la demande n’augmente pas d’autre part), donc l’inflation n’augmente pas. A contrario, lorsque le taux de chômage est faible, les salariés ont un pouvoir de négociation important. Ils vont donc demander des hausses de salaire, ce qui va se répercuter sur les prix.
Il y aurait donc un arbitrage à effectuer. En effet, les Banques Centrales ont pour mission de stabiliser l’inflation à un niveau faible. La validité de la courbe de Phillips implique donc que cette inflation faible entraîne mécaniquement un chômage structurel permanent et l’impossibilité du Plein-Emploi. Mais dans le même temps, les Etats (et les pouvoirs élus) ont un intérêt d’atteindre le Plein-Emploi (pour être réélus), donc à mener des politiques monétaire accomodantes qui entraînent une hausse de l’inflation sur le long-terme. Les sociétés doivent donc choisir, entre l’inflation et le chômage.
\[\begin{equation} \pi_t = a_0 + \beta U_t + \varepsilon_t \end{equation}\]où \(U_t\) est le taux de chômage (en %) de la France entre 1960 et 2020, \(\pi_t\) est l’inflation sur la même période en France, mesurée comme la différence de logarithme, et \(\varepsilon_t\) représentent les résidus.
Données
donnees <- read_excel("C:/Users/fkraus/Desktop/data_schularick.xlsx") %>%
filter(country == "France") %>%
select(year, unemp, cpi) %>%
mutate(inflation = (log(cpi) - lag(log(cpi)))*100)%>%
na.omit()%>%
filter(year >= 1950)
donneesggplot(donnees, aes(x=year, y=inflation))+
geom_line(lwd=1)+
geom_line(aes(y=unemp), col="red", lwd=1)+
theme_bw()
On voit que dans l’ensemble sur la période, inflation et chômage tendent à avoir une dynamique opposée : l’inflation est haute lorsque le chômage est bas, et inversement.
Relation Inflation-Chômage
Le meilleur moyen de visualiser la relation entre deux variables (en tout cas, pour les modèles univariés) est de tracerle nuage de point entre les deux variables :
ggplot(donnees, aes(x=unemp, y=inflation))+
geom_point()+
geom_smooth(method="lm")+
theme_bw()
Il semble y avoir une relation négative et significative entre le taux d’inflation et le taux de chômage. On peut le vérifier avec la régression directement :
regression <- feols(inflation ~ unemp, data=donnees)
etable(regression) regression
Dependent Var.: inflation
Constant 7.104*** (0.8730)
unemp -0.4323*** (0.1157)
_______________ ___________________
S.E. type IID
Observations 71
R2 0.16818
Adj. R2 0.15613
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1On voit effectivement l’impact négatif et significatif de l’inflation sur le taux de chômage. D’après les coefficients, une augmentation du taux de chômage d’une unité (un point de pourcentage en l’occurence), le taux d’inflation tend à diminuer de 0.43 unités (donc 0.43 points de pourcentages). Autrement dit, nos résultats montrent bien qu’il y a un arbitrage entre inflation et chômage.
En revanche, il se peut que cette relation soit biaisée : on n’a pas vérifié la stationnarité des variables inclues dans le modèle !
Vérification de la stationnarité
Inflation
summary(ur.df(donnees$inflation))
###############################################
# 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
-11.4646 -0.1584 0.3398 1.2137 9.6494
Coefficients:
Estimate Std. Error t value Pr(>|t|)
z.lag.1 -0.15403 0.05198 -2.964 0.00421 **
z.diff.lag 0.04228 0.11262 0.375 0.70854
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.448 on 67 degrees of freedom
Multiple R-squared: 0.1163, Adjusted R-squared: 0.08993
F-statistic: 4.409 on 2 and 67 DF, p-value: 0.01589
Value of test-statistic is: -2.9635
Critical values for test statistics:
1pct 5pct 10pct
tau1 -2.6 -1.95 -1.61La série d’inflation en niveau n’est pas stationnaire (la t-stat n’est pas inférieure à la valeur critique à 5%. Il faut donc explorer le type de non-stationnarité
summary(ur.df(donnees$inflation, type="trend"))
###############################################
# Augmented Dickey-Fuller Test Unit Root Test #
###############################################
Test regression trend
Call:
lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
Residuals:
Min 1Q Median 3Q Max
-10.7968 -0.8037 -0.1813 0.8070 8.3911
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.24104 0.90244 2.483 0.015606 *
z.lag.1 -0.36310 0.08906 -4.077 0.000127 ***
tt -0.02489 0.01651 -1.508 0.136378
z.diff.lag 0.15430 0.11484 1.344 0.183745
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.343 on 65 degrees of freedom
Multiple R-squared: 0.2087, Adjusted R-squared: 0.1722
F-statistic: 5.714 on 3 and 65 DF, p-value: 0.001557
Value of test-statistic is: -4.0771 5.897 8.5434
Critical values for test statistics:
1pct 5pct 10pct
tau3 -4.04 -3.45 -3.15
phi2 6.50 4.88 4.16
phi3 8.73 6.49 5.47La trend est significative, et le coefficient d’autocorrélation est également significatif. L’inflation suit donc un processus DS, qu’il faut passer en différence pour stationnariser
donnees <- donnees %>%
mutate(dinflation = inflation - lag(inflation))%>%
na.omit()
summary(ur.df(donnees$dinflation))
###############################################
# 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
-11.8436 -0.6425 -0.1740 0.5989 9.4231
Coefficients:
Estimate Std. Error t value Pr(>|t|)
z.lag.1 -1.2165 0.1642 -7.410 3.01e-10 ***
z.diff.lag 0.2411 0.1126 2.141 0.036 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.49 on 66 degrees of freedom
Multiple R-squared: 0.5265, Adjusted R-squared: 0.5121
F-statistic: 36.69 on 2 and 66 DF, p-value: 1.931e-11
Value of test-statistic is: -7.4096
Critical values for test statistics:
1pct 5pct 10pct
tau1 -2.6 -1.95 -1.61La série d’inflation en différence est stationnaire en différence première, elle est intégrée d’ordre 1: \(I(1)\)
Chômage
summary(ur.df(donnees$unemp, type="none"))
###############################################
# 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.73911 -0.31311 0.03339 0.32811 1.85398
Coefficients:
Estimate Std. Error t value Pr(>|t|)
z.lag.1 0.0003894 0.0098586 0.039 0.9686
z.diff.lag 0.3098390 0.1185854 2.613 0.0111 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.6146 on 66 degrees of freedom
Multiple R-squared: 0.09587, Adjusted R-squared: 0.06848
F-statistic: 3.499 on 2 and 66 DF, p-value: 0.03594
Value of test-statistic is: 0.0395
Critical values for test statistics:
1pct 5pct 10pct
tau1 -2.6 -1.95 -1.61La série de chômage n’est pas stationnaire en niveau. On regarde donc le type de non-stationnarité
summary(ur.df(donnees$unemp, type="trend"))
###############################################
# Augmented Dickey-Fuller Test Unit Root Test #
###############################################
Test regression trend
Call:
lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
Residuals:
Min 1Q Median 3Q Max
-1.62799 -0.38689 -0.09593 0.28019 1.75622
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.214037 0.159619 1.341 0.1847
z.lag.1 -0.037551 0.035179 -1.067 0.2898
tt 0.002646 0.006719 0.394 0.6950
z.diff.lag 0.312849 0.121616 2.572 0.0124 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.6123 on 64 degrees of freedom
Multiple R-squared: 0.1134, Adjusted R-squared: 0.07182
F-statistic: 2.728 on 3 and 64 DF, p-value: 0.05122
Value of test-statistic is: -1.0674 0.8323 0.9375
Critical values for test statistics:
1pct 5pct 10pct
tau3 -4.04 -3.45 -3.15
phi2 6.50 4.88 4.16
phi3 8.73 6.49 5.47la trend n’est pas significative. On teste donc avec le drift
summary(ur.df(donnees$unemp, type="drift"))
###############################################
# Augmented Dickey-Fuller Test Unit Root Test #
###############################################
Test regression drift
Call:
lm(formula = z.diff ~ z.lag.1 + 1 + z.diff.lag)
Residuals:
Min 1Q Median 3Q Max
-1.65504 -0.37445 -0.08824 0.29124 1.80355
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.23291 0.15126 1.540 0.1285
z.lag.1 -0.02613 0.01980 -1.320 0.1914
z.diff.lag 0.30168 0.11749 2.568 0.0125 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.6083 on 65 degrees of freedom
Multiple R-squared: 0.1112, Adjusted R-squared: 0.08389
F-statistic: 4.067 on 2 and 65 DF, p-value: 0.02166
Value of test-statistic is: -1.3201 1.1863
Critical values for test statistics:
1pct 5pct 10pct
tau2 -3.51 -2.89 -2.58
phi1 6.70 4.71 3.86La constante est significative, et la série n’est pas stationnaire. On en conclut donc que la série de chômage est un DS, qu’il faut passer en différence première pour stationnariser.
donnees <- donnees %>%
mutate(dunemp = unemp - lag(unemp))%>%
na.omit()
summary(ur.df(donnees$dunemp))
###############################################
# 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.70353 -0.31249 0.02489 0.32736 1.89772
Coefficients:
Estimate Std. Error t value Pr(>|t|)
z.lag.1 -0.64331 0.14524 -4.429 3.71e-05 ***
z.diff.lag -0.06672 0.12349 -0.540 0.591
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.6147 on 65 degrees of freedom
Multiple R-squared: 0.3485, Adjusted R-squared: 0.3285
F-statistic: 17.39 on 2 and 65 DF, p-value: 8.943e-07
Value of test-statistic is: -4.4292
Critical values for test statistics:
1pct 5pct 10pct
tau1 -2.6 -1.95 -1.61la série de chômage en différence est stationnaire en différence première, donc intégrée d’ordre 1 : \(I(1)\).
Régression linéaire univariée
Maintenant, on sait que les deux séries sont stationnaires lorsqu’elles sont passées en différence première. On peut donc mener la régression avec les séries stationnaires :
regression_stationnaire <- feols( dinflation ~ dunemp, data=donnees)
etable(regression_stationnaire) regression_sta..
Dependent Var.: dinflation
Constant -0.1538 (0.3101)
dunemp -0.7342 (0.4905)
_______________ ________________
S.E. type IID
Observations 69
R2 0.03236
Adj. R2 0.01791
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1On voit ici que le coefficient associé au chômage n’est pas significatif à 5%. Autrement dit, la courbe de Phillips n’est pas valide. Il n’y a donc pas d’arbitrage entre inflation et chômage.
En revanche, un problème apparaît ici : on a fait une régression avec deux variables qui sont intégrées du même ordre : \[\begin{align*} & dinflation \sim I(1) \\ & dunemp \sim I(1) \end{align*}\] Et la co-intégration (intégration du même ordre) peut biaiser le coefficient, il convient alors d’utiliser un modèle qui corrige cette co-intégration : le Modèle à Correction d’Erreurs (voir séance 7 sur la Cointégration et les Modèles à Correction d’Erreurs).