WANtaroHP (F90 STATISTICS)

ページ概要

Contents
単回帰分析
重回帰分析
主成分分析
標準正規分布における％点および超過確率
確率楕円
クラスター分析 (K-means++ 法)
多次元尺度法

単回帰分析

プログラム概要

単回帰分析を行うプログラムです．
回帰直線 y=aa*x+bb における係数 aa, 係数 bb および相関係数を得ます．

ソースコードと関連文書

Filename	Description
f90_KAIKI.txt	program for single regression analysis

重回帰分析

プログラム概要

重回帰分析を行うプログラムです．
重回帰式の形は y=b₀+b₁*x₁+b₂*x₂+b₃*x₃+....
ここに, y は目的変数, b_i は偏回帰係数， x_i は説明変数です．
連立一次方程式は，Gauss-Jordan の掃き出し法で解かれています．

入力データファイル書式

 001 | Sample data
 002 | 4,46
 003 | 44523,6364,90547,19.7,640.9
 004 | 10395,7135,12634,20.0,128.2
 005 |  8669,7266,20101,16.3,113.8
..........

001: コメント
002: 説明変数の数 (=4), 入力データの組数 (=46)
003: (この場合，y値と3つのx値)
004: (3 行目以降はデータ：y値と3つのx値)

実行用バッチコマンド

gfortran -o f90_MRA.exe f90_MRA.f90
f90_MRA inp_MRA.txt out_MRA.txt

バッチコマンド書式

f90_MRA fnameR fnameW

f90_MRA	コンパイルされた実行ファイル名
fnameR	入力ファイル名
fnameW	出力ファイル名

ソースコードと関連文書

Filename	Description
f90_MRA.txt	プログラムソースコード
inp_MRA.txt	入力データサンプル
out_MRA.txt	出力データサンプル
R_DRAWinp.txt	入力データ作図用Rスクリプト
R_DRAWout.txt	出力データ作図用Rスクリプト
fig_MRA_inp.pdf	入力データ作図事例
fig_MRA_out.pdf	出力データ作図事例

主成分分析

プログラム概要

主成分分析を行うプログラムです．
固有値および固有ベクトルは， Jacobi法にて計算されます．

入力データ書式

 001 | iris
 002 | 4,150
 003 | 5.1,3.5,1.4,0.2,setosa,1
 004 | 4.9,3.0,1.4,0.2,setosa,2
 005 | 4.7,3.2,1.3,0.2,setosa,3
..........

001: コメント
002: データ列数 (=4), 入力データ組数 (=150)
003: (この場合4列の入力データ：5列目および6列目は無効)
004: (3行目以降入力データが続く

実行用バッチコマンド

gfortran -o f90_PCA.exe f90_PCA.f90

f90_PCA 1 1 inp_iris.csv out_iris.csv

rscript R_DRAWinp.R
ren Rplots.pdf fig_R_PCA_inp.pdf

rscript R_DRAWout.R
ren Rplots.pdf fig_R_PCA_out.pdf

実行用バッチコマンド書式

f90_PCA knor kvok fnameR fnameW

f90_PCA	コンパイルされた実行ファイル
knor	正規化の有無 (0: 正規化する, 1: 正規化しない)
kvok	分析に用いる行列 (0: 分散共分散行列使用, 1: 相関行列使用)
fnameR	入力データファイル名
fnameW	出力データファイル名

ソースコードと関連文書

Filename	Description
f90_PCA.txt	f90プログラムソース
inp_iris.txt	入力データ事例
out_iris.txt	出力データ事例
R_DRAWinp.txt	入力データ作図用Rスクリプト
R_DRAWout.txt	出力データ作図用Rスクリプト
fig_PCA_inp.pdf	入力データ作図事例
fig_PCA_out.pdf	出力データ作図事例

標準正規分布における％点および超過確率

プログラム概要

標準正規分布における％点および超過確率を計算するプログラムです．
％点を得るためには，「戸田の近似」が使われています．
超過確率を得るためには，「Shentonの近似」が使われています．

ソースコードと関連文書

Filename	Description
f90_TODAPP.txt	プログラムソースコード
inp_TODAPP.txt	入力データサンプル
out_TODAPP.txt	出力データサンプル

確率楕円

プログラム概要

ここには，2つのプログラム f90_SREG.f90 と f90_ELLIP.f90 が紹介されています．
プログラム f90_SREG.f90 は，単回帰分析と等確率楕円の特性をみるためのものです．
プログラム f90_ELLIP.f90 は，GMTの楕円を描画するオプション 'Se' を使うためのデータを作成するものです． 'Se' はGMTのコマンド 'psxy' のオプションであり，大量のプロットデータを必要としないで楕円を描くことができます．このプログラムでは，既知の3点から楕円の方程式を求め，更に非線形連立方程式を解くことにより，楕円の回転角(傾斜)，長軸直径，短軸直径を求めています．

実行用バッチコマンド

gfortran -o f90_SREG.exe f90_SREG.f90
gfortran -o f90_ELLIP.exe f90_ELLIP.f90

f90_SREG inp_reg1.txt 0.95 > out_reg1.txt
gawk "15<=NR{print $2,$3}" out_reg1.txt > inp_elp1.txt
f90_ELLIP inp_elp1.txt 5 5 5 3 > _elp11.txt
f90_ELLIP inp_elp1.txt 5 5 5 5 > _elp12.txt

コマンドライン引数の並びは以下のとおりです．

f90_SREG fnameR pp > fnameW

f90_SREG	コンパイルされたF90実行ファイル
fnameR	入力ファイル名
pp	描画する等確率楕円の非超過確率
fnameW	出力ファイル名

f90_ELLIP fnameR LX LY dx,dy > fnameW

f90_ELLIP	コンパイルされたF90 実行ファイル
fnameR	入力ファイル名
LX	グラフにおけるx軸長さ(単位：cm)
LY	グラフにおけるy軸長さ(単位：cm)
dx	x軸範囲 (x軸最大値マイナスx軸最小値)
dy	y軸範囲 (y軸最大値マイナスy軸最小値)
fnameW	出力ファイル名

ソースコードと関連文書

Filename	Description
f90_SREG.txt	F90プログラムソース
f90_ELLIP.txt	F90プログラムソース
inp_reg1.txt	入力データサンプル
TeX_ellipse.pdf	プログラム解説と出力事例

クラスター分析 (K-means++ 法)

プログラム概要

クラスター分析のプログラムです．
K-means++ 法を用いています．
k 個の初期中心を得るためのアルゴリズムを以下に示します．

(0) クラスタ数 k を定め，3個の距離を格納するベクトル {D₁},{D₂},{D_R} を準備します．
(1) 1個目の中心点 c₁ を乱数を用いて設定します．
(2) ベクトル {D₁} に大きな数値をダミーとしてセットします．
(3) 中心点 c₁ から全点までの距離を計算し，それらをベクトル {D₂}に格納します．
(4) {D₁} と {D₂} の要素を比較し，小さいほうの値を {D_R}に格納します．
(5) 2個目の中心点 c₂ を {D_R}の中で最大の距離を有する点として設定します．
(6) {D₁}={D_R} とおいて，(3) から (6) までのステップをk個の中心点を得るまで繰り返します．

入力データ書式

入力データ書式は，主成分分析のものと同一です．

 001 | iris
 002 | 4,150
 003 | 5.1,3.5,1.4,0.2,setosa,1
 004 | 4.9,3.0,1.4,0.2,setosa,2
 005 | 4.7,3.2,1.3,0.2,setosa,3
..........

001: コメント
002: データ列数 (=4), 入力データ組数 (=150)
003: (この場合4列の入力データ：5列目および6列目は無効)
004: (3行目以降入力データが続く

実行用バッチファイル

gfortran -o f90_KMEANSPP.exe f90_KMEANSPP.f90

f90_KMEANSPP 0 3 30 inp_iris.csv out_iris0.txt
f90_KMEANSPP 2 3 30 inp_iris.csv out_iris2.txt

del fig_R_CLUSinp.pdf
del fig_R_CLUSpca_0.pdf
del fig_R_CLUSpca_2.pdf
del fig_R_CLUSmds_0.pdf
del fig_R_CLUSmds_2.pdf

rscript R_CLUSinp.R inp_iris.csv
ren Rplots.pdf fig_R_CLUSinp.pdf

rscript R_CLUSpca.R out_iris0.txt
ren Rplots.pdf fig_R_CLUSpca_0.pdf
rscript R_CLUSpca.R out_iris2.txt
ren Rplots.pdf fig_R_CLUSpca_2.pdf

rscript R_CLUSmds.R out_iris0.txt
ren Rplots.pdf fig_R_CLUSmds_0.pdf
rscript R_CLUSmds.R out_iris2.txt
ren Rplots.pdf fig_R_CLUSmds_2.pdf

実行用コマンド書式

f90_KMEANSPP imethod kk mds fnameR fnameW

f90_KMEANPP	コンパイルされた実行ファイル
imethod	距離計算方法 (0: マハラノビス, 1: 市街地距離, 2: ユークリッド距離, >=3: ミンコフスキー距離)
kk	クラスタ数
mds	クラスタ内最小要素数
fnameR	入力ファイル名
fnameW	出力ファイル名

ソースコードと関連文書

Filename	Description
f90_KMEANSPP.txt	f90プログラムソース
inp_iris.txt	入力データ事例
out_iris0.txt	出力データ事例 (Mahalanobis)
out_iris2.txt	出力データ事例 (Euclidean)
R_CLUSinp.txt	作図用Rスクリプト (input data)
R_CLUSpca.txt	作図用Rスクリプト (PCA)
R_CLUSmds.txt	作図用Rスクリプト (MDS)
fig_R_CLUSinp.pdf	出力画像事例 (input data)
fig_R_CLUSpca_0.pdf	出力画像事例 (Mahalanobis,PCA)
fig_R_CLUSpca_2.pdf	出力画像事例 (Euclidean, PCA)
fig_R_CLUSmds_0.pdf	出力画像事例 (Mahalanobis,MDS)
fig_R_CLUSmds_2.pdf	出力画像事例 (Euclidean, MDS)

多次元尺度法

プログラム概要

古典的多次元尺度法のプログラムです．
解析の手順は以下のとおりです．

(0) n x m の次元を持つ観測データを入力することができます．この場合，以下の手順に先立って次元 n x n の距離行列を作成する必要があります．
(1) 次元 n x n の距離行列 [d] の入力．
(2) 距離行列 [d] の各要素を自乗した行列 [D] の作成．
(3) 距離行列の要素を自乗した行列[D] より Young - Householder 変換により行列 [B] を作成．
ここに，[B]=-0.5*[J][D][J], 中心化行列 [J]=[I]-1/n*{1}{1}^T
{1}={1, 1, 1, ... ,1, 1}^T
(4) 行列 [B]の固有値 {λ_i} (λ₁>=λ₂>=λ₃>= ... >=λ_n) および固有ベクトル [E] の計算．
(5) 座標行列の計算：[X_m]=[E_m][Λ_m]^1/2.
ここに, [Λ_m] は固有値を並べた対角行列であり，その要素 λ_{i=1 to m} は 0 以上の値である必要があります．
(6) 得られた座標プロットの実行．

入力データ書式

ベクトルデータの入力 (inp_iris.csv)

この場合の入力データ書式は，主成分分析のものと同じです．

 001 | iris
 002 | 4,150
 003 | 5.1,3.5,1.4,0.2,setosa,1
 004 | 4.9,3.0,1.4,0.2,setosa,2
 005 | 4.7,3.2,1.3,0.2,setosa,3
..........

001: コメント
002: データ列数 (=4), 入力データ組数 (=150)
003: (この場合4列の入力データ：5列目および6列目は無効)
004: (3行目以降入力データが続く

距離行列の入力 (inp_map.csv)

 001 | Distance
 002 | 9
 003 | KL       , 0
 004 | Tokyo    , 5332.837,  0
 005 | Jakarta  , 1184.231,  5774.036,  0
 006 | Bangkok  , 1179.106,  4605.718,  2308.805,  0
 007 | Vientiane, 1644.533,  4151.984,  2710.010,  513.445,   0
 008 | Ha_Noi   , 2033.150,  3679.648,  3010.784,  979.078,   480.029,  0
 009 | Seoul    , 4612.781,  1162.077,  5275.036, 3715.257,  3223.101,  2743.357,   0
 010 | Beijing  , 4339.820,  2113.622,  5197.282, 3284.870,  2771.456,  2321.983,   964.494,  0
 011 | Maynila  , 2480.550,  2997.510,  2776.534, 2203.659,  2000.776,  1757.266,  2614.778,  2846.152, 0

001: コメント
002: アイテム数(データ行数)
003-011: この場合，各都市間の距離

実行用バッチファイル

ベクトルデータ入力時の例 (inp_iris.csv)

gfortran -o f90_MDS.exe f90_MDS.f90

f90_MDS 1 inp_iris.csv out_MDS_iris.csv

del fig_R_MDS_inp.pdf
del fig_R_MDS_out.pdf
rscript R_DRAWinp.R
ren Rplots.pdf fig_R_MDS_inp.pdf
rscript R_DRAWout_1.R
ren Rplots.pdf fig_R_MDS_out.pdf

距離行列入力時の例 (inp_map.csv)

gfortran -o f90_MDS.exe f90_MDS.f90

f90_MDS 0 inp_map.csv out_MDS_map.csv

del fig_R_MDS_map.pdf
rscript R_DRAWout_0.R
ren Rplots.pdf fig_R_MDS_map.pdf

実行用コマンド書式

f90_MDS knor fnameR fnameW

f90_MDS	コンパイルされた実行ファイル
knor	データ取り扱い (0: 距離行列入力, 1: ベクトルデータ入力 - 正規化実行, 2: ベクトルデータ入力 - 正規化しない
fnameR	入力データファイル
fnameW	出力データだいる

ソースコードと関連文書

Filename	Description
f90_MDS.txt	f90ソースコード
inp_iris.txt	入力データ事例 (ベクトルデータ 'iris')
out_MDS_iris.txt	出力データ事例 'iris'
inp_map.txt	入力データ事例 (距離行列)
out_MDS_map.txt	出力データ事例 'map'
R_DRAWinp.txt	作図用R スクリプト事例(ベクトルデータ描画)
R_DRAWout_0.txt	作図用R スクリプト事例 (距離行列：map)
R_DRAWout_1.txt	作図用R スクリプト事例 (ベクトルデータ：iris)
fig_R_MDS_inp.pdf	出力画像事例 (input vector data)
fig_R_MDS_out.pdf	出力画像事例 (output for vector data)
fig_R_MDS_map.pdf	出力画像事例 (2-dim. map). 下のGMT出力参照のこと
fig_G_Malay.png	GMT出力画像 (Kuala Lumpur中心)