今回のデータ包絡問題のお話についてwikiを書く必要ある?と思いつついつか誰かの役に立てばいいかなと思うので書きます.
pythonプログラミングにより複数のインジケーターの組み合わせについて包絡分析を行う.効率性の良いインジケーターの組み合わせを見つける.
Data Envelopment Analysis(DEA)の日本語略で,多入力多出力な複数のシステムの相対的評価を計算し,効率性を求め,それらを比較評価する手法.
データ包絡問題で解く分数計画問題は制約条件を増やすことで線形計画問題として解くことができるらしい.詳しくは下のpdfを呼んでみてほしい.
ということで,以下ではpythonで線形計画法を解いていく.
今回はpython version 3.8.5でやっていく.まったく同じ環境じゃないと出来ないわけではないがpythonのバージョン違いには色々苦しめられたのでなるべく同じバージョン推奨.
pythonのバージョンはコマンドプロンプトで
python -V
で確認できる.
pythonで線形計画問題を解くことのできるpulpパッケージが存在するのでそれをインポートする.
pip install pulp
コマンドプロンプトで上記を実行.さらに動かすプログラムの中でpandasを使用しているのでpandasもインポート.
pip install pandas
環境構築は以上.pulpで線形計画問題を解きたいとき,基本的にpulpさえインポートしてしまえば解けるので便利.
環境構築は済んだのでいよいよ実行に移ります.
実行するpythonファイルは以下.今回はcsvファイルからデータを取ってきているのでそのcsvファイルも置いておく.
クリックでダウンロードできる.
pythonファイルの中身を説明する.
import pulp import pandas as pd
pulpとpandasを使うためにインポート.
dea_df = pd.read_csv('C:/Users/takky/富山県立大学/3年後期/奥原研究室/データ包絡問題/hyouka.csv',encoding="shift-jis",index_col = 0)
in_df = dea_df[['Exposure Time [%]','Max. Drawdown [%]','Avg. Drawdown[%]','Trades','Worst Trade [%]']].copy()
out_df = dea_df[['Equity Final [$]','Equity Peak [$]','Return [%]','Buy & Hold Return [%]','Win Rate [%]','Best Trade [%]','Avg. Trade [%]']].copy()]
指定したcsvファイルからデータを取ってきてin_dfとout_dfに格納.自分のパソコンでpythonファイルを実行するときはcsvファイルのパスを自分のものに変更.
problem = pulp.LpProblem('DEA_problem',pulp.LpMaximize)
線形計画問題の宣言.線形計画問題に名前を「DEA_problem」とつけ,最終的な目標は目的関数を最大化することであるので「pulp.LpMaximize」と指定した.(最小化したいときはpulp.LpMiniimize)とする.
u_in = []
for j in range(num_in):
exec(f"u_{j} = pulp.LpVariable('u_{j}', 0, 50 ,'Continuous')")
u_in.append(u_0)
u_in.append(u_1)
u_in.append(u_2)
u_in.append(u_3)
u_in.append(u_4)
u_in_df = pd.DataFrame(u_in)
u_in_df = u_in_df.T
v_out = []
for j in range(num_out):
exec(f"v_{j} = pulp.LpVariable('v_{j}', 0, 50 ,'Continuous')")
v_out.append(v_0)
v_out.append(v_1)
v_out.append(v_2)
v_out.append(v_3)
v_out.append(v_4)
v_out.append(v_5)
v_out.append(v_6)
v_out_df = pd.DataFrame(v_out)
v_out_df = v_out_df.T
変数の宣言.ひとつの線形計画問題を解くにあたって,その中で同じ名前の変数は使えない(エラー出る)ようなので使いたい変数はすべて宣言する必要がある.わざわざデータフレーム化しているのはあとで変数を扱いやすいようにするためである.
problem += (out_df.iat[i,0]*v_out_df.iat[0,0])
+ (out_df.iat[i,1]*v_out_df.iat[0,1])
+ (out_df.iat[i,2]*v_out_df.iat[0,2])
+ (out_df.iat[i,3]*v_out_df.iat[0,3])
+ (out_df.iat[i,4]*v_out_df.iat[0,4])
+ (out_df.iat[i,5]*v_out_df.iat[0,5])
+ (out_df.iat[i,6]*v_out_df.iat[0,6])
目的関数を宣言.目的関数の宣言はproblem+=(式)でできるので楽.長く見えるのはデータフレームから指定して値を取ってきているだけでたいしたことはしていない.
ここの部分をfor文で回そうとしたがそうすると,目的関数が上書きされてうまくいかなかったのでいちいち指定している.うまいことfor文で回せそうだったら教えてください.
for j in range(num_sitenn):
problem += ((in_df.iat[j,0] * u_in_df.iat[0,0])
+ (in_df.iat[j,1] * u_in_df.iat[0,1])
+ (in_df.iat[j,2] * u_in_df.iat[0,2])
+ (in_df.iat[j,3] * u_in_df.iat[0,3])
+ (in_df.iat[j,4] * u_in_df.iat[0,4]))
- ((out_df.iat[j,0]*v_out_df.iat[0,0])
+ (out_df.iat[j,1]*v_out_df.iat[0,1])
+ (out_df.iat[j,2]*v_out_df.iat[0,2])
+ (out_df.iat[j,3]*v_out_df.iat[0,3])
+ (out_df.iat[j,4]*v_out_df.iat[0,4])
+ (out_df.iat[j,5]*v_out_df.iat[0,5])
+ (out_df.iat[j,6]*v_out_df.iat[0,6])) >= 0
problem += (in_df.iat[i,0] * u_in_df.iat[0,0])
+ (in_df.iat[i,1] * u_in_df.iat[0,1])
+ (in_df.iat[i,2] * u_in_df.iat[0,2])
+ (in_df.iat[i,3] * u_in_df.iat[0,3])
+ (in_df.iat[i,4] * u_in_df.iat[0,4]) == 1
for j in range(num_in):
problem += u_in_df.iat[0,j] >= 0
for j in range(num_out):
problem += v_out_df.iat[0,j] >= 0
制約条件の宣言.制約条件はproblem += (式) >= (式) のような書き方で宣言できる.
すべての宣言が完了したのであとは問題を解くだけ.
problem.solve()
問題を解くのは簡単でproblem.solve()と入力すればOK.
あとは計算した結果をresultリストに格納してprintしているだけ.コードの説明は以上になります.
実行する前にコードのcsvファイルのパスを指定している部分を自分のパスに書き換えてください.書き換えないor書き換えたパスが間違っていると
FileNotFoundError: [Errno 2] No such file or directory:
こんな感じのエラーが出ます.
コマンドプロンプトなどを開いて先ほどダウンロードしたpythonファイルがあるディレクトリまで移る.
cd C:\Users\~
ディレクトリを移動したら
python DEA3.py
で実行できる.
結果を載せます.
今回,データ包絡問題(実際やってることは線形計画問題)を実装してみて,実装自体はうまくいったのですが何分インジケーターの組み合わせの数が多すぎて(37000くらい)このファイルを動かそうと思ったら10日くらいかかることがわかりました。
さすがに今のままだと使い物にならないので今後はいかに実行時間を少なく出来るかといったところが課題になってきます.
いろいろ試してみることがありそうなので今後試してみたときはまたwikiに書き連ねて移行と思います.
データ包絡問題の考え方の中に,いままで解いてきた線形計画問題を主問題とすると変数をいじることで線形計画問題の双対系として解くことができます.詳しくは先ほどのpdfを呼んでみてください.
ということでpythonで双対系として解くコードを実装したので追加する.