金子邦彦研究室プログラミングOctave の活用Octave での平均、分散の例

Octave でのグラフ作成

Octave での平均、分散の例

前準備

必要となるソフトウェア

CSVファイルを読み込み,テーブルに格納

  1. (前準備) 使用する CSV ファイルの作成

    ※ ここでは Book1.csv をダウンロードし,分かりやすいディレクトリに置く

    (参考: 「外国為替データ(時系列データ)の情報源の紹介」の Web ページ

    以下の説明では、

    として説明を続ける.

  2. 使用する CSV ファイルの確認

    属性名として seq, date, USD, EUR, AUD の5つが, CSV ファイルの1行目に書かれていることを確認する.

    [image]
  3. (オプション) 1列目に通し番号を追加したい場合の操作

    ◆ コマンドの例

    echo '"seq", "date", "USD", "EUR", "AUD"' > Book1.csv  
    awk '{print NR, ",", $0}' hoge.csv >> Book1.csv 
    
  4. octave の起動

    octave
    
  5. octave におけるパッケージの読み込み

    次のコマンドを実行.

    pkg load io
    

    [image]
  6. read.table」を用いて,CSV ファイルを R のデータフレームに読み込み

    Octave で次のコマンドを実行

    Windows の場合のOctaveコマンドの例

    A = dlmread( "C:\\octave\\Book1.csv", ",/", 1, 0 );
    

    [image]

    ◆ Linux の場合のOctaveコマンドの例

    A = dlmread( "/tmp/Book1.csv", ",/", 1, 0 );
    

    [image]
  7. 変数 A の中身の確認

    例えば,次のコマンドを実行する.これは,Aの2行目と3行目と4行目と5行目と6行目と7行目である.

    A(:,[2 3 4 5 6 7])

    [image]

reshape の例

例えば、A の 5 列目を表示してみる

A(:,5)

[image]

◆ A の 5 列目の reshape の例

n 要素ごとに平均 (mean) を求める

n = 7
m = fix( size(A)(1) / n )
S = mean( reshape( A([1:(n * m)],5), n, size(A)(1) / n ) )
plot(S)

[image]

先ほどの結果では、長さが m になる. 長さを n * m に伸ばす.

n = 7
m = fix( size(A)(1) / n )
S = mean( reshape( A([1:(n * m)],5), n, size(A)(1) / n ) )
T = S;
for i = [2:n]
  T = [T; S];  
end
plot( reshape(T, n * m, 1) )

[image]

関数化する

function retval = vecmean(A, x, n)
  m = fix( size(A)(1) / n );
  S = mean( reshape( A([1:(n * m)],x), n, size(A)(1) / n ) ); 
  T = S;
  for i = [2:n]
    T = [T; S];  
  end
  retval = reshape(T, n * m, 1); 
endfunction

hold on
plot(  vecmean(A, 5, 2)  );
plot(  vecmean(A, 5, 8)  );
plot(  vecmean(A, 5, 30)  );
plot(  vecmean(A, 5, 60)  );

[image]

平滑化の例


■ 平均

x = [zeros(1,20), ones(1,10), zeros(1,20)]
pkg load signal
plot([1:50], x, [1:50], filtfilt([1 1 1 1 1]/5, 1, x))

[image]

■ Savitsky-Golay 平滑化 (smoothing) フィルタ

x = [zeros(1,20), ones(1,10), zeros(1,20)]
pkg load signal
plot([1:50], x, [1:50], sgoleyfilt(x), [1:50], filtfilt([1 1 1 1 1]/5, 1, x))

[image]

■ Spencer's 15 point moving average

x = [zeros(1,20), ones(1,10), zeros(1,20)]
pkg load signal
plot([1:50], x, [1:50], spencer(x), [1:50], filtfilt([1 1 1 1 1]/5, 1, x))

[image]

「外国為替データ(時系列データ)」の平滑化(例) 「filtfilt ...の部分は,周辺 20 点の平均.

pkg load signal
plot(A(:,1), A(:,5), A(:,1), filtfilt(ones(1,20)/20, 1, A(:,5)))

[image]

タイトル,ラベル,範囲指定,格子

次の関数を使う

凡例

plot( X[,c("seq")], X[,c("USD")] )
lines(lowess( X[,c("seq")], X[,c("USD")] ), col = "red")
lines(lowess( X[,c("seq")], X[,c("USD")], f=0.2 ), col = "green")
legend( 5, 105, c( "f = 2/3", "f = 0.2" ), lty = 1, col = c("red", "green") )

テキストを書くこともできる.

text(locator(1), labels = "ほげほげ")

ファイルへの保存

グラフィカル・デバイスの種類(例)

など.

例えば,PDF ファイルを作りたいときは,最初に「pdf()」を実行.最後に「dev.off()」で,ファイルを閉じる. すると,作業ディレクトリにファイルができる.作業ディレクトリは,getwd() で分かる.

※ pdf( file="C:\1.png" ) のように,ファイル名を陽に指定する方が分かりやすいでしょう.

[image]