ttdeamのブログ -4ページ目

ツインテールdeエンジェルモード!! で Deep Learning!! （活性化関数3）

ツインテールdeエンジェルモード!! で活性化関数のうち、

「ReLU関数」を記述してみます。

この関数も、数値だけじゃなくて配列も処理できるように記述してみます。

こんな感じでうまく行きそうです。

def relu(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0 ? x[k]:0
   }
   else
           ans    = x   >0 ? x   :0
   retn(ans)
}

数値で動作を確認してみます。

% cat relu1
def relu(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0 ? x[k]:0
   }
   else
           ans    = x   >0 ? x   :0
   retn(ans)
}

for(i=-3;i<3;i++)
print("%d %f\n",i,relu(i))

% tt relu1
-3 0.000000
-2 0.000000
-1 0.000000
0 0.000000
1 1.000000
2 2.000000

配列で動作を確認してみます。

% cat relu2
def relu(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0 ? x[k]:0
   }
   else
           ans    = x   >0 ? x   :0
   retn(ans)
}

tmp={ -1, +1, +2 }
ans=relu(tmp)
p(tmp,ans)

% tt relu2
tmp[0]=-1,tmp[1]=1,tmp[2]=2
ans[0]=0,ans[1]=1,ans[2]=2

次に、グラフ化してみます。

% cat plot-relu
def relu(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0 ? x[k]:0
   }
   else
           ans    = x   >0 ? x   :0
   retn(ans)
}

for( x=-5 ; x<5 ; x+=0.1 )
print("%f %f\n",x,relu(x))

% tt plot-relu > data

% gnuplot

G N U P L O T
Version 5.0 patchlevel 5 last modified 2016-10-02

Copyright (C) 1986-1993, 1998, 2004, 2007-2016
Thomas Williams, Colin Kelley and many others

   gnuplot home:     http://www.gnuplot.info
   faq, bugs, etc:   type "help FAQ"
   immediate help:   type "help" (plot window: hit 'h')

Terminal type set to 'x11'
gnuplot> plot [][-0.1:+5.1]"data" with line

最後に、3つの活性化関数を同時にグラフ化してみます。

スクリプトはこんな感じでOKです。

（ :r とあるのは、ソーススクリプトを読み込む命令です。）

% cat plot-3
:r prog/step
:r prog/sigmoid
:r prog/relu

for( x=-5 ; x<5 ; x+=0.1 )
print("%f %f %f %f\n",x,step(x),sigmoid(x),relu(x))

% tt plot-3 > data

% gnuplot

G N U P L O T
Version 5.0 patchlevel 5 last modified 2016-10-02

Copyright (C) 1986-1993, 1998, 2004, 2007-2016
Thomas Williams, Colin Kelley and many others

   gnuplot home:     http://www.gnuplot.info
   faq, bugs, etc:   type "help FAQ"
   immediate help:   type "help" (plot window: hit 'h')

Terminal type set to 'x11'
gnuplot> set grid
gnuplot> plot [][-0.1:3]"data" u 1:1 t "Step " w l, "" u 1:2 t "Sigmoid" w l, "" u 1:3 t "ReLU " w l

うまく表示できたようです。

ツインテールdeエンジェルモード!! で Deep Learning!! （活性化関数2）

ツインテールdeエンジェルモード!! で活性化関数のうち、

「シグモイド関数」を記述してみます。

せっかくなので、数値だけじゃなくて配列も処理できるように記述してみます。

こんな感じでうまく行きそうです。

def sigmoid(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = 1.0/(1.0+exp(-x[k]))
   }
   else
           ans    = 1.0/(1.0+exp(-x   ))
   retn(ans)
}

数値で動作を確認してみます。

% cat sigmoid1
def sigmoid(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = 1.0/(1.0+exp(-x[k]))
   }
   else
           ans    = 1.0/(1.0+exp(-x   ))
   retn(ans)
}

for(i=-3;i<3;i++)
print("%d %f\n",i,sigmoid(i))

% tt sigmoid1
-3 0.047426
-2 0.119203
-1 0.268941
0 0.500000
1 0.731059
2 0.880797

配列で動作を確認してみます。

% cat sigmoid2
def sigmoid(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = 1.0/(1.0+exp(-x[k]))
   }
   else
           ans    = 1.0/(1.0+exp(-x   ))
   retn(ans)
}

tmp={ -1, +1, +2 }
ans=sigmoid(tmp)
p(tmp,ans)

% tt sigmoid2
tmp[0]=-1,tmp[1]=1,tmp[2]=2
ans[0]=0.268941,ans[1]=0.731059,ans[2]=0.880797

次に、グラフ化してみます。

% cat plot-sigmoid
def sigmoid(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = 1.0/(1.0+exp(-x[k]))
   }
   else
           ans    = 1.0/(1.0+exp(-x   ))
   retn(ans)
}

for( x=-5 ; x<5 ; x+=0.1 )
print("%f %f\n",x,sigmoid(x))

% tt plot-sigmoid > data

% gnuplot

G N U P L O T
Version 5.0 patchlevel 5 last modified 2016-10-02

Copyright (C) 1986-1993, 1998, 2004, 2007-2016
Thomas Williams, Colin Kelley and many others

   gnuplot home:     http://www.gnuplot.info
   faq, bugs, etc:   type "help FAQ"
   immediate help:   type "help" (plot window: hit 'h')

Terminal type set to 'x11'
gnuplot> plot [][-0.1:+1.1]"data" with line

うまく表示できました。

ツインテールdeエンジェルモード!! で Deep Learning!! （活性化関数1）

ツインテールdeエンジェルモード!! で活性化関数のうち、一番簡単な

「ステップ関数」を記述してみます。

ステップ関数自体はとても簡単なのですが、せっかくなので、

数値だけじゃなくて配列も処理できるように関数を記述してみます。

こんな感じでうまく行きそうです。

def step(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0
   }
   else
           ans    = x   >0
   retn(ans)
}

動作を確認してみます。まずは、数値から、

% cat step1
def step(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0
   }
   else
           ans    = x   >0
   retn(ans)
}

for(i=-3;i<3;i++)
print("%d %d\n",i,step(i))

% tt step1
-3 0
-2 0
-1 0
0 0
1 1
2 1

特に、問題なく動いています。

次に、配列です。

% cat step2
def step(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0
   }
   else
           ans    = x   >0
   retn(ans)
}

tmp={-30,-20,-10,0,10,20,30}
ans=step(tmp)
p(tmp,ans)

% tt step2
tmp[0]=-30,tmp[1]=-20,tmp[2]=-10,tmp[3]=0,tmp[4]=10,tmp[5]=20,tmp[6]=30
ans[0]=0,ans[1]=0,ans[2]=0,ans[3]=0,ans[4]=1,ans[5]=1,ans[6]=1

こちらも、上手く動作してるみたいですね。

次に、グラフ化してみます。（ x=-5 〜 x=+5 までの範囲）

例によって、餅は餅屋でデータ計算をスクリプトで実行して、

グラフ表示は gnuplot を利用してみます。

% cat plot-step
def step(x){
   if( type(x)=='A' ){
       each( k=keys(x) )
           ans[k] = x[k]>0
   }
   else
           ans    = x   >0
   retn(ans)
}

for( x=-5 ; x<5 ; x+=0.1 )
print("%f %d\n",x,step(x))

% tt plot-step > data

% gnuplot

G N U P L O T
Version 5.0 patchlevel 5 last modified 2016-10-02

Copyright (C) 1986-1993, 1998, 2004, 2007-2016
Thomas Williams, Colin Kelley and many others

   gnuplot home:     http://www.gnuplot.info
   faq, bugs, etc:   type "help FAQ"
   immediate help:   type "help" (plot window: hit 'h')

Terminal type set to 'x11'
gnuplot> plot [][-0.1:+1.1]"data" with line

綺麗に、表示できました。

MNISTデータの画像化（その2）

追加。。。

一旦、コマンドスクリプト mk-mnistbmp が出来てしまえば、

全データを画像化するなんてことは、シェルを使えば簡単にできてしまいます。

これも、餅は餅屋ですね！

得意な分野は、得意なソフトにまかせましょう！

全部1つのソフトでやろうなんて、Windowsじゃあるまいし、

UNIXでは、ソフトウェアレバレッジを効かせてこそ真価が発揮されるのです！

では、そのやり方ですが、、、

（1）テストイメージの場合。

% for i in `seq 0 999`; do mk-mnistbmp test.image $i > L$i.bmp; done

の1行で、OKです。

（2）学習用イメージの場合

% for i in `seq 0 5999`; do mk-mnistbmp test.image $i > L$i.bmp; done

の1行で、OKです。

どちらも、大量の画像ファイルが生成されますのでご注意を！（笑）

まぁ、消すのも次の1コマンドですみますけどね！

% rm *.bmp

では、

MNISTデータの画像化

ツインテールdeエンジェルモード!! で読み込んだMNISTデータを画像化するスクリプトを書いてみます。

画像ファイルの形式ですが、とりあえずはシンプルなBMP形式としてみます。

なお、BMP形式のフォーマットは、ググれば沢山でてきますので、その説明はスキップして、まずは、24ビットのBMP形式画像ファイルのヘッダ（54バイト）を出力するスクリプトを書いてみます。

X-Pixel値とY-Pixel値に依存しますので、こんな感じに書けばいいと思います。

% cat mk-bmphead
#!/usr/local/bin/tt

# 説明：２４ビットＢＭＰファイルヘッダーを出力します。
# この出力に、ＲＧＢ画像データを追加出力すれば、ＢＭＰファイルが完成します。
# なお、追加するＲＧＢ画像データのフォーマットは次の通りです。
#（１）ｙ＝０行（最下行）〜ｙ＝最終行（最上行）までの画素データについて、、、
#（２）    ｘ＝０点（最左点）〜ｘ＝最終点（最右点）までの画素データについて、、、
#（３）        ｛青＋緑＋赤｝の順に各１バイトの輝度データを出力する。
#（４）    各行の終わりに、パッドデータ（0x00）を出力する。
#（５）終了。

# パラメーターの取得
   if( argc!=3 )
       dying("Usage: ${cmd} X-Pixel Y-Pixel\n")
   xpic = int(shift())
   ypic = int(shift())

# 画像データ領域＆BMPファイル全体のサイズ計算
rem = xpic*3 % 4 # 0, 1, 2, 3
pad = (rem==0 ? 0:(4-rem)) # 0, 3, 2, 1

i_size = (xpic*3+pad)*ypic # [byte] 各行が４の倍数のバイト数を持つ。
f_size = i_size+54 # [byte] BMPファイル＝ヘッダ＋画像データ

# ファイルヘッダの出力
   putc('B'); putc('M')
   loop(i<4){ putc((f_size>>i*8)&0xFF) }
   loop(i<4){ putc(0x00) }
   loop(i<4){ putc((54    >>i*8)&0xFF) }

# 情報ヘッダの出力
   loop(i<4){ putc((40    >>i*8)&0xFF) }
   loop(i<4){ putc((xpic >>i*8)&0xFF) }
   loop(i<4){ putc((ypic >>i*8)&0xFF) }
   loop(i<2){ putc(( 1    >>i*8)&0xFF) }
   loop(i<2){ putc((24    >>i*8)&0xFF) }
   loop(i<4){ putc(( 0    >>i*8)&0xFF) }
   loop(i<4){ putc((i_size>>i*8)&0xFF) }
   loop(i<4){ putc(( 1    >>i*8)&0xFF) }
   loop(i<4){ putc(( 1    >>i*8)&0xFF) }
   loop(i<4){ putc(( 0    >>i*8)&0xFF) }
   loop(i<4){ putc(( 0    >>i*8)&0xFF) }

# インフォメーションの表示
   eprint("\n")
   eprint("File=%d[Byte] (Head=%d/Data=%d)\n",f_size,54,i_size)
   eprint("> X=%d*3+%d\n",xpic,pad)
   eprint("> Y=%d\n",ypic)
   eprint("\n")

MNISTデータの手書きイメージファイルは、横28X縦28なので、

使い方は、こんな感じとなります。（実行属性をつけておいて下さい。）

% mk-bmphead 28 28 > ファイル名

次に、このスクリプトを利用して、MNISTイメージデータのN番目（0が最初）

の画像データをBMP化するスクリプトを書いてみます。

ただし、BMPデータは、先頭画素が左下の画素から始まるのと、

1画素あたりRGBの3要素が必要なのを注意すれば、

簡単に、こんな感じに書けると思います。

% cat mk-mnistbmp
#!/usr/local/bin/tt

# 説明：MNIST Image ファイルを読み込む。戻り値はイメージ値配列。(28*28)
def mnist_rimage(f){
   a=NULL                               # Error Value !!
   loop(i<4*4){ getc(f) }               # Skip Header
   for( i=0 ; TRUE ; i++ ){           # All Image(s)
       (p,c)=read(f,28*28)
       if( c==0 || c!=28*28 )
           retn(a)
       a[i]=p
   }
}

# パラメーターの取得
   if( argc!=3 )
       dying("Usage: ${cmd} (Image)FileName Index\n")
   fin=    shift()
   idx=int(shift())

# データの読み込み
a_img=mnist_rimage(fin)

# 画像データの作成
   system("mk-bmphead 28 28")
   for( y=28-1 ; y>=0 ; y-- ){
       for( x=0 ; x<28 ; x++ ){
           i=y*28+x
           putc(a_img[idx][i])   # B #
           putc(a_img[idx][i])   # G #
           putc(a_img[idx][i])   # R #
       }
   }

使い方は、こんな感じになります。

例えば、先頭（0番目）の画像データをBMPファイル化するには、、、

% mk-mnistbmp test.image 0 > x.bmp

File=2406[Byte] (Head=54/Data=2352)
> X=28*3+0
> Y=28

といった感じとなります。

なお、BMPファイル自体の表示は、ブラウザ等の他の専門家にまかせましょう。

餅は餅屋ですから。

ちなみに、テストデータの最初の手書き文字を表示したら、

数字の「7」のようでした。

とりあえず、ここまで。

MNISTデータの読み込み

ツインテールdeエンジェルモード!! でMNISTデータを読み込むスクリプトを

書いてみます。

書式自体はシンプルなので、こんな感じに書けばいいと思います。

（1）イメージファイルの読み込み関数

# 説明：MNIST Image ファイルを読み込む。戻り値はイメージ値配列。(28*28)
def mnist_rimage(f){
   a=NULL                               # Error Value !!
   loop(i<4*4){ getc(f) }               # Skip Header
   for( i=0 ; TRUE ; i++ ){           # All Image(s)
       (p,c)=read(f,28*28)
       if( c==0 || c!=28*28 )
           retn(a)
       a[i]=p
   }
}

（2）ラベルファイルの読み込み関数

# 説明：MNIST Label ファイルを読み込む。戻り値はラベル値配列。
def mnist_rlabel(f){
   a=NULL                               # Error Value !!
   loop(i<4*2){ getc(f) }               # Skip Header
   for( i=0 ; (d=getc(f))!=NULL ; i++ )# All Label(s)
       a[i]=d
   retn(a)
}

（3）2つ合わせると、こんな感じです。

% cat rd-mnist
#!/usr/local/bin/tt

# 説明：MNIST Image ファイルを読み込む。戻り値はイメージ値配列。(28*28)
def mnist_rimage(f){
   a=NULL                               # Error Value !!
   loop(i<4*4){ getc(f) }               # Skip Header
   for( i=0 ; TRUE ; i++ ){           # All Image(s)
       (p,c)=read(f,28*28)
       if( c==0 || c!=28*28 )
           retn(a)
       a[i]=p
   }
}

# 説明：MNIST Label ファイルを読み込む。戻り値はラベル値配列。
def mnist_rlabel(f){
   a=NULL                               # Error Value !!
   loop(i<4*2){ getc(f) }               # Skip Header
   for( i=0 ; (d=getc(f))!=NULL ; i++ )# All Label(s)
       a[i]=d
   retn(a)
}

# パラメーターの取得
   if( argc!=3 )
       dying("Usage: ${cmd} (Image)FileName (Label)FileName\n")
   fi=shift()
   fl=shift()

# データの読み込み
a_img=mnist_rimage(fi)
a_lbl=mnist_rlabel(fl)

このスクリプトファイルに実行属性をつけると、

% rd-mnist test.image test.label

みたいな感じでコマンドラインから直接実行できるようになります。

読み込んだデータはイメージとラベルの配列となります。

具体的には、 a_img[0] が最初のイメージデータの28*28バイト。

その各ピクセルの値にアクセスするには、

a_img[0][n] （ただし、 n=0 から 28*28-1 の整数）

みたいにします。

また、a_lbl[0] が最初のラベルの1バイト（値は7）です。

ツインテールdeエンジェルモード!! で Deep Learning!! （MNISTデータ）

画像認識でつかうMNISTデータセットファイルを取得してみます。

これは、手書きの数字（0〜9）を書いたイメージデータで、

「イメージデータ」と「その答え（ラベルデータ）」があります。

イメージデータは、縦28X横28の白黒画像（1ピクセルあたり、8ビット）

ラベルデータは、イメージデータの答え（0〜9）

となっています。

また、学習用（6万イメージ）とテスト用（1万イメージ）の2種類があります。

全てスクリプトを書いて取得してもいいのですが、せっかくいいコマンドがあるので、

ここでは適材適所で対応します。（UNIXの基本精神ですよね？）

（1）まずは、ファイルの取得

% wget http://yann.lecun.com/exdb/mnist/train-images-idx3-ubyte.gz
% wget http://yann.lecun.com/exdb/mnist/train-labels-idx1-ubyte.gz
% wget http://yann.lecun.com/exdb/mnist/t10k-images-idx3-ubyte.gz
% wget http://yann.lecun.com/exdb/mnist/t10k-labels-idx1-ubyte.gz

取得したファイルのサイズはこんな感じです。

% ls -l t*
-rw-r--r--. 1 root root 1648877 Jul 22 2000 t10k-images-idx3-ubyte.gz
-rw-r--r--. 1 root root 4542 Jul 22 2000 t10k-labels-idx1-ubyte.gz
-rw-r--r--. 1 root root 9912422 Jul 22 2000 train-images-idx3-ubyte.gz
-rw-r--r--. 1 root root 28881 Jul 22 2000 train-labels-idx1-ubyte.gz

（2）次に、ファイルの解凍

% gzip -d train-images-idx3-ubyte.gz
% gzip -d train-labels-idx1-ubyte.gz
% gzip -d t10k-images-idx3-ubyte.gz
% gzip -d t10k-labels-idx1-ubyte.gz

解凍したファイルのサイズはこんな感じです。

% ls -l t*
-rw-r--r--. 1 root root 7840016 Jul 22 2000 t10k-images-idx3-ubyte
-rw-r--r--. 1 root root 10008 Jul 22 2000 t10k-labels-idx1-ubyte
-rw-r--r--. 1 root root 47040016 Jul 22 2000 train-images-idx3-ubyte
-rw-r--r--. 1 root root 60008 Jul 22 2000 train-labels-idx1-ubyte

（3）名前の変更

ファイルの名前がわかりにくいので、変更します。

ls 出力の上から順に、テスト用イメージ、テスト用ラベル、学習用イメージ、

学習用ラベル、となってますので、それぞれ、わかりやすい名前に変更します。

% mv train-images-idx3-ubyte learning.image
% mv train-labels-idx1-ubyte learning.label
% mv t10k-images-idx3-ubyte test.image
% mv t10k-labels-idx1-ubyte test.label

% ls -l l* t*
-rw-r--r--. 1 root root 47040016 Jul 22 2000 learning.image
-rw-r--r--. 1 root root 60008 Jul 22 2000 learning.label
-rw-r--r--. 1 root root 7840016 Jul 22 2000 test.image
-rw-r--r--. 1 root root 10008 Jul 22 2000 test.label

learning.image が、学習用の手書きイメージファイルです。

（16バイトのヘッダ＋6万イメージ入っています。

1イメージにつき、28＊28＝784バイト。）

learning.label が、その答えファイルです。

（8バイトヘッダ＋6万個の答えが入っています。

1イメージについき、1バイト。）

test.image が、テスト用の手書きイメージファイルです。

（16バイトのヘッダ＋1万イメージ入っています。

1イメージにつき、28＊28＝784バイト。）

test.label が、その答えファイルです。

（8バイトヘッダ＋1万個の答えが入っています。

1イメージについき、1バイト。）

ツインテールdeエンジェルモード!! で Deep Learning!! （パーセプトロン）

ここでは、AND／NAND／OR／XORゲートを関数で作ります。

まずは、単純なANDゲートとその実行。

% cat and
def and(x1,x2){
   w1=0.5
   w2=0.5
   theta=0.7
   tmp=x1*w1+x2*w2
   return( tmp<=theta ? 0 : 1 )
}

p(and(0,0))
p(and(0,1))
p(and(1,0))
p(and(1,1))

% tt and
0
0
0
1

次に、改良したANDゲートと実行

% cat and2
def and(x1,x2){
   x={x1,x2}
   w={0.5,0.5}
   b=-0.7
   each( k=keys(x) )
       z[k]=x[k]*w[k]
   tmp=sum(z)+b
   return( tmp<=0 ? 0 : 1 )
}

p(and(0,0))
p(and(0,1))
p(and(1,0))
p(and(1,1))

% tt and2
0
0
0
1

次に、NANDゲート。（パラメータが違うだけです。）

% cat nand
def nand(x1,x2){
   x={x1,x2}
   w={-0.5,-0.5}
   b=+0.7
   each( k=keys(x) )
       z[k]=x[k]*w[k]
   tmp=sum(z)+b
   return( tmp<=0 ? 0 : 1 )
}

p(nand(0,0))
p(nand(0,1))
p(nand(1,0))
p(nand(1,1))

% tt nand
1
1
1
0

次に、ORゲート。（これも、パラメータが違うだけです。）

% cat or
def or(x1,x2){
   x={x1,x2}
   w={0.5,0.5}
   b=-0.2
   each( k=keys(x) )
       z[k]=x[k]*w[k]
   tmp=sum(z)+b
   return( tmp<=0 ? 0 : 1 )
}

p(or(0,0))
p(or(0,1))
p(or(1,0))
p(or(1,1))

% tt or
0
1
1
1

最後に、XORゲート。（いままでの関数を使います。）

% cat xor
def and(x1,x2){
   x={x1,x2}
   w={0.5,0.5}
   b=-0.7
   each( k=keys(x) )
       z[k]=x[k]*w[k]
   tmp=sum(z)+b
   return( tmp<=0 ? 0 : 1 )
}

def nand(x1,x2){
   x={x1,x2}
   w={-0.5,-0.5}
   b=+0.7
   each( k=keys(x) )
       z[k]=x[k]*w[k]
   tmp=sum(z)+b
   return( tmp<=0 ? 0 : 1 )
}

def or(x1,x2){
   x={x1,x2}
   w={0.5,0.5}
   b=-0.2
   each( k=keys(x) )
       z[k]=x[k]*w[k]
   tmp=sum(z)+b
   return( tmp<=0 ? 0 : 1 )
}

def xor(x1,x2){
   s1=nand(x1,x2)
   s2=or (x1,x2)
   y =and (s1,s2)
   return( y )
}

p(xor(0,0))
p(xor(0,1))
p(xor(1,0))
p(xor(1,1))

% tt xor
0
1
1
0

重みとバイアス値が違うだけで、いろんな挙動が実現できました。

とりあえず、ここまで。

ツインテールdeエンジェルモード!! で Deep Learning!! （その4）

1.6 外部ライブラリ（Matplotlib）

ツインテールdeエンジェルモード!! では、インタプリタ本体である /usr/local/bin/tt のみで完結しますので、これ以外の外部ライブラリ等はありません。

よって、グラフ描画は外部のグラフプログラムに描かせます。（餅は餅屋です。）

ここでは、 gnuplot を利用してみます。

http://www.gnuplot.info/

具体的なプロット方法ですが、

（1）ツインテールdeエンジェルモード!! でグラフデータを生成する。

（2）それを、 gnuplot で描画する。

という流れとなります。

では、簡単な sin() 関数のグラフから、、、

まずは、こんな感じのスクリプトファイル vvv を作成します。

% cat vvv
for( x=0 ; x<6 ; x+=0.1 ){
y=sin(x)
print("%f %f\n",x,y)
}

それを普通に実行して、結果をファイル data にリダイレクトして保存します。

% tt vvv > data

結果はこんな感じになります。

% cat data
0.000000   0.000000
0.100000   0.099833
0.200000   0.198669
0.300000   0.295520
0.400000   0.389418
0.500000   0.479426
0.600000   0.564642
0.700000   0.644218
0.800000   0.717356
0.900000   0.783327
1.000000   0.841471
1.100000   0.891207
1.200000   0.932039
1.300000   0.963558
1.400000   0.985450
1.500000   0.997495
1.600000   0.999574
1.700000   0.991665
1.800000   0.973848
1.900000   0.946300
2.000000   0.909297
2.100000   0.863209
2.200000   0.808496
2.300000   0.745705
2.400000   0.675463
2.500000   0.598472
2.600000   0.515501
2.700000   0.427380
2.800000   0.334988
2.900000   0.239249
3.000000   0.141120
3.100000   0.041581
3.200000   -0.058374
3.300000   -0.157746
3.400000   -0.255541
3.500000   -0.350783
3.600000   -0.442520
3.700000   -0.529836
3.800000   -0.611858
3.900000   -0.687766
4.000000   -0.756802
4.100000   -0.818277
4.200000   -0.871576
4.300000   -0.916166
4.400000   -0.951602
4.500000   -0.977530
4.600000   -0.993691
4.700000   -0.999923
4.800000   -0.996165
4.900000   -0.982453
5.000000   -0.958924
5.100000   -0.925815
5.200000   -0.883455
5.300000   -0.832267
5.400000   -0.772764
5.500000   -0.705540
5.600000   -0.631267
5.700000   -0.550686
5.800000   -0.464602
5.900000   -0.373877
6.000000   -0.279415

ここで、1つ注意。

x<6 だった筈なのに x=6.000000 が出力されています。

これは電子計算機で実数計算を行う場合に発生する、「誤差」によるものです。

ツインテールdeエンジェルモード!! は、生のコンピューターの演算結果を

そのまま返しますので、このようになります。

他の言語であっても、注意が必要です。

ちなみに、C言語で計算しても全く同じ結果となります。

これを回避するには、（1）ループ変数と終了判定を整数で行う。

（2）終了判定に余裕をもたせる。などの方法があります。

最初の方法だと、変数 t を導入して、

例えばこんな感じにすることで回避できます。

for( t=0 ; t<60 ; t+=1 ){
   x=t/10.0
   y=sin(x)
   print("%f   %f\n",x,y)
}

あるいは、こんな感じに終了判定を変更することで、

回避することもできます。

for( x=0 ; x<6-0.01 ; x+=0.1 ){
y=sin(x)
print("%f %f\n",x,y)
}

話がそれましたが、データファイル data のグラフをプロットしてみます。

gnuplot を立ち上げて次のコマンドを入力します。

% gnuplot

G N U P L O T
Version 5.0 patchlevel 5 last modified 2016-10-02

Copyright (C) 1986-1993, 1998, 2004, 2007-2016
Thomas Williams, Colin Kelley and many others

   gnuplot home:     http://www.gnuplot.info
   faq, bugs, etc:   type "help FAQ"
   immediate help:   type "help" (plot window: hit 'h')

Terminal type set to 'x11'
gnuplot> plot "data" with line
gnuplot>

同様に、 sin(x) と cos(x) のグラフも描いてみます。

こんな感じになります。

% cat vvv
for( x=0 ; x<6 ; x+=0.1 ){
   y=sin(x)
   z=cos(x)
   print("%f   %f   %f\n",x,y,z)
}

% tt vvv > data

% cat data
0.000000   0.000000   1.000000
0.100000   0.099833   0.995004
0.200000   0.198669   0.980067
0.300000   0.295520   0.955336
0.400000   0.389418   0.921061
0.500000   0.479426   0.877583
0.600000   0.564642   0.825336
0.700000   0.644218   0.764842
0.800000   0.717356   0.696707
0.900000   0.783327   0.621610
1.000000   0.841471   0.540302
1.100000   0.891207   0.453596
1.200000   0.932039   0.362358
1.300000   0.963558   0.267499
1.400000   0.985450   0.169967
1.500000   0.997495   0.070737
1.600000   0.999574   -0.029200
1.700000   0.991665   -0.128844
1.800000   0.973848   -0.227202
1.900000   0.946300   -0.323290
2.000000   0.909297   -0.416147
2.100000   0.863209   -0.504846
2.200000   0.808496   -0.588501
2.300000   0.745705   -0.666276
2.400000   0.675463   -0.737394
2.500000   0.598472   -0.801144
2.600000   0.515501   -0.856889
2.700000   0.427380   -0.904072
2.800000   0.334988   -0.942222
2.900000   0.239249   -0.970958
3.000000   0.141120   -0.989992
3.100000   0.041581   -0.999135
3.200000   -0.058374   -0.998295
3.300000   -0.157746   -0.987480
3.400000   -0.255541   -0.966798
3.500000   -0.350783   -0.936457
3.600000   -0.442520   -0.896758
3.700000   -0.529836   -0.848100
3.800000   -0.611858   -0.790968
3.900000   -0.687766   -0.725932
4.000000   -0.756802   -0.653644
4.100000   -0.818277   -0.574824
4.200000   -0.871576   -0.490261
4.300000   -0.916166   -0.400799
4.400000   -0.951602   -0.307333
4.500000   -0.977530   -0.210796
4.600000   -0.993691   -0.112153
4.700000   -0.999923   -0.012389
4.800000   -0.996165   0.087499
4.900000   -0.982453   0.186512
5.000000   -0.958924   0.283662
5.100000   -0.925815   0.377978
5.200000   -0.883455   0.468517
5.300000   -0.832267   0.554374
5.400000   -0.772764   0.634693
5.500000   -0.705540   0.708670
5.600000   -0.631267   0.775566
5.700000   -0.550686   0.834713
5.800000   -0.464602   0.885520
5.900000   -0.373877   0.927478
6.000000   -0.279415   0.960170

そしてプロットです。

% gnuplot

G N U P L O T
Version 5.0 patchlevel 5 last modified 2016-10-02

Copyright (C) 1986-1993, 1998, 2004, 2007-2016
Thomas Williams, Colin Kelley and many others

   gnuplot home:     http://www.gnuplot.info
   faq, bugs, etc:   type "help FAQ"
   immediate help:   type "help" (plot window: hit 'h')

Terminal type set to 'x11'
gnuplot> plot "data" , "data" using 1:3

using 1:3 とは、1カラム目（x）と3カラム目（z）のデータを用いるということです。

あと、画像の表示ですが、これは普通にブラウザで見た方が良さそうです。

餅は餅屋です。

とりあえず、ここまで。

2次元配列（追加）

a[0,0]=1, a[0,1]=2, a[1,0]=3, a[1,1]=4
p(a)

b[0,0]=3, b[0,1]=0, b[1,0]=0, b[1,1]=6
p(b)

for( i=0 ; i<2 ; i++ ){
   for( j=0 ; j<2 ; j++ ){
       z[i,j] = a[i,j] + b[i,j]
   }
}
p(z)

と書きましたが、さらに、

a = { {1,2} , {3,4} }

b = { {3,0} , {0,6} }

k = keys(a)

each( t=k ){ z[t]=a[t]+b[t] }

と書く方が、よりスマートでした。

なお、 k=keys(a) は、配列 a の添字を集めて配列化したものです。

配列 a が1次元配列の時は、そのままの添字値が集められますが、

配列 a が多次元配列の時は、内部表現での添字値が集められます。

今日は、ここまで。