バックナンバーはこちら。
https://www.simulationroom999.com/blog/stock-predict-matlabpython-backnumber/
はじめに
前々回、前回でMATLAB版の売却、買付タイミング時のVTI単価特定コードと、
それの実行結果を確認。
いい感じに取得出来ていそうではある。
今回からこれらのPthon側の作業となる。
登場人物
博識フクロウのフクさん
イラストACにて公開の「kino_k」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=iKciwKA9&area=1
エンジニア歴8年の太郎くん
イラストACにて公開の「しのみ」さんのイラストを使用しています。
https://www.ac-illust.com/main/profile.php?id=uCKphAW2&area=1
売却、買付タイミング時のVTI単価特定コード(Python版)
太郎くん
MATLABで売却、買付タイミング時のVTI単価特定はできたから、今度はPython版だね。
フクさん
作って来たコードを見せよう。
import numpy as np
import matplotlib.pyplot as plt
VTI=np.loadtxt('VTI3.csv',delimiter=',') # 変換用波形読み込み
N=len(VTI)
L=N/2
x=np.linspace(-L,L-1,N)
ft=VTI-np.mean(VTI)
fig = plt.figure()
# 変換前波形
ax1 = fig.add_subplot(5, 1, 1)
ax1.plot(x,ft)
ax1.set_title('f(t)')
ax1.grid()
# FFT後にローテーション
Fw=np.fft.fft(ft)
Fw_tmp = np.roll(Fw,int(L))
ax2 = fig.add_subplot(5, 1, 2)
mask = (-50<=x) & (x<=50)
ax2.plot(x[mask],np.abs(Fw_tmp[mask]))
ax2.set_title(r'$F(\omega)$')
ax2.grid()
# 拡大
ax3 = fig.add_subplot(5, 1, 3)
mask = (0<=x) & (x<=10)
ax3.plot(x[mask],np.abs(Fw_tmp[mask]))
ax3.set_title(r'$F(\omega)$ expansion')
ax3.grid()
# 特定周波数のみ抽出
Fw_Filter=Fw_tmp;
Hz = 4;
Low = Hz-2.1;
High = Hz+2.1;
Fw_Filter[ ((0 <= np.abs(x)) & (np.abs(x) < Low)) | (High < np.abs(x))]=0;
ax4 = fig.add_subplot(5, 1, 4)
mask = (-10<=x) & (x<=10);
ax4.plot(x[mask],np.abs(Fw_Filter[mask]))
ax4.set_title(r'$F(\omega)$ Filter')
ax4.grid()
# IFFT前にローテーション
fx=np.fft.ifft(np.roll(Fw_Filter,int(L)))
max_fx=max(fx);
max_ft=max(ft);
ax5 = fig.add_subplot(5, 1, 5)
ax5.plot(fx.real)
max_fx=np.max(fx.real)
max_ft=np.max(ft)
ax5.plot(ft*(max_fx/max_ft))
ax5.set_title('f(x)')
ax5.grid()
# 極大値、極小値特定
maxima = np.zeros(N)
minima = np.zeros(N)
fxTmp=fx[0]
mode = 0
for i in range(1,N):
if mode == 0:
if fx[i]>fxTmp:
minima[i]=ft[i]*(max_fx/max_ft)
mode = 1
fxTmp=fx[i]
if mode == 1:
if fx[i]<fxTmp:
maxima[i]=ft[i]*(max_fx/max_ft)
mode = 0
fxTmp=fx[i];
maxima[maxima==0]=np.nan # 極大値だけを残す
minima[minima==0]=np.nan # 極小値だけを残す
ax5.plot(maxima,'ro') # 極大値をplot
ax5.plot(minima,'bo') # 極大値をplot
maxima_index=np.where(~np.isnan(maxima))
minima_index=np.where(~np.isnan(minima))
print('maxima_index:', end='')
print(maxima_index)
print('minima_index:', end='')
print(minima_index)
VTImaxima=VTI[maxima_index] # 極大値要素のみ抽出
VTIminima=VTI[minima_index] # 極大値要素のみ抽出
print('VTI maxima value:', end='')
print(np.round(VTImaxima*127)) # $1=\127
print('VTI minima value:', end='')
print(np.round(VTIminima*127)) # $1=\127
plt.show()コード解説
太郎くん
流れとしてはMATLABの時と全く一緒だね。
太郎くん
あと、論理インデックス検索のところもMATLABの時と同じように修正したんだね。
フクさん
うん。
全く同じ理屈が通用したので修正しておいた。
Fw_Filter[ ( x < -High) | ((-Low < x) & (x <= 0)) | ((0 <= x) & (x < Low)) | (High < x)]=0;↓
Fw_Filter[ ((0 <= np.abs(x)) & (np.abs(x) < Low)) | (High < np.abs(x))]=0;太郎くん
こっちの方が楽そうだよね。
フクさん
これが出来るように複素共役を負の周波数側に持ってきたわけだしね。
太郎くん
まぁ、それを忘れて慌てて修正したわけなんだけど。
フクさん
(ぐうの音も出ない・・・。)
フクさん
よし！
次回は、Python(Numpy)版コードの実働確認だな。
太郎くん
(ごまかしたな。)
まとめ
フクさん
まとめだよ。
- 売却、買付タイミング時のVTI単価特定コードのPython版を作成。
- MATLABの時に実施した論理インデックス検索のやり方を修正。
- このために複素共役を負の周波数側に持ってきた。(忘れてたけど)
バックナンバーはこちら。
コメント