第4回(最終回) 関数のfitting (1/2)


1. gnuplotを使ったfitting

何かの実験データについての考察を行う上で、得られたデータを、理論的(あるいは現象論的)に予想される関数でフィットすることは非常に重要です。gnuplotでは、このfittingを比較的簡単に行うことができ、さらに誤差の見積もりも自動で行ってくれます。今回の授業では、この方法を紹介します。 (作業を始める前に、ターミナルを使って、Desktopに作成したcompphysのディレクトリに移りましょう。)

a) データファイルの用意

まず、フィットするためには、フィットを行う実験データのファイルが必要になります。新しく作るのは面倒なので、前回の授業の最後で使った「data.txt」を使いましょう。まず、data.txtの内容を確認してみましょう。

1.0 12.25 0.15
1.2 14.34 0,18
1.4 16.21 0.26
1.6 18.35 0.34
1.8 21.01 0.58
2.0 22.56 0.72

このような内容のファイルがありましたか?(もしなければemacsで作りましょう。)前回の授業ではgnuplotを使って、このデータからグラフを作成しました。gsを使ってそのepsファイル(横軸T, 縦軸Eとしたもの)も確認してみましょう。

1+1

このような図が確認できたでしょうか?このグラフを見ると、横軸の変数Tと縦軸の変数Eは、E= aT+b(a,bは定数)の相関関係があるように推測できます。そこで以下では、このデータをE= aT+bの線形関数でフィットすることを考えてみましょう。


b) fittingの実行

フィットをするためにはまずgnuplotを立ち上げて、以下のように線形関数f(x)=ax+bを定義します。

1+1

次に、フィットするために以下のコマンドを入力しましょう。

1+1

コマンドの意味はほぼ明らかだと思いますが、一つだけ説明を加えておきます。「using 1:2:3」という部分ですが、これは3列目にある誤差も考慮したフィットであることを示しています。詳細は省きますが、gnuplotではデータに誤差がある場合、誤差の大きいものはあまり信頼せず誤差の小さいものを重要視するというように、誤差をウエイトとしてフィットすることができます。using 1:2:3の部分はこのようなフィットを行うことを意味しています。(逆に単にusing 1:2とすると、誤差を考慮しないフィットが行われます。)

ともあれ、上のコマンドを書いて[enter]を押すと、f(x)の未定係数であったaとbがデータに最もよくフィットするように決定されます。[enter]を押して以下のような画面が表示されれば成功です。

1+1

この結果の見方を以下で説明します。


c) フィット結果の見方

フィットの結果がターミナル上(に立ち上げたgnuplot上)で表示されていると思いますが、ここではいったんquitでgnuplotを終了しましょう。先ほど表示されたフィット結果は、logファイルを見ることでいつでも参照することができます。gnuplotを終了したら、lsでカレントディレクトリにあるファイルを確認してみましょう。

1+1

ここに、fit.logというファイルが増えていますね。これが先ほどgnuplotで行ったfitの結果です。cat等のコマンドで、内容を確認してみましょう。

1+1

以下のように、先ほどgnuplot上で見た結果が再度表示されたでしょうか?

1+1

さて結果の見方ですが、最も重要な行は二重線が引いてある「Finial set of parameters」と「Asymptotic Standard Error」です。これらはそれぞれ、フィットで得られた係数a,bの値と、その誤差を表します。また、その上に表示されているreduced chisquareはフィットの良しあしを測るものなのですが、この意味は今回のレポートで調べてもらうことにします。


d) 結果のプロット

さて、フィットが出来たらその結果を図に書いてみて、フィットがいいかどうか目で確認してみましょう。先ほど得られたaとbの値を使ってf(x)をフィットしてもいいのですが、もっと楽な方法はやはりスクリプトを書くことです。以下のようなスクリプトファイル「gnu3.txt」をemacsで作成してみましょう。このファイルは前回作成した「gnu2.txt」とほぼ同じなので、それをコピーしてから作ると楽に作れると思います。

1+1

このファイルでは、gnu2.txtにおけるデータのプロットに加えて、関数のフィットまで同時に行っています。さらに最後の行では、フィットで得られた関数f(x)をreplotコマンドによって書き足しています。ここで一つ注意点ですが、「set output」の行が二つあることに注目してください。スクリプトによって関数をreplotする場合、初めに宣言するoutputファイルは単なるダミーであって、一番最後にプロットする関数の手前に真のoutputファイルを宣言する必要があります。今の場合、garbage.epsは文字通りダミーファイルで、本当のoutputファイルは「test4.eps」になっています。ともかくこれをgnuplot上で実行すれば、test4.epsに以下のような図が出力されるはずです。

1+1

これを見ると、確かにきちんとフィットされていることが分かりますね。




2. シェルスクリプトでgnuplot

gnuplotを実行するには、まずターミナルでgnuplotと打って、さらにそれからいくつかコマンドを打って、最後はquitするわけですが、これだけでも少なくとも3回コマンドを打たなければならないので大変です。そこで、前回紹介したシェルスクリプトを使ってこれを一度に実行する方法を紹介します。

まず以下の内容のシェルスクリプトをemacsで作成しましょう。ここで***.txtの部分には、実行したいgnuplotのスクリプトファイルを入れましょう。(例えば上で作成したgnu3.txt等)

1+1

後はこれを実行すると、シェルスクリプト中でgnuplotが立ち上がって、***.txtの内容を実行してくれるわけです。その際、chmodで実行権限を与えるのを忘れないようにしましょう。

このようにシェルスクリプト上でgnuplotを立ち上げることは、大規模な数値解析でしばしば有用になります。上のシェルスクリプトはとてもシンプルで、そんなに役に立ちそうにありませんが、このシェルスクリプトの前後に数行付け足すことで、さらに複雑な解析を実行させることができます。例えば、[1]データのフィット⇒[2]フィットで得られた係数を別ファイルに書き出す⇒[3]その結果を使って何か別の量を計算⇒[4]最後にそれらをプロットするといった一連の操作を、一つのシェルスクリプトによって行うことができます。さらに、一度シェルスクリプトを作成しておけば、もとの実験データがアップデートされても、即座に同じ解析を繰り返すことができます。


Home Top Back Next