こんにちは。農業経済学を勉強している大学院生です。
いきなりですが、経済学の研究をしていて尚且つプログラミングにも興味がある人の中で、
政策評価などで用いられるCGE（computable general equilibrium）やSAM（Social Accounting Matrix）を研究で使った人はいますか？
僕はそういったシミュレーション手法に憧れてはいたものの、いろいろな理由注1から研究で使えませんでした。
しかし、どうしても自分でもやってみたかったので、それならもう研究で使えなくてもいいやと思い直し、趣味として農業経済シミュレーションソフトを作ってみることにしました。
また、普通に作るだけだとつまらないので、最終的にゲームみたいな形にして、STEAMでインディーズゲームとして売れたらいいなと思っています。
最悪、youtubeやniconicoで動画にできたら面白いなぁ（biimシステムみたいな感じとかで、、）。
以上のモチベーションの下でいろいろ試行錯誤していたら、農家経済シミュレーションゲームの型となる農家の効用最大化シミュレーションができたので、モデルの説明なんかも適度に挟みつつ記録を残そうと思います。
環境情報
この記事は
python 3.8.3
scipy 1.5.0
をjupyter notebook上で使ってます。
（この書き方であってますか、、、？）
想定する農家経済モデルの説明
以下、シミュレーションの核となった効用最大化問題をどうやって解いたのか、および農家の行動をどうやって設定したのかを説明しました。
かったるい人や、農業経済的なモデル設定に興味がない人は次の記事に飛んでも大丈夫だと思います。次の回でも最小限の説明は入れますので。
言葉を用いたモデルの説明
農家はどんな農作物を作るのか？
本モデルの農家は、想定した地域の主食となる穀物Aだけを作ります。
そして、自分の食欲を満たすために一定量の穀物Aを食べますし、売ったり買ったりもします。
なお、想定した地域はメキシコ農村部で、穀物Aのイメージは飼料用トウモロコシです。
効用最大化過程で、農家は「主に売る農家」「自給農家」「主に買う農家」に分かれる
本モデルの農家は、自分の効用を最大化する過程で、「穀物Aの生産・農外労働・農地貸借」という3つの収入を得る手段に農地と労働力を配分します。
しかし、効用最大化した状態は、農家の置かれた条件によって「穀物Aを自給したうえで余った分を売る」「自分だけのために穀物Aを作る」「穀物Aを自給して、足りない分を他人から買う」など様々です。
したがって、Janvry(1991)に基づき、上記の3パターンを
net-seller＝「穀物Aを自給したうえで余った分を売る」
自給農家＝「自分だけのために穀物Aを作る」
net-buyer＝「穀物Aを自給して、足りない分を他人から買う」
に分けました。
以上の分類から明らかですが、上記3パターンの農家は「農家」なのに、農業収入を得られるのはnet-sellerだけです。
農家は、消費者として何を消費するの？
本モデルでは、最適化の目的関数として「効用関数」を想定しています。
効用は「余暇」と「その他の財」を消費することによって得られます。
厳密には、「余暇」は自分の所有する時間を配分することで得られ、「その他の財」は購入することで得られます。
なお、穀物Aの消費は効用関数から外し、一定量消費しなければいけないものとして制約を作りました（満腹制約と命名）。
農業以外の収入はどうやって得るの？
net-sellerは自分の農地で穀物Aを作って売るしか収入を得る手段がありませんが、
net-buyerと自給農家は農外労働と農地貸借で収入を得られます。
取引費用も入れました
本モデルでは、net-sellerとnet-buyerが穀物Aを売買するときに、取引費用を課しました注2。
net-sellerは、穀物Aを売る時に (価格 - 取引費用) が価格となり、
net-buyerは、穀物Aを買うときに (価格 + 取引費用) が価格となります。
同様に、net-seller、自給農家、net-buyerが労働力の売買を行う際にも上記と同じように取引費用を課しました。
なんでこれを入れたかはJanvry(1991)を見てほしいのですが、
ざっくりいうと、売買する市場を探したり、穀物Aを運んだりするのにコストがかかるからです。
本来は、市場がちゃんと機能していれば理論的には全ての農家は「net-seller」か「net-buyer」になるはずですが、
こういった取引費用が発生する条件下で「自給的農家」が生まれます。
最適解はどうやって導いた？
ぶっちゃけ、最適解は「計算」してないです。
scipy.optimize の minimize関数を使って、SLSQP法でそれっぽい解を見つけ出しました。
「おいおい、そりゃないだろ」と思った人がいると思うので、ちょっと言い訳させてください。
実は、もともと、効用関数と生産関数についてコブ＝ダグラス型関数を設定して、複数作物を生産する農家モデルを考えてました。
初めは手計算で最適解を求めて、その結果をコードに写そうと思っていたのですが、面倒＋モデル型によっては無理そうだったのでminimize関数で解を探すことにしました。
最適解を探す方法としては、目的関数と制約の両方が非線形になるので、SLSPQ法を用いました。
しかし、どうも最適解が見つけられなかったので、コブ＝ダグラスをやめて、もっと簡単な関数型を使って現在に至ります。
それでは、以下、最適化問題の式を貼りました。
効用最大化問題の式
net-seller
Max
U = a_oC_o+a_ll+{1/2}*(a_{o,o}C_o^2+a_{l,l}l^2+a_{o,l}*C_o*l)...(1)
s.t.
(p_a-t_a)*(Yield*Equip*La-appe)>(w+t_w)Ln+r*An+p_oC_o+MilkMeat...(2:流動性制約) \\ Yield*Equip*La>appe...(3:満腹制約) \\ (Aa+An)/(La+Ln)=Equip...(4:土地装備率の制約) \\ L>La+l...(5:労働力賦存制約) \\ A>Aa...(6:農地賦存制約) \\ 全変数>0
$U$: 効用（満足感）
$a_o, a_l, a_{o,o}, a_{ll}, a_{o,l}$: 2次式モデルのパラメータ注3。
$C_o$: その他の財の消費量（イメージはNetflixや課金スマホゲー）
$L, La, l$: 所有労働力（24-6=18時間）、穀物A生産に用いた労働力、余暇
$A, Aa$: 所有農地面積（2ha）、穀物A生産に用いた農地
$Ln, An$: 雇った労働力、借りた農地
$w, r, t_W, p_a, t_a, p_o$: 賃金（時給）、地代（1シーズン）、農外労働に参加する際の取引費用（取引量当たり）、穀物Aの価格、穀物Aの購入時にかかる取引費用、その他の財の価格
$appe, MilkMeat, Yield, Equip$: appetite（定数、1シーズンの穀物Aの消費量）、ミルクと肉の消費額（定数、1シーズン）、単収（定数、kg/ha）、土地装備率（定数、ha/hour）
自給農家
max
U = a_oC_o+a_ll+{1/2}*(a_{o,o}C_o^2+a_{l,l}l^2+a_{o,l}*C_o*l)...(1)
s.t.
(w-t_w)Ln+rAn>p_oC_o+MilkMeat...(2:流動性制約) \\ Yield*Equip*La>appe...(3:満腹制約) \\ Aa/La=Equip...(4:土地装備率の制約) \\ L>La+Ln+l...(5:労働力賦存制約) \\ A>Aa+An...(6:農地賦存制約) \\ 全変数>0
$U$: 効用（満足感）
$a_o, a_l, a_{o,o}, a_{ll}, a_{o,l}$: 2次式モデルのパラメータ注3。
$C_o$: その他の財の消費量（イメージはNetflixや課金スマホゲー）
$L, La, Ln, l$: 所有労働力（24-6=18時間）、穀物A生産に用いた労働力、農外労働、余暇
$A, Aa, An$: 所有農地面積（2ha）、穀物A生産に用いた農地、貸した農地
$w, r, t_W, p_o$: 賃金（時給）、地代（1シーズン）、農外労働に参加する際の取引費用（取引量当たり）、その他の財の価格
$appe, MilkMeat, Yield, Equip$: appetite（定数、1シーズンの穀物Aの消費量）、ミルクと肉の消費額（定数、1シーズン）、単収（定数、kg/ha）、土地装備率（定数、ha/hour）
net-buyer
max
U = a_oC_o+a_ll+{1/2}*(a_{o,o}C_o^2+a_{l,l}l^2+a_{o,l}*C_o*l)...(1)
s.t.
(w-t_w)Ln+rAn>(p_a+t_a)C_a+p_oC_o+MilkMeat...(2:流動性制約) \\ Yield*Equip*La+C_a>appe...(3:満腹制約) \\ Aa/La=Equip...(4:土地装備率の制約) \\ L>La+Ln+l...(5:労働力賦存制約) \\ A>Aa+An...(6:農地賦存制約) \\ 全変数>0
$U$: 効用（満足感）
$a_o, a_l, a_{o,o}, a_{ll}, a_{o,l}$: 2次式モデルのパラメータ注3。
$C_a, C_o$: 穀物Aの購入量、その他の財の消費量（イメージはNetflixや課金スマホゲー）
$L, La, Ln, l$: 所有労働力（24-6=18時間）、穀物A生産に用いた労働力、農外労働、余暇
$A, Aa, An$: 所有農地面積（2ha）、穀物A生産に用いた農地、貸した農地
$w, r, t_W, p_a, t_a, p_o$: 賃金（時給）、地代（1シーズン）、農外労働に参加する際の取引費用（取引量当たり）、穀物Aの価格、穀物Aの購入時にかかる取引費用、その他の財の価格
$appe, MilkMeat, Yield, Equip$: appetite（定数、1シーズンの穀物Aの消費量）、ミルクと肉の消費額（定数、1シーズン）、単収（定数、kg/ha）、土地装備率（定数、ha/hour）
効用関数について
注3を参考にして、2次式モデルを用いました。
lとCoに掛かっている係数は、それぞれの財単体を消費したときに得られる効用です。
l×Coに掛かっている係数は、（名前なんだっけ、、、）2つを同時に消費したときに得られる効用です。
なんで2つ同時に消費したときの効用があるのかというと、余暇だけあれば幸せな人や、Netflixやスマホゲーだけあれば幸せな人は存在せず、遊ぶ時間と遊ぶ物が同時にないとダメだからです。
生産関数
おそらく、一番問題が大きいのが生産関数でしょう。
一般的に、経済学では、各投入要素の限界効果逓減を仮定します。したがって、その仮定を表現できる関数型がよく使われます。
だから、多分、コブ＝ダグラス型が一番無難だと思ったのですが、上述の通り、その場合は最適解を見つけられなかったです。
コブ＝ダグラス型は凹関数だから問題ないはずですが、多分たくさんぶっこんだ制約式が悪さをしているのでしょう。
ということで、2か月悩んだ結果、上記の感じにしちゃいました。
でも、本当はコブ＝ダグラス型が使いたかった、、、。
パラメータとして使った値はどうやって決めたの？
効用最大化で用いた各パラメータは、実際のデータや参考資料などから抜き出しました。
以下、説明です。
効用関数の係数
効用関数の係数とかは本気でやったらたぶん構造型推定やらないといけないので（やったことない）、
メキシコにおける世帯消費調査結果の消費額シェアを用いました。2財を同時に消費した際の効用は適当にそれっぽい値を入れました。
実は、コブ＝ダグラス型効用関数＋収穫一定を仮定すると理論的にこのやり方でいいのですが、今回は2次式モデルを使ったので厳密には正しい方法ではないはずです。
2次式モデルだと、消費額シェアなどからパラメータは計算できないのかな？ 理論難しいから分かんないです。
具体的には以下のように効用関数の係数を求めました。
メキシコ国家地理統計院（INEGI）が公開する2018年世帯消費調査結果（ENIGH）注4の"trabajos_2018.csv"、"poblacion_2018.csv"、"concentradohogar_2018.csv"を使いました。
求めた値は、①田舎における自営業者以外の平均賃金、②余暇を機会費用として計算した消費額シェア、③生活必需品でない財の消費額シェアです。
Equip（土地装備率）、1シーズンを何日とするか？
Equip（土地装備率）に関するメキシコのデータが見つからなかったので、日本における飼料用トウモロコシの資料を用いました注5。
〇土地装備率（p24）
飼料用トウモロコシ生産の全工程に要する都道府県の10a当たり労働時間
播種関連：1.2時間、追肥除草：0.8時間、収穫：1.6時間、裁断詰め込み：1.6時間
→5.2時間/10a=52時間/ha
→土地装備率：1/52（ha/時間）
〇1シーズンを何日間と仮定するか？（p32）
ここでは、相対熟度が110-120日の品種が仮定されているので、
→栽培日数：4か月＝120日
地代（メキシコpeso）、1シーズン当たり賦存労働力
地代は探しても本当に見つからなかったので、メキシコで行った個人的なインタビュー結果を参考にしました。
インタビューした農家さんは0.1haのグリーンハウスを2000pesoで借りていた。
→1haの農地も2000pesoで借りられると仮定。
1シーズン＝4か月中ずっと農地を借りると仮定して、
4か月あたり地代：2000 * 4 = 8000peso
また、
4か月あたり賦存労働力：(24 - 6) * 30 * 4
Yield（単収）, 穀物A（飼料用トウモロコシ）の価格
メキシコ農牧漁業情報局（SIAP）の公開データを参考にしました注6 。
〇2018年の1年間、天水農業（temporal）、全国
作付面積：396397.92 ha, 収穫量：8324248200 kg, 生産額：5024207380 peso
→価格：5024207380 / 8324248200 = 0.6036 (peso/kg)
→単収：5024207380 / 396397.92 = 12675 (kg/ha)
なお、1か月あたり単位面積当たり生産額：1912 (peso/ha), 1日あたり単位面積当たり生産額：64 (peso/ha)
農外収入時の取引費用、食費、トウモロコシの一人当たり消費量
これもわからなかったので、個人的に行ったメキシコでのインタビュー結果を用いました。
〇1日の一般的な世帯の食費（父、母、息子2人分＝4人の食費）
肉：130 peso/kg → 130/4 peso/0.25kg/人
トルティージャ：15 peso/kg → 15/4 peso/0.25kg/人
ミルク：25peso/L → 25/4 peso/0.25L/人
面倒くさいから、1kgのトルティージャ＝1kgのメイズと仮定。
→1日当たりの肉（0.25kg）とミルク（0.25L）の消費額：130/4 + 25/4 （peso/人）
→1日当たりのトルティージャの消費額：15/4 →一人の人が一日に食べるトルティージャは0.25kgである。と仮定
同じ聞き取り調査の結果、職場までの交通費は大体20pesoだと言ってました。
だから、外部労働の際の取引費用は20peso（てきとう）。
→今更気づきましたが、本モデルの中では20pesoを投入労働力当たりのコストとしましたが、これは固定取引費用とするべきでした。、、、ま、いいか。
その他の消費財の価格
一日当たりの普通の世帯の消費額が200pesoくらいらしいので、
適当に、その他の消費財を一日当たり1つ消費するとして、30pesoとしておく。
さいごに
ということで、条件設定は以上です。
次回、コードをお見せしますが、2点だけ断っておきます。
まず、本シミュレーションはこれが完成形ではありません。
あくまでも、多数の農家が協力/競争する農家経済「ゲーム」を目指しているため、最後はマルチエージェントシミュレーションの形に持っていく予定です。
ですので、コードはマルチエージェントに拡張可能な形でまとめました。
単なるシミュレーションとしては見にくいかもしれないですが、その点ご了承ください。
次に、このシミュレーションでは最適解を計算せずに、Pythonでscipy.optimizeを使ってそれらしい解を無理やり見つけ出してます。
だから、まさかいないと思いますが、このコードを使用して厳密な議論が必要な資料を作る人は気を付けてくださいね。これ、大域的な最適解じゃないかもしれないので。
あと、生産関数の関数型と効用関数のパラメータ推定が多分ダメなやり方です。
最後に、もしも「経済モデルに詳しい方」や「プログラミング言語を使った最適化問題の解き方に詳しい方」がいらっしゃいましたら、
何卒うまいシミュレーションモデルの作り方についてアドバイスを頂ければ幸いです。
それでは、また次の回にお会いしましょう。
注釈
注1 ①よく使われるソフト（GAMSとか）の価格が高い、②研究の意義が見つけられないから手が出しづらい、③専用ソフトの使い方を勉強するのが面倒くさい。
注2 今後、この取引費用が農家の条件によって異なる感じにしようと思ってます。例えば、市場から遠い農家は距離分だけ輸送費用が高くなり、取引費用が高くなるとか。
注3 北村行伸「1999年比較統計システム論、講義資料(?)」http://www.ier.hit-u.ac.jp/~kitamura/lecture/Hit/99Statsys4_c.pdf
注4 https://www.inegi.org.mx/programas/enigh/nc/2018/
注5 菅野勉「飼料用トウモロコシの栽培の現状と今後の作付け拡大に向けた技術開発」https://www.maff.go.jp/kanto/seisan/tikusan/kihon_zyoho/pdf/corn_gizyutu.pdf