最近になってようやく、Maxima というフリーソフトがあることに気づきました。無料のソフトでここまでできれば十分じゃないでしょうか?(といっても使いこなせていないのでよくわかりません。また、プログラムに関してはまったくの素人なので評価もできません。)
img_JIRO

 ぜいたく言うな、なんといってもタダだ



 Maxima のページ(http://maxima.sourceforge.net)からダウンロードするだけで使えます(このページは英語ですが、日本語のマニュアルがインターネット上にいろいろあります)。

 前回の記事で、ラムゼイ・モデルの場合、r (資本収益率、資本からのリターン)の変化や、g (経済成長率)の変化を取り出して示すことははできない(できない?)、と書きました。そこで、今回は、Maxima を使って、ラムゼイモデルの r と g の変化を示してみようと思います(といっても、やはり r と g の変化を方程式のかたちで明示的に示すことはできません)。また、Maxima を使ってラムゼイモデルにおける均衡にいたる経路(サドルパス)を描いてみます。

 しかし、連続時間モデルのままでは解けません。そこで、最初に、ラムゼイモデルを離散時間モデルに直すことが必要です。

 経済の人口増加率は n,知識・技術の増加率を g とします。t期の知識・技術のレベルをA,労働力人口を Lとすると、t+1 期の知識・技術のレベルは A t+1=(1+g)A,t+1 期の労働力人口は Lt+1=(1+n)Lです。
 生産関数は収穫一定とします(つまり、生産要素を t 倍すると、生産量も t 倍になる)。次の関係が成り立ちます。
150214 (1)   (1)

 個人の最適化は次の目的関数を最大化することです(ρは主観的時間割引率です)。また、このモデルの個人は全員同質であると想定します。
150214 (2)
   (2)

 以下、1人当たりの消費 c と資本ストック k には、後から出てくる効率労働当たりの消費と資本と区別するためにハット( ^ )をつけて表記します。
img_JIRO





  
 以下、1人当たりで計算していきますが、1人当たりで考えるのは、個人の最適化(効用最大化)は、個人が1人当たりで考えて行われる、と想定しているからです(ローマーや Barro and Sala-i-Martin の教科書に出てくるラムゼイモデルは「家計当たり」です)。
 個人の予算制約は次のものです。
150214 (3)   (3)

 wは t 期の時間当たり賃金、hは t 期の労働時間、rは t 期の資本1単位からのリターン、iは t 期の1人当たりの投資です。(3)式は、個人は、労働所得と資本所得を所得とし(左辺)、消費と投資に支出する(右辺)、ということを表しています。この(3)式を制約条件として最適化問題を解こうとすると、変数が多すぎます、そこで。
 経済全体では、次の関係が成り立っています(政府支出と税金はなしです)。
150214 (4)
   (4)
 生産(左辺)が消費と投資(右辺)になる、ということを表しています。この(4)式を労働力人口 Lで割ります。生産関数は収穫一定と想定しているので、
150214 (5)



 が成り立ちます。そこで、Lで(4)式をで割ると、
150214 (6)
   (5)
 となります。この(5)式と(3)式を比べると、
150214 (7)

 とならなければならない、とわかります。そこで、(3)式の左辺を f(k^ ) に変えれば、
150214 (6)
   (6)
 です。少し変形すれば、
150214 (8)
   (7)
 となります。これが制約条件になるのですが、最適化問題を解こうとするとまだ変数が多いので、iを消したい(つまり、c と k だけの式にしたい)、そこで。
 経済全体の資本蓄積過程は、次の式で表されます。
150214 (9)

 Kは t 期の経済全体の資本ストックです。Lt+1で割り、1人当たりに直します。
150214 (10)



 Lt+1=(1+n)Lなので、
150214 (12)


 iに(7)式から iを求め代入すれば、
150214 (13)
   (8)

 となります。この(8)式を最適化問題の制約条件とします。ラグランジュ方程式を設定します。
150214 (14) 



 一階の条件として、
150214 (15)
   (9)


150214 (16)
   (10)


150214 (17)
   (11)


 が得られます。ここで、u(・), f(・) の添え字は、その変数で微分した関数であることを表しています。
(9)式、(10)式、(11)式から、
150214 (18)

   (12)

 が得られます。これがこのモデルのオイラー方程式(t 期の消費と t+1 期の消費の関係を表す式)になります。効用関数、生産関数をそれぞれ次のような対数型効用関数、コブ・ダグラス型生産関数とすると、
150213 (9)
   (13)
 (12)式は、
150214 (19)

   (14)

 となります。これを効率労働当たりに直したいです( 効率労働当たり= AL1単位当たりにする)。なぜなら、技術・知識の増減、人口の増減を想定している場合、効率労働当たりに直したほうが、定常状態を考えるときに簡単だからです。技術・知識の増減、人口の増減を想定している場合、定常状態でも、1人当たりのcやkは変化しています。効率労働当たりでは一定と想定できます。
 効率労働当たりの消費を cとすると、(14)式の左辺は、次のようになります。
150214 (20)



 この関係を使って、(14)式の左辺を効率労働当たりの消費に変えると、
150214 (21)

   (15)

 となります。右辺の分子にある、生産関数の1人当たりの資本ストックによる微分は、効率労働当たりに直しても同じです。なぜなら (効率労働当たりの資本ストックを k とします。k =K/(A) です)、
150214 (22)




150214 (23)



 が成り立つからです。したがって、(15)式の右辺を効率労働当たりに直すと、
150213 (8)



 という関係が得られます。生産関数が f(k)=kα ((13)式)なら、
150213 (10)

   


 となります。cを右辺に移項すれば、
150213 (11)

   (16)

 となります。

 既出の(8)式も効率労働当たりに直したいです。(8)式のk^ ,ck^ を k ,c に変える場合、k^ =A, c^ =Aなので、
150214 (24)


 となります。A t+1=(1+g)Aなので、
150214 (25)


 これを整理すれば、
150213 (2)


 となります。生産関数が f(k)=kα ((13)式)の場合、
150213 (12)
   (17)

 となります。この(16)式と(17)式が、効率労働当たりの資本ストックと効率労働当たりの消費の時間変化を表す式です(これを使ってプログラムで計算します)。適当な初期値を与えれば、その後の変化がわかります。

 定常状態を求めておくと、定常状態では、ct+1 = c=c,kt+1 = kt+1 =k になります( は定常値を表します。先に述べたように、効率労働当たりでないとこの関係が使えません)。 定常状態では、(10)式の左辺が 1 になります。
150213 (13)
   (18)

 この関係から、kの定常値 kが求まります。
150213 (14)
   (19)

 また、(17)式から、定常状態では次の関係が成り立ちます。
150213 (15)


 cの式に直せば、
150213 (16)
   (20)
 となります。これに、kを代入すれば、cが求まります。

(1)α=0.36, g=0.02, n=0.01, ρ=0.02, δ=0.1 の場合
 定常状態における k と c の値は、(19)式、(20)式から、
k* = 3.894507656030513
c* = 1.124344360296009
  になります。その均衡から外れた場所から、どのように均衡にたどり着くのかを見てみます。

[a:0.36,n:0.01,g:0.02,delta:0.1,rho:0.02]$
k1[0]:1$
c1[0]:0.5$ 
k1[i]:=(k1[i-1]^a-c1[i-1]+(1-delta)*k1[i-1])/((1+g)*(1+n))$
c1[i]:=(a*k1[i]^(a-1)+1-delta)*c1[i-1]/((1+g)*(1+n)*(1+rho))$

 パラメーターと初期値を与え、(16)式、(17)式から、k と c の漸化式を上記のようにつくります。k と c の値を組み合わせて、リストをつくります。

R1:makelist([k1[i],c1[i]],i,0,30)$

 (R1はリストの名前。Ramsey のRのつもりですが) このリストの離散値をグラフで表示させれば(plot2d(discrete データ名)  という形で離散値データをグラフで表示させる)、 初期値以後の経路が示されます。また、上記の$記号をセミコロン(;)に変えて、プログラムを実行させれば、それぞれの数値も示してくれます。ただし、上記の初期値では、均衡にたどりつけません。

  適当に値を入れながら、探していくと・・・  最初の効率労働当たりの資本ストック k が 1 のとき、均衡にたどり着くのは、消費の初期値が 0.54815 のときです。そこで、それより低い c=0.5 のときと、それより高い c=0.58 のときも合わせてグラフに描くと次のようになります。

k1[0]:1$
c1[0]:0.5$ 
k1[i]:=(k1[i-1]^a-c1[i-1]+(1-delta)*k1[i-1])/((1+g)*(1+n))$
c1[i]:=(a*k1[i]^(a-1)+1-delta)*c1[i-1]/((1+g)*(1+n)*(1+rho))$

k2[0]:1$
c2[0]:0.54815$ 
k2[i]:=(k2[i-1]^a-c2[i-1]+(1-delta)*k2[i-1])/((1+g)*(1+n))$
c2[i]:=(a*k2[i]^(a-1)+1-delta)*c2[i-1]/((1+g)*(1+n)*(1+rho))$

k3[0]:1$
c3[0]:0.58$ 
k3[i]:=(k3[i-1]^a-c3[i-1]+(1-delta)*k3[i-1])/((1+g)*(1+n))$
c3[i]:=(a*k3[i]^(a-1)+1-delta)*c3[i-1]/((1+g)*(1+n)*(1+rho))$

[a:0.36,n:0.01,g:0.02,delta:0.1,rho:0.02]$
R1:makelist([k1[i],c1[i]],i,0,30)$
R2:makelist([k2[i],c2[i]],i,0,30)$
R3:makelist([k3[i],c3[i]],i,0,13)$

plot2d([[discrete,R1],[discrete,R2],[discrete,R3],x^a-(n+g+g*n+delta)*x,[parametric, 3.894507656030513, t,[t,0,1.7]],[parametric, 1, t,[t,0,0.6]]],[x,0,14],[y,0,1.7],[color,blue,red,green,magenta,gray,gray],[legend, "c(0):0.5", "0.54815","0.58","","",""],[xlabel, "k"],[ylabel, "c"]);

15021201

















 ピンクの曲線は、定常状態 kt+1 = kt+1= k  が成り立つ軌跡(つまり、(20)式)を表しています。k=3.894507656030513 のグレーの直線は、定常状態 ct+1 = c= cが成り立つ軌跡(つまり、(18)式、(19)式が成り立つ k の値です)を表しています。その両者の交点が均衡です。最初の効率労働当たりの資本ストック k が 1 のとき、均衡にたどり着くのは、c の初期値が 0.54815 のときだけ(赤いライン)とわかります。この赤いラインがいわゆるサドルパス(鞍点経路)です。

 それぞれの変数の時間変化も簡単に示せます。離散時間 ( i=0,1,2, ・・・ )とそれぞれの変数で新たなリストをつくって表示させれば、その変数の時間変化を示すことができます。サドルパスを通って均衡に近づいていく場合、c と k の両方が増加することは上の図からもわかりますが、例えば、c の場合、次のプログラムを上記のプログラムに続けて実行させるだけです。

Cdata:makelist([i, c2[i]], i, 0, 30)$
plot2d([discrete,Cdata],[x,-5,30],[y,0,1.5],[xlabel, "t"],[ylabel, ""],[legend, "c"]);

15021203

















 r (資本の収益率、資本からのリターン) を示す場合には、次の r の式にもとづいて、それぞれの k の値から r を求め、それをリストにして表示します。
150213 (17)
 (21)

 上記のプログラムに続いて、次のプログラムを実行させます。
rdata:makelist([i, a*k2[i]^(a-1)-delta], i, 0, 30)$
plot2d([discrete,rdata],[x,-5,30],[y,0,0.4],[xlabel, "t"],[ylabel, ""],[legend, "r"]);
15021204


















(2) n (人口増加率)と g (知識・技術の増加率)が、それぞれ2%と3%から、1%と2%に低下した場合
 つまり、経済成長が低下する場合です。他のパラメーターの値は、上記の場合と同じとします。
 n (人口増加率) と g(知識・技術の増加率)が、それぞれ2%と3%のとき、定常状態における k と c の値は、(19)式、(20)式から、
k*=3.182306091427404
c*=1.037750016414477
 です。しかし、n と g が変化した後に、その初期値から出発すると、新たな均衡にだどりつきません。それは、初期値をその値にして、上記のプログラムを実行させればわかります(下の図の青いライン)。
 
k1[0]:3.182306091427404$
c1[0]:1.037750016414477$
k1[i]:=(k1[i-1]^a-c1[i-1]+(1-delta)*k1[i-1])/((1+g)*(1+g))$
c1[i]:=(a*k1[i]^(a-1)+1-delta)*c1[i-1]/((1+g)*(1+n)*(1+rho))$

k2[0]:3.182306091427404$
c2[0]:1.00365$ 
k2[i]:=(k2[i-1]^a-c2[i-1]+(1-delta)*k2[i-1])/((1+g)*(1+n))$
c2[i]:=(a*k2[i]^(a-1)+1-delta)*c2[i-1]/((1+g)*(1+n)*(1+rho))$

[a:0.36,n:0.01,g:0.02,delta:0.1,rho:0.02]$
R1:makelist([k1[i],c1[i]],i,0,15)$
R2:makelist([k2[i],c2[i]],i,0,30)$

plot2d([[discrete,R1],[discrete,R2],x^a-(n+g+delta+n*g)*x,x^0.36-0.1506*x,[parametric, 3.894507656030513, t,[t,0,1.7]],[parametric, 3.182306091427404, t,[t,0,1.04]]],[x,0,14],[y,0,1.6],[color,blue,red,magenta,green,gray,gray],[legend, "c(0): 1.03775","1.00365","","","",""],[xlabel, "k"],[ylabel, "c"])$

 (k1,c1 が、n と g が低下する前の均衡におけるk と c の値を初期値として出発した場合、k2,c2 が均衡に到達できる経路、つまり、 c を 1.037750016414477 から1 .00365 に低下させて出発した場合)
15021205

















 n=0.02,g=0.03 のときの均衡は緑の曲線とグレーの直線の交点です。そして、n=0.01,g=0.02 に低下したとき、その均衡から出発すると、青いラインの経路にしたがって、左上に移動していきます。経済が新たな均衡に到達するためには、c を低下させなければなりません。まず瞬時的に c を 1.037750016414477 から1.00365 に低下させ、サドルパス(赤いライン)に乗るようにし、そこから均衡に向かっていく、ということがわかります(赤いライン)。
 (この図を見ると、g (知識・技術の増加率) が低下しているのに消費が増えていて変だと思われるかもしれませんが、これは「効率労働当たり」の消費です。「1人当たり」の消費成長率は低下しています。)


 この場合、 r と g (経済成長率)はどのように変化するでしょうか? r (資本収益率)は(21)式を使い、上記の方法で示すことができます。g (経済成長率)は連続時間モデルの場合、
150213 (18)


 となります(知識・技術の増加率 g と区別するために、ここでは gと表記しています)。離散時間モデルでは時間微分が使えないので、最後の項を次のように変えます。
150213 (19)
   (22)


 この式をもとにして、gのリストをつくって表示させます。(ただし、(22)式のかたちからわかるように、i=0 の時点の gのデータが欠けます。グラフに示したとき、その部分が欠けてしまうので、以下のプログラムでは、i=-1 のときの k,c 値 (k(-1),c(-1)) を(19)式、(20)式から逆算して求め、その値をもとにして、i=0 のときの gを求め、リストに加えています。)

rdata:makelist([i, a*k2[i]^(a-1)-delta], i, 0, 30)$
G1data:makelist([i, g+n+a*(k2[i]-k2[i-1])/k2[i-1]], i, 1, 30)$
Gdata:cons([0,0.0429],G1data)$

plot2d([[discrete,rdata],[parametric,t,0.071612,[t,-5,0]],[discrete,Gdata],[parametric,t,0.05,[t,-5,0]],[parametric, 0, t,[t,0.0429,0.05]]],[x,-5,30],[y,0,0.09],[xlabel, "t"],[ylabel, ""],[color,red,red,blue,blue,blue],[legend, "r","","g", "",""]);

15021206

















 前回の記事で紹介したソローモデルの場合と同様に、知識・技術の増加率( g ) と人口増加率( n )が低下すると、r と g のギャップが広がるということがわかります。ソローモデルの場合と同様に、g は瞬時的に下がりますが、r はゆっくりとしか低下しないからです。

 しかし、ラムゼイモデルの場合、ソローモデルに比べて、収束していく速度(r と g の差が一定の値に収束する速度)がソローモデルに比べて速いということもわかります(前回のソローモデルが g+n=8% から 4% への低下で、今回のラムゼイモデルが g+n=5% から 3% への低下で、比較しにくいですが)。

 こうなる理由は、ソローモデルが貯蓄率一定を想定しているのに対して、ラムゼイモデルでは、個人が将来を考え、消費と貯蓄(投資)を調整しているからです。そのために調整が速く行われます。ラムゼイモデルでは、定常状態でも r>g となりますが(ソローモデルでは、定常状態では r=g )、逆にラムゼイモデルでは、定常状態への収束は速くなります(いっぽう、ソローモデルでは、定常状態では r=g となりますが、いったん差ができると、r>g の差がかなりの期間持続する、という結論が得られます。ただし、ソローモデルには個人の最適化が組み込まれていない、という問題がありますが)。

 ピケティは、ソローモデルもラムゼイモデルも、格差を想定していないモデルは現実を説明するときにはあまり信用できないと考えているので、どっちもどっちということでしょう(例えば、こちら)。