From da843333489db3727c46220c498ea6fccabe5c13 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?St=C3=A9phane=20Adjemian=20=28Charybdis=29?= Date: Wed, 23 Feb 2022 15:20:18 +0100 Subject: [PATCH] Various fixes in blog post. --- Makefile | 2 + blog/simulation-du-modele-de-solow/index.org | 69 ++++++++++++-------- 2 files changed, 44 insertions(+), 27 deletions(-) diff --git a/Makefile b/Makefile index 4dda93b..a75916a 100644 --- a/Makefile +++ b/Makefile @@ -11,6 +11,8 @@ publish: sed -i '/‎<\/title>/d' ./output/teaching/economic-calculus/index.html sed -i '/<title>‎<\/title>/d' ./output/teaching/growth/index.html sed -i '/<title>‎<\/title>/d' ./output/teaching/time-series/index.html + sed -i 's/<pre class="src src-python">/<pre><code class="language-python">/g' ./output/posts/simulation-du-modele-de-solow/index.html + sed -i 's/<\/pre>/<\/code><\/pre>/g' ./output/posts/simulation-du-modele-de-solow/index.html assets: @rsync --recursive -avz assets/fonts output diff --git a/blog/simulation-du-modele-de-solow/index.org b/blog/simulation-du-modele-de-solow/index.org index b9a4ea8..ef33cbe 100644 --- a/blog/simulation-du-modele-de-solow/index.org +++ b/blog/simulation-du-modele-de-solow/index.org @@ -1,6 +1,9 @@ #+OPTIONS: H:3 num:nil toc:nil \n:nil @:t ::t |:t ^:nil -:t f:t *:t TeX:t LaTeX:t skip:t d:t tags:not-in-toc creator:t timestamp:nil author:nil title:nil #+HTML_HEAD: <link rel="stylesheet" type="text/css" href="../../css/stylesheet-blog.css"> #+HTML_HEAD: <link rel="stylesheet" type="text/css" href="../../fontawesome/css/all.css" /> +#+HTML_HEAD: <link rel="stylesheet" href="../../highlight/styles/base16/solarized-light.min.css"> +#+HTML_HEAD: <script src="../../highlight/highlight.min.js"></script> +#+HTML_HEAD: <script>hljs.highlightAll();</script> #+LANGUAGE: fr #+TITLE: Simulation du modèle de Solow #+DATE: Février 2022 @@ -24,19 +27,6 @@ physique agrégé, $K$, est donc : \dot K(t) = s K(t)^{\alpha}\bigl(A(t)L(t)\bigr)^{1-\alpha} - \delta K(t) \] -On aimerait, pour une condition initiale arbitraire, déterminer la trajectoire -du capital ou de la production à partir de cette équation. Il s'agit d'une -équation différentielle non linéaire que nous pourrions résoudre, au moins -numériquement, mais en l'état l'analyse de la dynamique serait malaisée car -cette équation n'est pas autonome : le lien entre la variation de $K$ et son -niveau n'est pas stable dans le temps, à cause de la population et de -l'efficacité du travail qui augmentent à chaque instant. Au préalable on -transforme le modèle en éliminant ces deux sources exogènes de croissance. Si on note $n$ -et $x$ les taux de croissance (supposés constants) de la population $L$ et de -l'efficacité du travail $A$, on peut facilement montrer que la dynamique du -stock de capital par tête efficace $\hat k(t) = \frac{K(t)}{A(t)L(t)}$ est -donnée par : - où $s\in]0,1[$ est le taux d'épargne, $\delta\in]0,1[$ le taux de dépréciation, $L$ la population, $A$ l'efficacité du travail, et on notera $Y=K^{\alpha}\bigl(AL\bigr)^{1-\alpha}$ la production (identique au revenu dans @@ -44,6 +34,21 @@ ici). Cette équation nous dit que le stock de capital physique agrégé de l'économie augmente si et seulement si l'investissement est supérieur à la dépréciation du capital. +\\ + +On aimerait, pour une condition initiale arbitraire, déterminer la trajectoire +du capital ou de la production à partir de cette équation. Il s'agit d'une +équation différentielle non linéaire que nous pourrions résoudre, au moins +numériquement, mais en l'état l'analyse de la dynamique serait malaisée car +cette équation n'est pas autonome : le lien entre la variation de $K$ et son +niveau n'est pas stable dans le temps, car la population et +l'efficacité du travail augmentent à chaque instant. Au préalable on +transforme le modèle en éliminant ces deux sources exogènes de croissance. Si on note $n$ +et $x$ les taux de croissance (supposés constants) de la population $L$ et de +l'efficacité du travail $A$, on peut facilement montrer que la dynamique du +stock de capital par tête efficace $\hat k(t) = \frac{K(t)}{A(t)L(t)}$ est +donnée par : + \[ \dot {\hat k}(t) = s \hat k(t)^{\alpha} - (n+x+\delta)\hat k(t) \] @@ -55,6 +60,8 @@ le niveau invariante dans le temps, il sera plus simple d'étudier la dynamique partir de cette équation. On pourra toujours revenir aux variables agrégées, en notant que par définition on a $K(t)=\hat k(t)A(t)L(t)$. +\\ + Cette équation différentielle admet un unique état stationnaire strictement positif : @@ -77,12 +84,14 @@ analytiquement l'équation différentielle, mais nous pourrions établir de faç plus générale ce résultat (cela fera l'objet d'un autre article) sans identifier explicitement la solution de l'équation différentielle. +\\ + Nous utilisons Python, ainsi que les librairies =numpy= et =scipy= pour résoudre numériquement l'équation différentielle. Nous devons d'abord définir le problème à résoudre. La fonction suivante retourne la variation du stock de capital par tête efficace associée à un niveau de capital par tête efficace : -#+begin_src python :session :exports code +#+begin_src python :session :exports code def dotk(k, t, alpha, s, n, x, delta): return s*k**alpha-(n+x+delta)*k #+end_src @@ -106,6 +115,8 @@ initiale en ayant les idées précises sur la distance initiale à l'état stationnaire (qui en principe doit aussi être le niveau de long terme du stock de capital par tête efficace). +\\ + Finalement le bloc de code suivant utilise la fonction =odeint= de la librairie =scipy= pour résoudre l'équation différentielle : @@ -120,7 +131,7 @@ Finalement le bloc de code suivant utilise la fonction =odeint= de la librairie n = .02 x = .02 delta = .02 - # Condition initiale (50% de l'état stationnaire) + # Conditions initiales (en déviation à l'état stationnaire) kinit0 = 0.05*kstar(alpha, s, n, x, delta) kinit1 = 1.95*kstar(alpha, s, n, x, delta) # Résolution de l'équation différentielle (transition vers l'état stationnaire) @@ -186,9 +197,11 @@ $\hat k<\hat k^{\star}$ et que le taux de croissance est une fonction monotone décroissante du niveau de stock de capital par tête efficace (ces deux propriétés sont des conséquences de l'hypothèse $\alpha<1$). +\\ + La figure suivante représente la fonction $g(\hat k)$ (la courbe bleue), elle -passe bien par 0 à l'état stationnaire (intersection des droites rouges) et -conformément à ce que nous pouvions voir plus haut elle est monotone +passe bien par 0 à l'état stationnaire (intersection des droites rouges) et, +conformément à ce que nous pouvions voir plus haut, elle est monotone décroissante (on devine l'asymptote verticale en 0 et l'asymptote horizonale en $-(n+x+\delta)$). #+begin_src python :session :exports none @@ -213,6 +226,8 @@ n'est pas forcément très intéressant de simuler ce modèle, mais l'approche numérique peut se révéler intéressante, voire indispensable, si on complexifie le modèle. +\\ + Il convient néanmoins de s'interroger sur la pertinence de cette solution numérique. D'un point de vue qualitatif nous venons de voir qu'il n'y avait a priori pas de problème ici (nous retrouvons les propriétés attendues). Dans la @@ -267,7 +282,7 @@ second membre ($Ae^{-(1-\alpha)(n+x+\delta)t}$) : z(t) = \frac{s}{n+x+\delta} + A e^{-(1-\alpha)(n+x+\delta)t} \] -Si on dispose d'une condition initiale z(0), alors la solution qui satisfait cette condition est : +Si on dispose d'une condition initiale $z(0)$, alors la solution qui satisfait cette condition est : \[ z(t) = \frac{s}{n+x+\delta} + \left(z(0)-\frac{s}{n+x+\delta}\right) e^{-(1-\alpha)(n+x+\delta)t} @@ -293,6 +308,8 @@ Nous avons donc bien que, dans ce modèle, l'état stationnaire est asymptotiquement globalement stable. En passant, on reconnaît sous l'exponentielle la vitesse de convergence $\beta = (1-\alpha)(n+x+\delta)$. +\\ + Nous souhaitons maintenant comparer cette solution, par construction exacte, avec la solution numérique. La fonction suivante permet de calculer $\hat k$ a un instant $t$ quelconque étant donnée une condition initiale $\hat k(0)$. @@ -392,13 +409,13 @@ courbe bleue par sa tangente en $\hat k^{\star}$ (la droite verte dans la figure #+RESULTS: [[file:None]] -#+CAPTION: *Taux de croissance du stock de capital physique par tête efficace (approximation).* +#+CAPTION: *Taux de croissance approximé du stock de capital physique par tête efficace.* #+LABEL: fig:k-growth-2 [[file:img/k-growth-2.svg]] On devine sur ce graphique que les erreurs d'approximation (la différence entre la courbe bleue et la droite verte) sont ici beaucoup plus importantes que les -erreurs numériques calculées dans la section précédante. On remarque que les erreurs +erreurs numériques calculées dans la section précédante. On remarque que les erreurs : - sont d'autant plus importantes que l'on s'éloigne de l'état stationnaire, @@ -410,12 +427,14 @@ erreurs numériques calculées dans la section précédante. On remarque que les convexe, le taux de croissance approximé est toujours inférieur au taux de croissance théorique. +\\ + Ces différences, que nous observons ici graphiquement, ont évidemment des conséquences sur les trajectoires que nous pouvons calculer. Avec l'approximation d'ordre 1, nous avons : \[ -g(\hat k) \approx -(1-\alpha)s \hat k^{\star}^{\alpha-1}\frac{\hat k - \hat k^{\star}}{\hat k^{\star}} +g(\hat k) \approx -(1-\alpha)s \left. \hat k^{\star}\right. ^{\alpha-1}\frac{\hat k - \hat k^{\star}}{\hat k^{\star}} \] En substituant lé définition de l'état stationnaire, il vient : @@ -469,8 +488,8 @@ initiale à l'état stationnaire est donnée par : La dernière équation, issue de la log-linéarisation de la dynamique, doit être comparée à l'équation exacte obtenue plus haut pour $\hat k(t)$. La figure -[[fig:k-approx-1]] compare les dynamiques approximées (tirets) et exactes (traits -pleins). Comme nous pouvions l'anticiper la différence est plus importante quand +[[fig:k-approx-1]] compare les dynamiques approximées (tirets) et exactes (trait +plein). Comme nous pouvions l'anticiper la différence est plus importante quand on considère une condition initiale inférieure à l'état stationnaire, car c'est quand le stock de capital se rapproche de zéro que la courbure du modèle devient plus importante. @@ -514,7 +533,3 @@ semble aussi significativement plus lente à rejoindre l'état stationnaire. Cec suggère que la mesure de la vitesse d'ajustement que nous avons l'habitude d'utiliser ($\beta = (1-\alpha)(n+x+\delta)$, obtenue en log-linéarisant la dynamique), sous estime la vitesse d'ajustement vers l'état stationnaire. - -# Local Variables: -# org-confirm-babel-evaluate: nil -# End: