Memorandum: 4月 2016

オープンソースである"R"をMac OS X (Yosemite)にインストールした際のメモです．

CRANにアクセスし，トップページのリンクから"Download R for (Mac) OS X"を選択してクリックします．

OS Xであればpkg形式，Windowsであればexe形式のファイルがダウンロードできます．
ダウンロード後は他のソフトウェアと同じようにインストールを実行します．

pkgファイルをダブルクリックすると，以下のようにインストーラーが立ち上がります．特に設定を変更したい場合を除いて，"Continue" をクリックし続ければOKです．

Rを起動すると，以下のようなコンソール画面が立ち上がります．

gnuplotの機能の一つに，関数のフィッティングがあります．
パラメータを含んだ関数を用意すると，データファイルの数値に対して最適なパラメータの数値を自動的に探索します．
関数は，gnuplotで表現できるものなら，線形・非線形を問いません．

具体例を示します．
以下のようなデータ（圧密試験結果，data.csv）があったとします．

# Consolidation Pressure p(kN/m2),Void Ratio e
0,0.593
9.8,0.584
19.6,0.58
39.2,0.573
78.5,0.564
157,0.551
314,0.536
628,0.515
1256,0.484

このデータに関して，関数
f(x)=c*(x-a)**b
をフィッティングして，変数a, b, cの値を求めます．

まずは，以下のようなスクリプトを実行してデータをプロットしてみます．

set logscale x 10
set xrange [1:10000]
set tics font "Times New Roman,25"

set yrange [0.35:0.65]
set ytics 0.05
set mytics 5
set xlabel "{/=30 Consolidation stress p(kN/m^{2})}"
set ylabel "{/=30 Void Ratio e}"
set size square
set key font "Times New Roman,25"
set key right top

plot "./data.csv" using 1:2 w p ps 1.5 title "Sample 1"

次に変数a,b,cの初期値を求めます．
gnuplotのフィッティングでは，初期値を与えないと全て"1"として計算を始めます．
しかし，反復回数が多くなったり，結果が得られないこともあるようなので，当たりをつけて初期値を決めた方が良いようです（以下の例では初期値を与えない場合を示します）．

最初に，フィッテイングする関数を定義します．
探索するパラメータの変数名は，xやt等のgnuplotが使う特別な変数でなければ何でも大丈夫です．
関数フィッティングには，fitコマンドを使います．
オプションに，探索したいパラメータ名を via に続けて与えます．

f(x)=c*(x-a)**b
fit f(x) "./data.csv" using 1:2 via a,b,c

フィットを行うデータファイルに対しては，using, index, every等のオプションが使うことができ，データファイル中の任意の数値をフィッティングの対象にできます．
また，データの誤差を与える事で，各点での重みを変えることができます．
データの重みは，誤差σの2乗の逆数(σ^{-2})となります．
重み付き最小自乗法を行うためには，データファイルの3カラム目に誤差を用意し， fit f(x) "exp.dat" using 1:2:3 via a,b のようにusing を使って誤差のカラムを指定します．誤差を与えない場合は，全データの重みは1になります．

以下のスクリプト（test.gp）を実行すると，下図が得られます．

set datafile separator ","

set logscale x 10
set xrange [1:10000]
set tics font "Times New Roman,25"

set yrange [0.35:0.65]
set ytics 0.05
set mytics 5
set xlabel "{/=30 Consolidation stress p(kN/m^{2})}"
set ylabel "{/=30 Void Ratio e}"
set size square
set key font "Times New Roman,25"
set key right top

plot "./Dc90-2.csv" using 1:2 w p ps 1.5 title "Sample 1"

f(x)=c*(x-a)**b
fit f(x) "./Dc90-2.csv" using 1:2 via a,b,c

replot f(x), "./Dc90-2.csv" using 1:2 w yerrorbars

ちなみに，Mac OS Xのターミナルから，上記のスクリプトを実行した場合，ターミナルには以下のように表示され，イタレーションの回数（今回の例では15回）や，フィッテングされた変数の値（a, b, c )などがわかります．

gnuplot> load 'test.gp'

iter      chisq       delta/lim  lambda   a             b             c            

   0 2.0958191547e+06   0.00e+00  2.07e+03    1.000000e+00   1.000000e+00   1.000000e+00

   1 3.1401236259e+05  -5.67e+05  2.07e+02    1.000059e+00   8.652118e-01   9.785288e-01

   2 4.3661665492e+04  -6.19e+05  2.07e+01    1.000868e+00   7.302030e-01   9.158026e-01

   3 3.9742713760e+03  -9.99e+05  2.07e+00    1.008882e+00   6.584271e-01   4.525282e-01

   4 7.1730628717e+01  -5.44e+06  2.07e-01    9.923941e-01   6.351242e-01   7.812126e-02

   5 8.6890175128e+00  -7.26e+05  2.07e-02    6.980505e-01   4.985268e-01   7.998705e-02

   6 5.2167471432e-01  -1.57e+06  2.07e-03   -7.175606e+00   2.711390e-01   1.407226e-01

   7 3.0807202872e-01  -6.93e+04  2.07e-04   -2.955954e+01   2.491715e-02   3.138176e-01

   8 1.1567525775e-01  -1.66e+05  2.07e-05   -1.338551e+01  -9.192671e-02   6.734814e-01

   9 5.8730090513e-03  -1.87e+06  2.07e-06   -3.335614e+01  -4.813472e-02   7.368142e-01

  10 9.9947209776e-05  -5.78e+06  2.07e-07   -8.772884e+01  -6.842443e-02   7.982868e-01

  11 8.0106216817e-05  -2.48e+04  2.07e-08   -1.408457e+02  -8.119137e-02   8.727729e-01

  12 3.6356267939e-05  -1.20e+05  2.07e-09   -1.692429e+02  -8.667967e-02   9.118675e-01

  13 3.3880317042e-05  -7.31e+03  2.07e-10   -1.775257e+02  -8.825277e-02   9.234001e-01

  14 3.3854660121e-05  -7.58e+01  2.07e-11   -1.791170e+02  -8.857467e-02   9.256578e-01

  15 3.3854399746e-05  -7.69e-01  2.07e-12   -1.793814e+02  -8.862971e-02   9.260331e-01

iter      chisq       delta/lim  lambda   a             b             c            

After 15 iterations the fit converged.

final sum of squares of residuals : 3.38544e-05

rel. change during last iteration : -7.69104e-06

degrees of freedom    (FIT_NDF)                        : 5

rms of residuals      (FIT_STDFIT) = sqrt(WSSR/ndf)    : 0.00260209

variance of residuals (reduced chisquare) = WSSR/ndf   : 6.77088e-06

Final set of parameters            Asymptotic Standard Error

=======================            ==========================

a               = -179.381         +/- 47.19        (26.3%)

b               = -0.0886297       +/- 0.01025      (11.57%)

c               = 0.926033         +/- 0.06838      (7.384%)

correlation matrix of the fit parameters:

                a      b      c      

a               1.000 

b               0.968  1.000 

c              -0.979 -0.999  1.000 

gnuplot> 

Memorandum

2016年4月20日水曜日

[TeX] 目次をsubsubsectionまで記載する方法

2016年4月7日木曜日

[Mac] RをMac OS X(Yosemite)にインストール

2016年4月5日火曜日

[gnuplot] 関数のフィッテイング