ソローモデルの各パラメーターをいじった時にどうなるか的なやつ

労働節約的なコブダグラス型生産関数：F(K,AL)=(K^a) * (AL)^(1-a)
資本蓄積式は⊿K=sY-δK
政府や海外は存在しない
プログラミングでは微分を扱えない（ことはないだろうけど簡単のため）差分を扱うことに注意すると

# -*- coding: utf-8 -*-
    
import numpy as np
import matplotlib.pyplot as plt

time = 100 #何期までやるか

Y = np.zeros(time) #空き配列の作成
K = np.zeros(time)
L = np.zeros(time)
A = np.zeros(time)
y = np.zeros(time)
k = np.zeros(time)

K[0] = 10 #初期の水準
A[0] = 5
L[0] = 10

a = 0.5 #コブダグラス型生産関数の資本分配率
d = 0.2 #資本減耗率。減価償却率と読み替えても良いかもしれない
s = 0.2 #貯蓄率
g = 0.02 #技術進歩率
n = 0.02 #人口成長率

for t in range(time-1):
    Y[t] = (K[t]**a)*((A[t]*L[t])**(1-a)) #コブダグラス型生産関数
    K[t+1] = s*Y[t] + (1-d)*K[t] #資本蓄積式
    A[t+1] = (1+g)*A[t] #技術進歩
    L[t+1] = (1+n)*L[t] #人口成長
    y[t] = Y[t] / (A[t]*L[t]) #効率労働あたりGDP
    k[t] = K[t] / (A[t]*L[t]) #効率労働あたり資本

plt.subplot(2,2,1) #(縦分割数、横分割数、ポジション)
plt.plot(y,"red")
plt.plot(k,"green")

plt.subplot(2,2,2)
plt.plot(Y,"red")
plt.plot(K,"green")

plt.subplot(2,2,3)
plt.plot(A,"red")
plt.plot(L,"blue")

左上のグラフは効率労働あたりのGDPと資本、右上のグラフは経済全体のGDPと資本
人口成長と技術進歩の天井がない $a_{t+1} = (1+n)a_t$
ので当然の帰結として経済全体も際限なく成長していく
ただ効率労働あたりの数字は早い段階(25期目あたり）で収束しているのが分かる
このパラメータをいじってみる

定常状態では⊿k=0のため、k*とy*は（生産関数と資本蓄積式を連立して）
$\displaystyle k^*={s/(n+g+δ)}^{1/(1-a)}$
$\displaystyle y^*={s/(n+g+δ)}^{a/(1-a)}$
と解ける
見て分かる通り、sの上昇、またはnとgとδの減少はkとyの上昇につながる
逆も然り
これを見てみる

# -*- coding: utf-8 -*-
    
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as ani

time = 100 #何期までやるか

Y = np.zeros(time) #空き配列の作成
K = np.zeros(time)
L = np.zeros(time)
A = np.zeros(time)
y = np.zeros(time)
k = np.zeros(time)

K[0] = 10 #初期の水準
A[0] = 5
L[0] = 10

a = 0.5 #コブダグラス型生産関数の資本分配率
d = 0.1 #資本減耗率。（財が1つしかないので）減価償却率と読み替えても良い
s = 0.1 #貯蓄率
g = 0.01 #技術進歩率
n = 0.01 #人口成長率

def solow(a,d,s,g,n):
    for t in range(time-1):
        Y[t] = (K[t]**a)*((A[t]*L[t])**(1-a)) #コブダグラス型生産関数
        K[t+1] = s*Y[t] + (1-d)*K[t] #資本蓄積式
        A[t+1] = (1+g)*A[t] #技術進歩
        L[t+1] = (1+n)*L[t] #人口成長
        y[t] = Y[t] / (A[t]*L[t]) #効率労働あたりGDP
        k[t] = K[t] / (A[t]*L[t]) #効率労働あたり資本

fig = plt.figure()
imgset = []

for i in range(10):
    new = s + i/10
    solow(a,d,new,g,n)
    imgk, = plt.plot(k,color="red") #変数のあとに,つけること
    imgy, = plt.plot(y,color="blue")
    imgset.append([imgk,imgy])
    
animation = ani.ArtistAnimation(fig,imgset)
plt.show()
animation.save("solow_s_yk.mp4", writer="ffmpeg")
#animation.save("solow_s_yk.gif", writer="imagemagick")

gifで保存するなら要imagemagick、mp4で保存するなら要ffmpeg
vimeo.com
↑こんなかんじ、sの上昇とともに効率労働あたりのGDPと資本が底上げされてるのがわかる
vimeo.com
↑こちらはsの上昇による経済全体への影響
以下同様に他のパラメータも動かしてみる

for i in range(10):
    new = d + i/50 #資本減耗率0.1から0.02刻み
    solow(a,new,s,g,n)
    imgk, = plt.plot(k,color="red")
    imgy, = plt.plot(y,color="blue")
    imgset.append([imgk,imgy])
animation = ani.ArtistAnimation(fig,imgset)
plt.show()
animation.save("solow_d_yk.mp4", writer="ffmpeg")
#animation.save("solow_s_yk.gif", writer="imagemagick")

vimeo.com
↑資本減耗率を0.1から0.02刻みで増やしていった場合のyとk
減耗率が高いと資本が減耗（減少）していく様が見て取れる
vimeo.com
↑経済全体でも同様

for i in range(30):
    new = n + i/1000
    solow(a,d,s,g,new)
    imgk, = plt.plot(Y,color="red")
    imgy, = plt.plot(K,color="blue")
    imgset.append([imgk,imgy])
animation = ani.ArtistAnimation(fig,imgset)
plt.show()
animation.save("solow_n_yk.mp4", writer="ffmpeg")

↑人口成長率を0.01から0.001刻みで0.04まで（人口が成長していくのではなく人口成長率そのものが変化していることに注意）
vimeo.com
↑経済全体では急激に（指数的に）成長しているのがわかる、が
↓効率労働単位では逆に減少している（効率労働あたりのGDP資本の成長率＜人口成長率　なので）
もちろん人口成長率は技術進歩率に読み替えても同じ
vimeo.com

ところでspyderだとmatplotlib.animationがデフォで動かないので以下参照
python - Animation from matplotlib not working in spyder - Stack Overflow

無機野

ここはお前の日記帳

ソローモデルの各パラメーターをいじった時にどうなるか的なやつ