hokan 関数を使い１→100までの数列を作るとき、ちょっと注意することがあります。それは得られた数列が見かけ上、整数に見えても実は誤差を 含んでいるため注意する必要があります（mc で取り扱うのは実数）。その誤差の大きさは大体 10^-14 程度で非常に小さいのですが、この値をピクセルの位置座標に使用する tensya 関数では大きな問題となります。

xx = |1,100|';
xx = hokan(xx,100,1);
yy = |1,100|'+0.5;
yy = int(hokan(yy,100,1));
zz = selectIf(xx-yy,(xx-yy)<>0);
plot(zz);

１〜100の整数のはずが44個も誤差を含んでいるものがある

　この問題を解決するには上記プログラムの yy の項目にも書いているよう、四捨五入を利用します（全体に 0.5 を加え int 関数を使って整数化する）。それでは tensya 関数を使った解決例を紹介しましょう。

xx = |1,256||1,256|; ~２行２列の行列を定義
xx = hokan(xx,256,256); ~256×256に補間する
yy = xx'; ~xxの転置行列をy軸データとする
zz = size(256,256)+1; ~転写元データ（全て１）
xx = changeMN(xx,1); ~１列のデータに変換：これはx軸データ
yy = changeMN(yy,1); ~１列のデータに変換：これはy軸データ
zz = changeMN(zz,1); ~１列のデータに変換：これはz軸データ（濃淡）
cc = addP(xx,yy); ~xx,yy,xxを並べ３列のデータにする
cc = addP(cc,zz);
e = tensya(cc[#,1],cc[#,2],cc[#,3],|256,256|,0); ~転写データを作成
densityPlot(e,|1,0|);

　本来であれば転写後の行列eはすべて１になるはずで、得られた濃淡図は白一色のはずです。ところが座標データとして使用するx軸・y軸の補間処理で誤差 が発生しています。簡単な例で示すと、本来は２なのだが実際のデータは1.9999999999であり整数化すると目的の２ではなく１になってしまうとい う問題です。このような問題から上記図のよう tensya 関数で転写できなかった領域が残り黒 い線として残ってしまいます。

xx = |1,256||1,256|+0.5; ~各要素一律に 0.5 を加える
xx = int(hokan(xx,256,256)); ~補間後に整数化する
yy = xx'; ~以下は同じ
zz = size(256,256)+1;
xx = changeMN(xx,1);
yy = changeMN(yy,1);
zz = changeMN(zz,1);
cc = addP(xx,yy);
cc = addP(cc,zz);
e = tensya(cc[#,1],cc[#,2],cc[#,3],|256,256|,0);
densityPlot(e,|1,0|);

何も見えないが！　正常な画像

　整数化が四捨五入で問題なく実行できたので転写後の行列eの要素全てが１になり得られた画像は白一色となっている。本来であれば tensya 関数の方で直すべきなのかもしれない（どうしようかな〜）。

　ここで示す例は重ね合わせ位置調整済みの SPECT 画像データと MRI データを使いピクセル・バイ・ピクセルのSPECT値とMRI値の関連を分析するものです。

　文字列変数 f1$ にはSPECTデータのファイル名が、 f2$ には MRI T1 強調画像データのファイル名が保存されています。まずは画像データの読み込みと分析する画像を表示します。

s = 7;
p = loadVolume(f1$,1);
p = loadVolume(f2$,2);
changeCurrentVolumeData(1);
a = getVolumeData(s);
RGB = popRGBColor(0);
a# = densityPlot(a,|500,0|);
changeCurrentVolumeData(2);
b = getVolumeData(s);
RGB = popRGBColor(1);
b# = densityPlot(b,|500,0|);
c# = a#*b#*(a#<*>b#);
changePictSize(c#,1/2);

左から SPECT 、MRI、合成画像

　次に SPECT 値と MRI 値の関係を見るためピクセル・バイ・ピクセルで散布図を描いてみましょう。

aa = changeMN(a,1);
bb = changeMN(b,1);
pointMode(1);
setPlotWindowSize(400,300);
graphTitle(1,"SPECT",|180,283|);
graphTitle(2,"MRI",|50,35|);
parametricPointPlot(aa,bb)+0;

この散布図の中に 65536個の点が描かれている

　プログラム最後の作図関数に+0という項目が追加されています。これは得られた画像がDrawデータなのでそのまま使うと表示に時間がかかるためです（65536個の点を表示 の度に描く）。表示を高速化させる目的でゼロを加えDrawデータではなくビットマップデータに変換しています。このように画像をビットマップデータにすると表示が速くなります。

　さてSPECTとMRIの関係ですが散布図によると L 字型に分布しているのが判りました。それでは特定の領域にあるピクセルが画像のどの部位に相当するか調べてみましょう。　まずはピクセルデータに座標情報を付加する前処理です。

xx = |1,256||1,256|;
xx = hokan(xx,256,256);
yy = xx';
xx = changeMN(xx,1);
yy = changeMN(yy,1);
cc = addP(aa,bb);
cc = addP(cc,xx);
cc = addP(cc,yy);
cc = selectIf(cc,cc[#,1]>50 or cc[#,2]>100); ~小さい値の要素を削除
setPlotWindowSize(400,300);
graphTitle(1,"SPECT",|180,283|);
graphTitle(2,"MRI",|50,35|);
parametricPointPlot(cc[#,1],cc[#,2])+0;

値が小さいピクセルを selectIf 関数で取り除いたのでピクセル数は数千程度

　これで前処理は終わりです。次は grEdit 関数（この関数は未だ未公開）を使い散布図の特定の領域のデータだけを選択してみます。

dd = grEdit(cc,1,2);

散布図の領域をマウス操作で囲うと、囲んだ領域にあるデータが選択される。
ここではMRI 値が大きいが SPECT値はそう大きくない領域を選択し行列ddに代入する。

e = tensya(dd[#,3],dd[#,4],dd[#,2],|256,256|,0);
densityPlot(e,|500,0|);

　散布図の選択した領域にあるピクセルを元の座標に戻してやると脳の外側が見えてきます。しかし、画像に黒い線が入り汚く表示されています。これが位置座標を作るときの補間で生じた誤差の影 響です。それでは四捨五入を使い修正してみましょう。

　四捨五入により位置座標に誤差が無くなったので黒い線が無くなりました。

dd = grEdit(cc,1,2);

今度は別の領域を選択

e = tensya(dd[#,3],dd[#,4],dd[#,2],|256,256|,0);
densityPlot(e,|500,0|);