2022-05-02

現場に近いデータサイエンスしたい話

漁業体験＠天洋丸

１ヶ月以上前のことですが、大学院の先生から紹介を受けて天洋丸さんの「一年漁師」に参加してきました。データサイエンスと経営の研究をする身として、実際の現場で何が必要とされているのかとても興味があり最高の機会でした。漁業体験を振り返りながら、参加後に考えたことを残しておきます。

「一年漁師」は、寮に滞在しながら天洋丸さんの実際の業務に参加できるプログラムです。煮干しの原材料のカタクチイワシを採る「まき網漁」などを体験できます。僕は短期の滞在でしたが、数週間から数ヶ月の滞在も珍しくないそうです。詳しい内容についてはホームページの方をぜひご覧ください。

tenyo-maru.com

天洋丸さんは、料理の出汁に使われる「煮干し」や、煮干しで育てた鯖を使った「ニボサバ茶漬け」などなど出荷・販売し、漁獲や養殖から食品加工まで一貫して手がける会社です。自宅で「ニボサバ」を漬けにして食べたのですがとても美味しかったです。旨味の強さに脂のノリ具合もちょうど良くて、港町育ちで頻繁に魚を食べて育ってきた僕としも抜群でした。帰宅後1日経って母が捌いたのですが、身が割れる気配が全くないことに驚いてました。お茶漬けは福岡市大濠公園にある&LOCALSで食べれるそうです(^ ^)

参加した理由

大学院で研究しているものは大雑把にいうと経営課題の数理的な解決法¹です。大学・大学院で勉強を進めるにつれて、「実際に役に立つ手法を作ってみたい」とずっと思ってきました。しかし、現場を知らずに作ってしまうと、所詮はただの思いつきです。実際、機会があるときに現場を知る人達に話してみても「かけた時間の割に効果が小さそう」、「一番改善したいのはそこじゃない」という反応ばかりでした。

とはいえ、実際に現場を知る機会はなかなかありません。そういった点で、綺麗な寮があって長期間参加させてもらえた天洋丸さんの環境は最高でした。

参加したもの

初日は漁業体験の参加者向けに実施しているイカかご漁・蛸壺漁に行きました。湾内で揺れも少なく、獲れたイカとタコを加工室で捌いて、その日の夜に寮のみんなと食べました。めちゃくちゃ美味しかったです。「はしり蛸」²と言われるだけあります。

午後には、煮干しで育てる「ニボサバ」の餌やりに同行しました。バッシャバッシャと煮干しを食べるサバも迫力ありましたが、何より横のタイやアジがかなり立派でそちらの迫力もすごかったです。

初日の残りの時間は、「漁網たわし」³を作りました。漁で長年使ってきた網をタワシに再利用する商品です。僕は手先が不器用でこの手の作業にやや気後れしていましたが、みんな丁寧に教えてくれて徐々に覚えることができました。

2日目はロープワークを習いました。縄の太さ、縄にかかる力の大きさ、結ぶ対象によっていろいろな種類があるそうです。もともと知りたかった結び方も教えてもらえて終始興味が尽きない時間でした。

イカ・タコ漁とニボサバの餌付けは漁業体験として半日から参加できるそうです。今回参加した時期は、まき網漁がどちらかというとオフシーズンとのことで、「漁網たわし」作りのほかに、倉庫の整理などをお手伝いさせてもらいました。まき網漁は6月以降は頻繁にあるとのことで、次はそこを狙って参加したいです。まき網漁の様子はくるみ@漁師の旅人さんの投稿がとてもわかりやすいと思います。（「まき網漁が続いた時の私のスケジュール」）

データサイエンスと繋げてみる

天洋丸さんの漁業体験に参加しつつ、聞いたお話を参考にして何ができそうか考えてみました。ただ、冒頭でも言った通り、素人の妄想に過ぎないのは変わりないかもしれませんがとにかく考えてみました。

何ができそう？

僕が大学院で勉強している分野は、「何らかの計画」をうまくやるための手法を構築してきました⁴。具体的には、輸送や在庫などの配分を調整するときに応用されています。

例えば、燃料の状況を加味した運用計画を自動で策定した関西電力の例などがあります⁵。火力発電にはさまざまな種類の燃料が必要でそれぞれ使うべきタイミングがあるそうです。ベストなタイミングで使えるように在庫量や輸送の計画を調整する役目は、長年の経験を積んだ熟練者しか担えませんでした。それを自動で調整するソフトを開発したのがこの例です。もちろん導入するには開発コストがかかるというデメリットが生じますが、

熟練者以外でも計画の策定ができるようになる
計画の策定が早くなり、よく検討できるようになる

という点でメリットが優ったようです。しかし、事業の規模や質によって、必ずしもメリットが優るとも限らないので、そのあたりはよく検討する必要があると思います。さて、今回は原材料など「モノ」ではなく、人員配置のフレームワークを考えてみたいと思います。

非定型業務の記録とスケジューリング問題

「その日の仕事を数字として記録したいけど、毎日いろいろな種類の仕事をこなしていく必要があり、それらを一つ一つ記録するのは難しい」というお話を天洋丸さんで聞きました。記録さえできれば、どういう仕事を、どのくらいの時間で出来たか検討し、スケジューリングや人員配置の改善につながるはずです。今回はざっくりと、「より良い仕事の割り振り方を得る」というのを目標にして、記録の方法を自分なりに検討して、メリットとデメリットを整理してみたいと思います。

問題設定

ベテランのAさん、中堅のBさん、新人のCさんの3人が毎日さまざまな業務をこなしていきます。ある日の業務は以下であるとします。

倉庫の片付け
出荷
雑貨製作

従業員ごとに時給は変わるので、割り振り方を工夫してできるだけ安く業務が終わるようにしたいです。まず、1人で業務を担当した場合、それぞれの業務がどのくらいで終わりそうか（ふわっと）予測します。今回の例はAさん、Bさん、Cさんになるにつれて時間がかかるイメージです。

さて、少し考えやすくするため、上記から「1人だと1時間あたりだと何%終わるか」に捉え直します。もし1人で2時間かかる仕事があれば、1時間で50%すすむというわけです。2人いれば1時間で終わそうです。

以下の表がAさん、Bさん、Cさんの値です。Aさん、Bさん、Cさんの順で1時間あたりに終わる割合が小さくなります。

次に、仮に時給が以下のようだとします。仕事が早さに対応して、Cさん、Bさん、Aさんになるにつれて時給が高いです。

人件費を最小化するように割り振る

ここで一旦、話を整理します。人件費を最小化する仕事の割り振り方を見つけたいです。ただし、以下のような条件がつきます。

それぞれの労働時間は8時間以内
その日のうちにそれぞれの仕事を100%終わらせる

この時、その日の人件費を最小にする仕事の割り振り方がエクセルの分析ツール「ソルバー」⁶で（ある程度）わかります。結果は以下のようになります。

Aさんは出荷に3時間、Bさんは雑貨製作に2時間、Cさんが倉庫の片付けと雑貨製作に6時間と2時間ずつ割り振られれば、人件費の総額が最小化されるということです。ちなみに人件費は以下のようになります。

手作業でしようとするとそこそこ複雑ですが、エクセルで瞬時に解決できてしまいました。

そんなうまくいく...？

現実がこのような教科書通りの問題設定なら、ここまでの分析によって人件費が最小化することが出来ます。しかし、今回は少なくとも以下が成り立っていると仮定しています。

人数が増えるのに比例して仕事が終わる時間が減少する。
仕事場の間に距離がない。
仕事をどの順序でしても時間が変わらない。
有資格者の参加が必要なタイプの仕事がない。

(1.)は例えば、1時間で30%進められる人と、1時間で50%進められる人が一緒にやれば、1時間で80%進むという仮定です。しかし、現実は必ずしもそうではなく、3人入れる倉庫を片付けるときに100人投入しても、そこまで早くならなそうです。(2.)は、仕事場の間に距離があると移動に時間がかかるため、できるだけ1人が1つの仕事を担当すると時間がかからないという点を本来は考慮するべきです。(3.)は、体力を使う仕事の場合、その日の前半でやるのと後半でやるのでかかる時間が変わりそうです。(4.)は、操船など資格が必要で一部の人しかできない仕事を仮定していないということです。

逆にいうと、これらを分析に組み込むことで、より現実に沿った結果を得られそうです。というより、これらの検討事項に至るまでが一番難しい気がします。上記のほかに検討するべき仮定があるかもしれません。むしろ、まだまだたくさんある気がします。

何を記録するべきか

今回の例で必要なデータはAさん、Bさん、Cさんが各仕事にかける時間です。とはいえ、日ごとに違う仕事を細かく分類して記録するのは手間なので、ざっくり分けて記録するのがいいかもしれません。例えば以下のように毎日の仕事を分類して3人がどれだけ時間をかけたか記録すると手間が少ないです。

場所を問わずに「片付け」
対象を問わずに「輸送」
商品を問わずに「加工・製作」

あとはこれらの記録と、経験者の勘を頼りにそれぞれが各仕事にかかる時間を決めるといいかなと感じています⁷。あるいは、追加の仮定を検討すると、必要なデータが増えるかもしれません。その辺の塩梅は、実際に試しながら検証しつつやるしかない気がします。ただし、これはあくまで「人件費を最小化する」を目的としたときの1つのやり方に過ぎません。目的が違えば、手段も変わって必要なデータも変わるかもしれません。何を目的にするかというのも、現場でのコミュニケーションを経て考察するべきかなと感じています。

まとめ

天洋丸さんで業務に参加させてもらったのを振り返りつつ、データサイエンスで何ができそうか考えてみました。やはり、付加価値を生む施策は「現場でのコミュニケーションなしには辿り着かない」と改めて感じました。話を聞くだけでなく、実際に経験しないと何が問題かわかりませんし、業務中に聞く方が詳細に知れることが多いです。一次産業のDXに興味ある院生の皆さん、天洋丸さんまでぜひお問い合わせください。

www.instagram.com

「問題設定」の文章冒頭を「人件費を最小にする」という趣旨に修正しました。5/8, 2022

雑すぎてすみません。具体的にはオペレーションズ・リサーチで「応用」してます。↩
天洋丸走り蛸セット ↩
天洋丸の網エコたわし ↩
詳しくはこちらをご覧ください：オペレーションズ・リサーチとは｜公益社団法人日本オペレーションズ・リサーチ学会 ↩
石炭火力発電所の燃料運用最適化 | DeNA×AI WORKS ↩
分析したエクセルのリンク ↩
必要な精度が出れるなら勘が大部分でも良い気がします。↩

2022-01-18

生産・サービスの効率性を測る包絡分析法（DEA）：使い方を具体例で確認してみる

オペレーションズ・リサーチ python

生産・サービスの効率性を測りたい

包絡分析法（DEA）を使うと、以下のような利点があります。

包絡分析法を使えば、生産活動やサービス提供の総合的な効率性を「効率値」という1つの指標で判断できる。
現状で、無駄になっている材料の量や、増やせそうな生産・サービスの量も数字で把握できる。

箇条書きの最初は、言い換えると、材料や設備を複数つかったとしても「効率値」という1つの指標で効率性を表すことができるということです。なぜこれが可能なのかというのを理解するのは少々時間を要しますが、ひとまず結果を解釈するのはそこまで難しくありません。それでは分析の流れを見ていきます。

分析例

分析用のデータを準備

3つの材料と2つの設備、2人の人員で完成品を作る例を考えてみます。

各回で費やされた材料や設備使用時間、作業時間が横向きに並んでいます。No.1の行を横に見ていくと、材料１を10、材料２を1・・・とあります。その行の最後には完成品が3とあります。これは、材料１を10、材料２を1・・・を使って、完成品を３作ったということです。これが回ごとに下向きに記載されています。

包絡分析法の結果

これに対して、包絡分析法を実施した結果が以下です。

この出力結果を使って各生産の効率性の優劣を比較していきます。解釈の仕方は以下です。

上図に従って、まずは効率値が1の生産から確認していきます。出力結果によると第1, 2, 3, 4, 8, 10回の生産は効率値が1です。

次にラムダ変数の値を確認すると、それぞれの生産のラムダ変数は自らのラムダ変数（下の図は赤丸）以外は０です。つまり、第1, 2, 3, 4, 8, 10回の生産は効率的で改善の余地はない、という解釈になります。

一方で、効率値が1未満の第5, 6, 7, 9回の生産は非効率的です。

効率値の値が小さい順に第6, 7, 9, 5回となっており、第6回の生産が最も非効率的、という解釈になります。ちなみに、包絡分析法では「非効率な第5, 6, 7, 9回の生産については、材料と設備使用時間と作業時間をそれぞれ(効率値)倍の値まで削減できたはず」という解釈も得られます。つまり、第5, 6, 7, 9回の生産の材料と設備使用時間と作業時間は、各回ごとに0.875倍, 0.6倍, 0.769倍, 0.846倍までに本来は抑えられたはず、ということです。

注意点

包絡分析法の使い方を簡単に触れてきました。最初に述べたとおり、この分析法は生産・サービスの総合的な効率性を一つの指標で比較できる利点があります。しかし、いつでも使えるわけではありません。例えば、今回のモデルは「材料や設備、作業時間に比例して完成品の数が増えるはず」という仮定のもとで成り立つ点に注意する必要があります。このモデルを特にCCRモデルと言いますが、「材料や設備、作業時間を増やせば増やすほど、完成品の増加率が減少する」という仮定を採用したBCCモデルなど別のモデルがあり、生産・サービスごとに適切なモデルを選ぶべきです。

また、今回の計算例ではありませんでしたが、効率値が１でも非効率なこともあります。（その場合は、スラック変数の分だけ削減可能です。）

まとめ

包絡分析法（DEA）を使って、効率値を使ってまとめて効率性を比較する例を示しました。非効率な活動を特定するだけでなく、その改善に対して具体的な目標値を得ることができます。さらには、効率値を悪化させる要因を調べることで、製品・サービスの製造や提供により長期的な安定性が期待できます。

参考文献

刀根薫（1993）『経営効率性の測定と改善ー包絡分析法DEAによるー』日科技連.

データ

ブログ_DEA_サンプル_構想 (1).xlsx - Google スプレッドシート

ブログ_DEA_16Jan.xlsx - Google スプレッドシート

計算用コード

Google Colab

2022-01-14

生産・サービスの効率性を測る包絡分析法（DEA）：利点と仕組みを具体例で考えてみる

オペレーションズ・リサーチ

効率性の分析で用いられる通常の手法

何らかの生産やサービスにおける効率性を比較するとき、材料の量や作業時間と、最終的な製品やサービスの量を比にした指標で検討することが多いと思います。例えば以下のような指標です。

業種	効率性の指標例
木造アパート	面積あたりの木材量面積あたりの作業時間面積あたりの道具・設備の使用時間
自動車生産	１台あたりの直接作業時間工場あたりの直接作業時間面積あたりの道具・設備の使用時間
コーヒーショップ	一杯あたりのコーヒー豆(g) 一杯あたりの作業時間家賃あたりの提供数

通常の分析：うまくいきそうなケース

木造アパートを例にとって考えてみます。アパートを建てるとき、より少ない木材で完成させるのが望ましいはずです。しかし、木材を破損させたり、加工に失敗したりすると想定より多くの木材が必要です。このとき、「面積あたりの木材量」という指標を使えば、各物件ごとに必要だった木材の量に関して効率性を比較できそうです。同じように「面積あたりの作業時間」、「面積あたりの道具・設備の使用時間」も、物件ごとに必要だった作業の時間や、道具・設備を使用した時間を比較できそうです。このように、生産やサービス提供の効率性を測ろうとすると、通常は複数の指標で判断することになります。ここで、4つのアパートを建て、以下のような結果が得られたとします。

	物件A	物件B	物件C	物件D
面積あたりの木材量	100	200	140	180
面積あたりの作業時間	10	20	18	14
面積あたりの道具・設備の使用時間	6	12	10	10

上の表で、各指標の優劣を物件ごとに比較してみます。すると下図のようになります。

f:id:Noooooz:20220114014158p:plain

すべての指標で、物件Aが他の物件より優れています。つまり、物件Aが最も効率的であると判断できます。新たにアパートを建てるときは、物件Aの指標を達成できるように対策を考えれば良さそうです。

通常の分析：うまくいかないケース

ところが、4つの物件を建てた結果が以下のような場合はうまくいきません。

	物件A	物件B	物件C	物件D
面積あたりの木材量	100	140	180	200
面積あたりの作業時間	20	18	10	12
面積あたりの道具・設備の使用時間	12	16	14	6

前の結果と同じように、各指標の優劣を物件ごとに比較してみます。

f:id:Noooooz:20220114014203p:plain

指標を一つずつ確認していきます。まず、「面積あたりの木材量」で最も小さい値なのは物件Aです。同様に「面積あたりの作業時間」と「面積あたりの道具・設備の使用時間」で最も小さい値はそれぞれ物件Cと物件Dです。つまり、指標ごとに最も優れる物件が異なっており、どの物件が最も効率的か一概には言えない状況です。

通常の分析：妥協案と問題点

このままでは、今後の新規物件で目指すべき最も効率的な値がわかりません。とりあえず上記の結果に基づいて分析を進めるため、試しに解決策を講じてみます。

それぞれの指標で最も良い値を目標にする。

差し当たって、表中にある一番良い値を各指標の目標にしておくという考え方です。しかし、よく考えると、この目標はそもそも達成不可能かもしれません。なぜなら、どれかの1つの指標を最も良い指標に近づけようとすると、他の指標が悪化するからです。木造アパートの例では、加工でミスしないように十分に注意すると「面積あたりの木材量」を改善できそうです。一方で、ミスしないように注意すると加工に時間がかかります。すると「面積あたりの作業時間」と「面積あたりの道具・設備の使用時間」が増加してしまいます。要は、ある１つの指標で最も良い値を達成しようとすると、他の指標を犠牲にする必要があり、上記の目標はそもそも達成可能ではないということです。それにもかかわらず改善のために予算や時間を無駄に費やしてしまうかもしれません。

指標の順位が平均的に最も優れた物件を、最も効率的だとする。

このとき、3つの指標の平均順位は各物件ごとに以下のようになります。

	物件A	物件B	物件C	物件D
指標の平均順位	2.33	3	2.33	2.33

同率一位が3つ出てきました。このとき、新規の物件で具体的にどの物件の指標を目指すべきかわからず、分析としてはあまり意味がありません。

木造アパートの例に限らず、冒頭の表にあった自動車生産やコーヒショップなど他の業種でも同じです。複数の指標を使って効率性を分析しようとする限り、以上のような状況と問題点を回避は困難です。

包絡分析法（DEA）

包絡分析法の利点

ここで有用なのが包絡分析法（DEA）です。包絡分析法は「効率値」と呼ばれる指標が効率性の程度を示します。この分析法の利点は簡単にいうと以下です。

「効率値」という1つの指標のみを使うため、効率性の優劣が容易に判断できる。
ある材料や設備の量で、生み出せる製品やサービスの最大の量がわかる。
ある製品やサービスの量を生み出すために必要な材料や設備の最小の量がわかる。

効率性を比較できる仕組み

それでは、なぜ上記のようなことがわかるか、その概要を説明します。まず、以下の図をみてみます。

f:id:Noooooz:20220114010401p:plain

左の図は、木材の量が増えるほどそれに応じて部屋数も比例して増えています。これはどの物件も想定通りの木材の量で済み、木材の量が多い物件ほど、部屋数もしっかり多くなっている理想的な状況です。一方、右の図は木材の量が増えても必ずしも部屋数が増えていません。物件Bや物件Cは左図と比べて下側に点があり、加工のミスなどによって木材を浪費し、想定通りの部屋数に到達しなかったケースです。

上の左図では、物件Bや物件Cは本当は点線上にあるべきで、点線上にあるとき最も効率的に木材を使えています。逆にいうと、点線までの距離が大きければ大きいほど非効率です。そこで、「点線までの距離」を効率性の指標にすると、距離の大小がそのまま各物件の効率性の優劣として判断できそうです。この考え方に基づき、さまざまな材料や人員、設備を用いて製品やサービスを生み出す活動の、総合的な効率性を示す指標を「効率値」として定めます（下図）。

f:id:Noooooz:20220114014522p:plain

また、図からは以下のようなこともわかります。下図の左で、物件B、物件Cから真上に伸ばした矢印と点線が交わる点は「同量の木材で到達可能な部屋数」を示しています。また下図の右では物件B、物件Cから左へ水平に伸ばした矢印と点線が交わる点は「同じの部屋数で最小限の木材の量」を示します。

f:id:Noooooz:20220114014556p:plain

ここまでの例では、部屋を作るのに必要な材料が木材のみのケースですが、実際はさまざまな種類の材料や、人員、道具や設備が必要になります。そのような場合でも、包絡分析法は効率値というただ1つの指標のみで効率性を比較することができます。これが、この分析手法の大きな利点です。実は、「点線までの距離」で効率性を測ることでそれが可能となっています。ただその仕組みは、やや技術的すぎるのでここでは省きます。

まとめ

振り返り

当記事は以下のような流れでした。

効率性を分析するとき、複数の指標を使って行われることが多い。しかし、具体的な改善案がはっきりしない。
包絡分析法を使えば、総合的な効率性を「効率値」という1つの指標で判断できる。
さらに、無駄になっている材料の量や、増やせそうな生産・サービスの量も数字で把握できる。

包絡分析法（DEA）の限界

包絡分析法にも限界があります。たとえば、効率値だけで、「効率が悪くなる原因」までは特定できません。人為的なミスの他に、外的な要因によっても生産・サービスの効率性は悪くなります。つまり、たとえ包絡分析法が示す効率値が悪かったとしても、その原因が必ずしも改善できるとは限りません。

その生産・サービスの効率性に関係する要素が少なければ、経験則等でそれらの関係がはっきりするかもしれません。しかし、影響を与える要因が多かったり、そもそもノウハウがなかったりすると、効率値を出したところで「何が原因か」まではわかりません。解決策としては、「何が原因か」を調べる別の統計的手法を使ったりします。

最後に

今回は木造アパートを例として説明しましたが、他に例として挙げられそうな題材がありましたらぜひお声掛けください。

参考文献

刀根薫（1993）『経営効率性の測定と改善ー包絡分析法DEAによるー』日科技連.

より具体的な使い方はこちら↓ https://noooooz.hatenablog.com/entry/2022/01/18/020341

2021-11-14

残余利益とお友達になりたい

会計学ファイナンス python

株式価値評価における３つのアプローチ

配当割引モデル
割引キャッシュフロー・モデル
残余利益モデル

　残余利益がしっくりこないので、かなりざっくりですが手を動かして追ってみたいと思います。具体的には、企業の財務諸表を拾ってきて実際に計算してみるだけです。財務諸表は大変ありがたいことに伊藤忠商事株式会社*1が時系列データをエクセルで公開してくれてますのでそちらを使いました。

残余利益モデル*2

　株式の価値が、現在の自己資本額と将来期間の利益の現在価値を合計したものになるという考え方です。以下のように表されます。

$\text{将来の残余利益}=[\text{１株当たり予測利益}-\text{１株自己資本}\times\text{自己資本コスト}$ ]

　当期純利益がすべて配当に回すとき、以下のようになります。

$\text{自己資本の価値}V=\text{自己資本簿価}C_0 + \frac{\text{当期純利益}A-\text{自己資本コスト}r_e\times \text{自己資本簿価}C_0}{\text{自己資本コスト}r_e} \\ (\text{ただし、自己資本簿価}C_0=\text{株主資本}+\text{その他の包括利益累計額})$

実装

　それではあくまでざっくりですが計算してみます。基本的には、前節の式を反映させるだけですが自己資本コストだけは、前節で引用した文献にある例で使われていた値を流用します。したがって、関係ない他企業の自己資本コストを使っていることになるので、推定結果自体は無意味です。ご注意ください。

import pandas as pd
import numpy as np

df_BS = pd.read_excel('財務諸表エクセルファイルのパス').fillna(0)
df_PL = pd.read_excel('財務諸表エクセルファイルのパス').fillna(0)


#変数を定義
nv = df_BS['株主資本合計']
op = df_BS['その他の資本の構成要素合計']
C_0 = nv + op

A = df_PL['当社株主に帰属する当期純利益']
r_e = 0.0877
V = C_0 + (A- r_e * C_0) / r_e

df_rslt = pd.DataFrame(V.values,columns=['V'], \
                       index=['2020',\
                              '2019',\
                              '2018',\
                              '2017',\
                              '2016',\
                              '2015',\
                              '2014',\
                              '2013',\
                              '2012'])

# 財務諸表の単位が百万円なので戻す
df_rslt['V'] = df_rslt['V']*1000000

# 2012-2020の発行済総株式数を有価証券報告書から取ってくる,千株単位なので戻す
df_rslt['shareNum'] = [1584889, \
                       1584889, \
                       1584889, \
                       1662889, \
                       1662889, \
                       1662889, \
                       1584889, \
                       1584889, \
                       1584889]
df_rslt['shareNum'] = 1000*df_rslt['shareNum']

#一株当たりの株式価値
df_rslt['stockValue'] = df_rslt['V'] / df_rslt['shareNum']
df_rslt.to_excel('保存先のパス')

f:id:Noooooz:20211113235807p:plain — 年ごと・株価推定（最右列）

*1:https://www.itochu.co.jp/ja/ir/finance/financial_data/index.html

*2:桜井久勝（2020）「財務諸表分析」中央経済社. 第14章

2021-11-10

倒産確率推定の話②

ファイナンス python

アプローチ②：オプション・アプローチ

森平（2000）*1を主に参考にして勉強しました。将来資産が将来負債額を下回る「確率」を推定する、とざっくり捉えてます。これが、コールオプションの価格決定理論に基づいているため「オプション・アプローチ」と呼ばれたりします。同論文では以下において、 $E_t,\sigma_E,\mu_D,D_T$ は所与の値で、 $\mu_A, \sigma_A, A_t$ をガウス・ザイデル法で求めるよう書いています。後述する実装例では反復回数50のガウス・ザイデル法を用いました。

$\mu_A = (\frac{E_t}{A_t})\mu_E+(1-\frac{E_t}{A_t})\mu_D\\ d_1^* \equiv \frac{ln(A_t/D_T)+(\mu_A+\sigma^2_A/2)T}{\sigma\sqrt{T}} \\ \sigma_A = \sigma_E\frac{E_t}{A_t N(d_1^*)} \\ d_2^* \equiv d_1^* - \sigma_A\sqrt{(T)} \\ A_t = \frac{E_t+D_T\cdot N(d_2^*)}{N(d_1^*)} \\$

各変数が意味するところは以下通り。

$\mu_A$ は株式価格の期待値
$\mu_E$ は負債成長率の期待値
$E_t$ は株式（自己）資本総額（株価 $\times$ 発行済み株式数）
$A_T$ は資産価値
$D_T$ は負債価値

1994年から2016年までのシャープの推計を行うと以下の１枚目の画像のようになりました。実装例には含めていませんが、ある程度の幅をとってシミュレーションした結果が２枚目の画像です。

f:id:Noooooz:20211110005246p:plain — 年ごとの倒産確率の推移

f:id:Noooooz:20211110005322p:plain — シミュレーション例

シャープの財務データを使ったPythonの実装例は以下です。1994年から2016年の各変数に必要な財務データを有価証券報告書から見つけて、おそらく手作業になるけど、csvファイルにすると使えるはず。（あるいは所属大学や組織の商用データベースで。）

#This code calucurating 'Option-based Bunkruptcy risk' from time t to T=1.
#'T = 1' means 1 year later from time t.
#The Expected growth-rate of asset(= marketcapitalization) debt
#are calucurated with datas for last five years.

#Through the process below, calucurating the parameters; mu_A, sigma_A, A_t.
#Other parameters are given; E_t, sigma_E, mu_E, mu_D, D_T, T.
#The process uses Gauss-Seidel method.

#Company: SHARP

import pandas as pd
import numpy as np
from scipy.stats import norm
import math
import matplotlib.pyplot as plt


#Reading time-series data of debt and equity and maeket capitalization
#as Dataframe of Pandas
sharp_debt_pd = pd.read_csv('シャープの負債額csvの絶対パス', header = None)
sharp_equity_marketcap_pd = pd.read_csv('シャープの株式資本総額csvの絶対パス', header = None)

#Chaging datatypes of number as 'float'
sharp_debt_pd[2] = sharp_debt_pd[2].astype(float)

#Debt is given from the company's B/S, 
#so the number of figure should be fixed.
sharp_debt_pd[2] = sharp_debt_pd[2] * 1000000
print('Debt of sharp 1994-2016')
print(sharp_debt_pd)
print('-------------------------')


sharp_equity_marketcap_pd[1] = sharp_equity_marketcap_pd[1].astype(float)
sharp_equity_marketcap_pd[2] = sharp_equity_marketcap_pd[2].astype(float)
print('Market capitalization (=0) and equity price (=1) 1994-2016')
print(sharp_equity_marketcap_pd)
print('-------------------------')


#Calucurating change rates per yaer
#Changing the names of Dataframes
#Be careful of not forgetting first line would be NaN
#on 'changes_..._pd'
changes_debt_pd = pd.DataFrame()
changes_equity_marketcap_pd = pd.DataFrame()

changes_debt_pd[0] = sharp_debt_pd[2].pct_change()

changes_equity_marketcap_pd[0] = sharp_equity_marketcap_pd[1].pct_change()
changes_equity_marketcap_pd[1] = sharp_equity_marketcap_pd[2].pct_change()

final_result = list()



print('Changes of debt')
print(changes_debt_pd)
print('-------------------------')
print('Changes of market capitalization (=0) and equity price (=1)')
print(changes_equity_marketcap_pd)
print('-------------------------')



#Defining parameters gotten given from the time-series data upper at t.
for t in range(8,23):

    #'iloc' should be chosen to specify component of dataframe with number.
    #There is a explaination for it in Japanses.
    # https://yolo.love/pandas/loc-iloc-at-iat/

    #E_t is the market capitalaization at time t.
    E_t = sharp_equity_marketcap_pd.iloc[t, 1]

    #Calucurating a sigma of equity for last 5 years.
    sigma_E = changes_equity_marketcap_pd.iloc[t-5:t, 1].std()

    #Calucurating a mu of equity price change for last 5 years.
    mu_E = changes_equity_marketcap_pd.iloc[t-5:t, 1].mean()

    #Calucurating a mu of debt change for last 5 years.
    mu_D = changes_debt_pd.iloc[t-5:t, 0].mean()
    

    #D_T is amount of debt after 1 year.
    #[let a debt on book (at t = T) be regard as effectivness for 1 year]
    D_T = sharp_debt_pd.iloc[t, 2]

    #T=1 means measuring a risk of bankruptcy of 1 years later.
    T = 1
    print('E_t is ' + str(E_t))
    print('sigma_E is ' + str(sigma_E))
    print('mu_E is '+ str(mu_E))
    print('mu_D is ' + str(mu_D))
    print('D_T is ' + str(D_T))


    #Creating equation based Gauss-Seidel method 
    # with the parameters defined upper.

    #Number of times of roop to get convergence value(=answer at t) is n.
    #First, initial value should be created.

    #The initial values are be calucurated with historical data to
    #avoid a result of error caused by limitations.
    n = 50
    mu_A = np.array([0.25493648])
    sigma_A = np.array([1.24735172])
    A_t= np.array([1.348580e+12])

    #Creating a empty pd.DataFrame to show final results at each t.
    final_results = pd.DataFrame()

    
    #roop start
    for k in range(0,n+1):
        
        d_1 = (np.log(A_t[k] / D_T) + (mu_A[k] + (sigma_A[k])**2/2 ) * T) / (sigma_A[k] * np.sqrt(T))
        sigma_A_ = sigma_E * E_t / (A_t[k] * norm.cdf(x=d_1, loc=0.5, scale=1) )
        sigma_A = np.append(sigma_A, sigma_A_)

        d_2 = d_1 - sigma_A[k+1] * np.sqrt(T)
        A_t_ = (E_t + D_T * norm.cdf(x=d_2, loc=0.5, scale=1)) / norm.cdf(x=d_1, loc=0.5, scale=1)
        A_t = np.append(A_t, A_t_)

        mu_A_ = (E_t/A_t[k+1]) * mu_E + (1 - E_t/A_t[k+1]) * mu_D
        mu_A =np.append(mu_A, mu_A_)
    
    print(sigma_A[n])
    print(A_t[n])
    print(mu_A[n])
        
    

    d = (np.log(A_t[n] / D_T) + (mu_A[n] + (sigma_A[n])**2/2 ) * T) / (sigma_A[n] * np.sqrt(T))
    bankruptcy_risk_t = 1 - norm.cdf(x=d, loc=0.5, scale=1)

    print(d)
    print(bankruptcy_risk_t)
    final_result.append(bankruptcy_risk_t)
        
print(final_result)

y = final_result
x = list(range(2001, 2016))

plt.plot(x, y)
plt.show()

*1:森平爽一郎（2000）「信用リスクの測定と制御」計測と制御. （2.3 オプション・アプローチ）信用リスクの測定と制御

2021-11-10

倒産確率推定の話①

ファイナンス

倒産確率推定、まとめておきたい

　学部３年、ファイナンス・リスクマネジメントのゼミで半年かけて倒産確率についてグループで発表しました。最近になって、会計方面から改めて倒産リスクを勉強したこともあり全体的に振り返ってみたくなりました。 AltmanのZスコアとMertonモデルについてまとめておきます。

Zスコア

　前者は財務データのうち倒産の予測に寄与しそうなものを説明変数として厳選し、倒産予測を行おうとするものです。倒産の予測に寄与しそうなものを選定する過程で、財務データそれ自体への理解が深くなければならないので会計学ベースのアプローチという印象を持っています。

Mertonモデル（オプション・アプローチ）

　一方、後者はブラック・ショールズモデルを基礎とするオプション評価式に基づいた予測モデルです。数理ファイナンス寄りのアプローチという印象です。それでは、Zスコアから見てみます。 f:id:Noooooz:20211124200246p:plain

アプローチ①：Zスコア

　例えば、財務情報から５つの説明変数を構成し判別分析を行います。具体的には以下のような説明変数です。*1

$X_1 = \frac{\text{working capital (運転資本)}}{\text{total assets (総資産)}}$
$X_2 = \frac{\text{retained earnings (留保利益)}}{\text{total assets (総資産)}}$
$X_3 = \frac{\text{earnings before interest and taxes (利子・税金の控除前利益)}}{\text{total assets (総資産)}}$
$X_4 = \frac{\text{market value equity (株式時価総額)}}{\text{book value of total liabilities (総負債)}}$
$X_5 = \frac{\text{sales (売上)}}{\text{total assets (総資産)}}$

　このもとで、Altmanは以下のような判別関数を推定しました。

$Z = 0.012X_1 + 0.014X_2 + 0.033X_3 + 0.006X_4 +0.999X_5$

　ここで、 $Z=2.675$ を閾値として、これを下回ると倒産と判断する。ただし、検索で出てきた論文を見ると判別分析のモデルによって切片があったりなかったりすると述べています。それと、財務分析の入門書に平易な例による解説があります。*2。その例でも、ある程度高水準の判別力が得られるという説明がなされています。ただ、判別力の評価手法についてはいろいろ議論があるそうです。実際、あくまで印象ですが、2010年代の倒産確率の格付け予測の論文では少なくともAUCでの評価でモデルの性能を測っていると思います。

*1:E. I. Altman, https://pages.stern.nyu.edu/~ealtman/Zscores.pdf

*2:桜井（2021）「財務諸表分析」中央経済社. 第15章