Introduction
R を用いて、推定された B-Spline 回帰曲線の可視化などが可能です。
※現時点で R による回帰曲線の可視化は連続値変数のみのデータの場合だけに対応しています。
準備
TXT 形式によるネットワーク出力
R 解析用のスクリプトは TXT 形式のネットワークファイルを使います。ING 形式や SGN3 形式で保存している場合は TXT 形式に変換してください。TXT 形式での出力の際には full オプションを指定してください（--output-args full など）。可視化・解析に必要な情報が全て出力されます。枝がついてないノードも子が NA として出力されます。
–output-data オプションの指定
INGOR に --output-data file_prefix オプションを指定すると、file_prefix.X、file_prefix.Y、file_prefix.PR.Y の三つのファイルが出力されます。それぞれ R で読み込み可能なタブ区切りテキストによる数値データで、各行が TXT ファイル形式の各行（各枝）に対応する、説明変数の値、目的変数、部分残差のデータで、内部的に使用された値になります。したがって、行数は枝数、および列数がサンプル数に対応します。
ダイナミックモデルを使用した場合、始点および終点は目的変数および説明変数としては使用されませんので、入力ファイルのサンプル数より少ない値になります。
目的変数データは親が複数ついた場合は、回帰曲線とは合わなくなります。したがって、回帰モデルの可視化の際には説明変数データと部分残差データのみを使用します。
ブートストラップ結果を用いた B-Spline 回帰曲線可視化、尤度計算について
ブートストラップ結果の回帰曲線を可視化する場合、(1) ブートストラップ結果の枝に対して回帰しなおす、(2) B-spline による回帰モデルの範囲を固定し、推定された B-spline 係数自体を平均化する、(3) ブートストラップ結果の回帰をブートストラップで行う、の三通り方法がります。
注意：通常のブートストラップまとめ処理結果のネットワークファイルの B-spline モデルの情報は後の解析には使えません。
(1) B-spline 回帰モデルの再フィット
ScoreFilter (--score) フィルターを用いると、与えられたネットワーク構造とデータに対して、スコアを再計算しますが、その際に B-spline によるノンパラメトリック回帰を行い、結果をネットワークに残します。
実行例：
$ ingor --read input.sgn3 --score file=data.txt -o result.sgn3
この例では input.sgn3 を読み込んで、そのネットワーク構造に対し、B-spline ノンパラメトリック回帰モデルの再フィッティングを行い、その結果を result.sgn3 として SGN3 形式で保存します。
(2) ブートストラップ推定時の B-spline 係数を用いた回帰曲線可視化
ブートストラップ推定時に --fix-range オプションを推定時に指定してください。 推定時は --output-data オプションはつけないでください。 まとめ処理結果は SGN3 形式で保存してください。
--fix-range オプションにより、推定される B-Spline の定義域が固定されるので、係数の平均値を取ることが可能になります。
実行例：
$ ingor -B 1 -N 1000 --fix-range input.txt -o result -t sgn3
$ ingor --bs prefix=result,type=sgn3,ed=1 -o result.sgn3
次に、--pr オプション (PRFilter) を用いて、部分残差を計算させます。
実行例：
$ ingor --read file=result.sgn3 --pr file=input.txt,output=out -o result.txt -O full
この例では out.X, out.Y, out.PR.Y の３つのファイルが出力されます。
(3) B-spline 係数のみをブートストラップ法で求める
B-spline 係数のみブートストラップで推定する場合、-A オプションに "skip" を指定し、入力の（構造を固定する）ネットワークを --fixed オプションで指定します。 その後、まとめ処理を行う際に BSFilter に "da_mode=median" を指定すればブートストラップで推定した係数群の中央値が最終結果になります。 指定しなければ平均値を利用します。
実行例：
# ネットワーク構造の推定 (--fix-range は不要)
$ ingor -B 1 -N 1000 input.txt -o result -t sgn3
# まとめ処理（ネットワーク構造の決定）
$ ingor --bs prefix=result,type=sgn3,ed=1,th=0.1 -o result.sgn3
# B-spline 係数のブートストラップ推定 （--fixed, --fix-range 及び -A skip が必要）
$ ingor -B 1 -N 1000 --fixed result.sgn3 --fix-range input.txt -o result2 -t sgn3 -A skip
# まとめ処理２ （以下の例ではb-spline係数の中央値を求める。da_mode を指定しなければ平均値を求める）
$ ingor --bs prefix=result2,type=sgn3,ed=1,da_mode=median -o result2.sgn3
# 部分残差及び可視化用ファイルの出力
$ ingor --read file=result2.sgn3 --pr file=input.txt,output=out -o result3.txt -O full
注意事項
- ダイナミックモデルでは
--bs 及び --pr オプションに dynamic を指定してください。
- ブートストラップ結果による B-Spline 回帰曲線は実際のデータに対して適用した結果ではないため、 下で説明する部分残差に回帰曲線を重ねても、散布図とはずれて表示されます。特に、ブートストラップ確率の低い枝や、 あるいはブートストラップまとめ処理時に大きめの閾値を指定した場合、大きくずれます。 これは期待された結果ではないかもしれません。部分残差ではなく目的変数データに対して回帰曲線を描画したほうが、 両者は一致するかもしれません。
R スクリプトによる回帰曲線可視化
以下の例は TXT フォーマットのネットワークファイル result.txt の解析例です。適宜読み替えてください。また、回帰データファイル名は out.X および out.PR.Y と仮定しています。
#
# スクリプトの読み込み
source("ingor.r")
# ネットワークの読み込み
net <- read.TXT("result.txt")
# 回帰データの読み込み
X <- data.matrix(read.table("out.X"))
Y <- data.matrix(read.table("out.Y"))
R <- data.matrix(read.table("out.PR.Y"))
# １番目の枝の取り出し
edge <- net[1,]
# 以下、枝やのインデックスが含まれるので枝番号を指揮する必要なし
# B-spline の計算
B <- bspline.values.edge(edge)
# 可視化
plot(X[edge$id,],R[edge$id,],xlim=minmax(X[edge$id,]),ylim=minmax(R[edge$id,]))
par(new=T)
plot(B$x,B$y,col=2,type="l",xlim=minmax(X[edge$id,]),ylim=minmax(R[edge$id,]),ann=F)
# 部分残差ではなく目的変数データと重ねるには R を Y に置き換えます。
plot(X[edge$id,],Y[edge$id,],xlim=minmax(X[edge$id,]),ylim=minmax(Y[edge$id,]))
par(new=T)
plot(B$x,B$y,col=2,type="l",xlim=minmax(X[edge$id,]),ylim=minmax(Y[edge$id,]),ann=F)
#
尤度計算
#
# 尤度計算用の情報の抽出
model <- extract.TXT(net)
# 尤度計算用サンプルデータの作成。
x <- bspline.midvalues(model)
y <- x
# 尤度計算
bspline.loglikelihood(x,y,model)
#
その他の解析
R での解析例
# 特定の親・子の抽出
# 子が g1 である枝の抽出
P <- net[net$Child=="g1",]
#
P$Parent
#
INGOR 日本語ドキュメント
