[VBMEG-users:00162] current variance estimationにおけるArea IDの選択について
shibusawa
ご担当者様
お世話になっております。
慶應義塾大学大学院、理工学研究科修士1年の澁澤と申します。
お伺いしたいことが二点ございます。
まず一点目は、current variance estimationにおける、”Area ID”の選択についてです。
推定に用いる自身のEEGデータが、Sample frequency=250 Hz, Nsample=150,000（10min）と大きいため、”.area.mat”に少数のダイポールから成るROIを作成し、そのエリア内で”Dipole reduction ratio=1.00”で推定を行っています。
そこで、current variance estimationの際に”Area ID”にROIを指定している場合は、最終的な”.curr.mat”におけるZactでROIのダイポールに推定されていることが確認されました。一方で、”Area ID”に”Cortex”を指定すると、ZactでROIよりも広い範囲のダイポールが推定されることが確認されました。また、それら二つでは同じダイポールにおける推定信号が異なることが確認されました。
電流分散を計算する際のエリア指定は、どのような意味になるのでしょうか。
二点目は、”.eeg.mat”における”Pretrigger”についてです。
“trigger onset”を登録するということは、二状態のデータに対して推定を行うということでしょうか。その場合、その区別における数学的意味を教えて頂きたく思います。
また、一状態で計測されたデータに対して推定を行うことは可能でしょうか。
以上の点につきまして、ご教示頂けたら幸いです。
宜しくお願い致します。
------------------------------------------------
慶應義塾大学大学院理工学研究科
基礎理工学専攻 牛場研究室
修士1年 澁澤柊花
shib...@brain.bio.keio.ac.jp
------------------------------------------------
_________________________________________________________
VBMEG users mailing list
vbmeg...@lists.osdn.me
Options:
http://lists.osdn.me/mailman/listinfo/vbmeg-users
Archives:
http://groups.google.com/group/vbmeg-users
Yamashita Okito
ATRの山下です。
> まず一点目は、current variance estimationにおける、”Area ID”の選択についてです。
> 推定に用いる自身のEEGデータが、Sample frequency=250 Hz, Nsample=150,000（10min）と大きいため、”.area.mat”に少数のダイポールから成るROIを作成し、そのエリア内で”Dipole reduction ratio=1.00”で推定を行っています。
> そこで、current variance estimationの際に”Area ID”にROIを指定している場合は、最終的な”.curr.mat”におけるZactでROIのダイポールに推定されていることが確認されました。一方で、”Area ID”に”Cortex”を指定すると、ZactでROIよりも広い範囲のダイポールが推定されることが確認されました。また、それら二つでは同じダイポールにおける推定信号が異なることが確認されました。
> 電流分散を計算する際のエリア指定は、どのような意味になるのでしょうか。
す。電流源の活動位置に関する信頼できる事前情報があるときに、場所を限定し
（未知パラメータを減らすことによって）、false positiveを抑制することがで
きます。
プログラム仕様としましては、Area IDによって指定された頂点集合のみ電流分
散が推定され、それ以外の頂点では0となります。
Area ID = "Cortex" は、デフォルトで登録される ID で、皮質上の全頂点集合
として登録されてます。
従って、電流分散推定の時に、Area ID = "Cortex"を指定すると、全皮質頂点に
おいて、電流分散値が推定されます。
> 二点目は、”.eeg.mat”における”Pretrigger”についてです。
> “trigger onset”を登録するということは、二状態のデータに対して推定を行うということでしょうか。その場合、その区別における数学的意味を教えて頂きたく思います。
算する信号時間窓を指定する必要があります。
課題時の脳活動計測では、一般にノイズ時間窓として刺激提示や課題時直前の特
に課題をしていない時間帯を指定するため、その時間帯を識別するために、
trigger情報は必要です。
> また、一状態で計測されたデータに対して推定を行うことは可能でしょうか。
レスティング課題のようなケースを想定されていますか？
観測ノイズとして適切な情報を含んだファイル(MEGであれば被験者が入っていな
い時のノイズファイルなど)を、megfile_baseline に指定し、twin_noise に時
間インデックスを指定すれば可能です。
山下
shibusawa
ご多忙中のところ、ご回答ありがとうございます。
一点目につきまして、回りくどい表現となってしまい誤解をお招きしてしまいました。申し訳ございません。今一度、現状を述べさせていただきます。
現状で確認したのは以下の二点です。
目標：全皮質上20004点のダイポール中の、ROI99点について信号を推定したい
1）current variance estimationにてArea IDに"ROI"を選択(dipole reduction
ratio=1.00)→ROIの99点の電流分散が計算された
estimate currentにてArea IDに"ROI"を選択→ROIの99点の信号が推定された
2) current variance estimationにてArea IDに"Cortex"を選択(dipole reduction
ratio=1.00)→263点の電流分散が計算された
esimate currentにてArea IDに"ROI"を選択→上と同じ263点の信号が推定された
ここで、2)の際は J=WZにおいて、Z(Zact): 263×T、W(Wact): 263×99となっております。
2点目でご指摘頂いた通り、この度はレスティング課題のデータ用いて推定したいと考えております。そのため、電流源の活動位置に関する信頼できる事前情報がございません。
よってArea IDは"cortex"を選択して電流分散を推定するべきだと考えておりますが、上記のような状況となり、ご教示頂きたくご質問させて頂きました。
二点目につきましては、ご教示頂いた方法を試してみたく思います。
お忙しいところ重ね重ねお伺いしてしまい恐縮ですが、何卒よろしくお願い致します。
------------------------------------------------
慶應義塾大学大学院理工学研究科
基礎理工学専攻 牛場研究室
修士1年 澁澤柊花
shib...@brain.bio.keio.ac.jp
------------------------------------------------
Yamashita Okito
> 2) current variance estimationにてArea IDに"Cortex"を選択(dipole reduction
> ratio=1.00)→263点の電流分散が計算された
current variance
estimationで263点のみで電流分散が推定されることは考えづらいです。
current variance
estimation以前の処理工程に何か間違いがある可能性があります。
以下の点を確認してみてください。
(1) brain model の頂点数
V= vb_load_cortex(brainfile);
でVを読み出す。Vの次元はいくつですか？
(2) "Cortex"に登録された頂点情報
Area = vb_get_area(areafile, 'Cortex');
でAreaを読み出す。 Area.Iextract はどうなってますか?
山下
shibusawa
大変失礼致しました、電流分散が263点というのは誤りでした。お騒がせして申し訳ございません。
(1) brain model の頂点数
V= vb_load_cortex(brainfile);
でVを読み出す。Vの次元はいくつですか？
→20004×3となりました。
(2) "Cortex"に登録された頂点情報
Area = vb_get_area(areafile, 'Cortex');
でAreaを読み出す。 Area.Iextract はどうなってますか?
→20004×1となりました。
正しい状況は以下となります。
1）current variance estimationにてArea IDに"ROI"を選択(dipole reduction
ratio=1.00)→ROIの99点の電流分散が計算された
estimate currentにてArea IDに"ROI"を選択→ROIの99点の信号が推定された
2) current variance estimationにてArea IDに"Cortex"を選択(dipole reduction
ratio=1.00)→全脳皮質20004点の電流分散が計算された
esimate currentにてArea IDに"ROI"を選択→263点の信号が推定された
ここで、2)の際は J=WZにおいて、Z(Zact): 263×T、W(Wact): 263×99となっております。
Wactは、 ".curr.mat"の"Jinfo.Wact"を確認しております。
何度もお手を煩わせてしまい申し訳ございません。
ご教示よろしくお願い致します。
Yamashita Okito
> 2) current variance estimationにてArea IDに"Cortex"を選択(dipole reduction
> ratio=1.00)→全脳皮質20004点の電流分散が計算された
> esimate currentにてArea IDに"ROI"を選択→263点の信号が推定された
> ここで、2)の際は J=WZにおいて、Z(Zact): 263×T、W(Wact): 263×99となっております。
> Wactは、 ".curr.mat"の"Jinfo.Wact"を確認しております。
Wact
は空間平滑化フィルタです （平滑化の大きさはbayes_parm.Rfiltで指定されます)。
Jは "ROI"で指定された99個の頂点における電流Zを、空間平滑化フィルタ
をかけてぼやかしたものとなります。
ちなみに、Jの263点の脳モデル上の頂点インデックスは
Jinfo.ix_act_ex に収められています。
”高解像度の頂点で推定された電流から空間的にダウンサンプルした電流を計算する”ときに、現在渋沢さんが行っているような
Area ID
で抜き出す方法はお勧めしません。この方法は、空間的にサブサンプルしているだけであり、その他の頂点の電流値は
全て無視されます。解決策としては、ROI分割情報を利用することが考えられます。
Area ID = "ROI"
の99個の頂点は解剖学的または機能的に大脳皮質を分割した領野の代表点だと推察します。
そうであれば、99個の領野に属する頂点集合が定義できるので、それをArea file
に area_key をつけて登録します。
そして、estimate current でarea_keyごとに
その領野に属する電流を取り出し、平均値また第1主成分などで代表時系列を
計算します。
最後の代表時系列の取り出すところは、決定版がなくて難しいところですが、だいたい以上のような手続きで計算すると良いと
思います。参考まで。
山下
shibusawa
分かりづらい説明となってしまい申し訳ございません。
現状と疑問点について、再度説明させて頂きたく思います。
まず、特定のROIに該当する全頂点の各電流値を求めたいと考えております。これが最終的な目標となります。
そこで、他のソフト(spm)にて皮質上に機能的なROIを定義し、".brain.mat"と".spm.mat"を用いて".area.mat"".act.mat"を作成する工程と同様の方法で、領域を定義しています（このうちのact情報は不要で、area情報を用いる方針です）。ここで定義される頂点集合は、代表点集合ではなくROIに含まれる全頂点集合だと解釈しており、"view estimated current"でvertex indexに数値を打ち込みなが
ら確認をしました。
次に、前回お伝えした問題点について、時系列順にの記述します。
以下で言うROIは上で述べた方法で定義したものです。
＊解析1
current variance estimationにてArea IDに"Cortex"を選択(dipole reduction
ratio=1.00)→全脳皮質20004点の電流分散が計算された
esimate currentにてArea IDに"ROI"を選択→263点の信号が推定された
ここで、J=WZにおいて、Jはその領域に含まれる全頂点の各電流値、Wは空間平滑化フィルタ、Zは間引きして推定した代表点の各電流値と理解しています。
よって、サンプル数をTとすると、Jは99×T、Wは99×99、Zは99×Tとなることを期待していました。
しかしながら実際は、Z(Zact): 263×T、W(Wact): 99×263となっております。
Zactは、".curr.mat"の"Zact"を、Wactは、 ".curr.mat"の"Jinfo.Wact"を確認しております。
また、"Jinfo.ix_act"は263×1、"Jinfo.ix_act_ex"は99×1となっております。
ここで、ZactがROIより拡張した範囲の頂点集合になっていたこと、空間平滑化フィルタが頂点数を減らすサイズであったことを不思議に思い、次の方法を試しました。
＊解析2
今回はレスティング課題のデータを解析しているため、電流源の活動位置に関する信頼できる事前情報がなく、電流分散を推定する際にROIのみに推定することは、推定精度を下げる要因になると考えています。
しかし、計算処理の工程を知りたく、データ構造を確認するため、以下のように行いました。
current variance estimationにてArea IDに"ROI"を選択(dipole reduction
ratio=1.00)→ROIの99点の電流分散が計算された
stimate currentにてArea IDに"ROI"を選択→ROIの99点の信号が推定された
私の予想していた結果は、解析1の方法で、Zact,Jactの次元がともに99×Tとなることでした。
この現象が起こる原因について、ご意見を頂けたら幸いです。
お手数おかけしておりますが、よろしくお願い申し上げます。
Yamashita Okito
> ここで、ZactがROIより拡張した範囲の頂点集合になっていたこと、空間平滑化フィルタが頂点数を減らすサイズであったことを不思議に思い、次の方法を試しました。
この解釈が完全に間違っています。
estimate current で Area IDを指定した時に起こることは、
１． ROIで指定された頂点の電流を抽出 Z 99*T
２. それを空間平滑化フィルタで拡張 J = WZ 263*T
です。
山下
shibusawa
ご指摘ありがとうございます。
".curr.mat"の"Zact"は J = WZにおける J を意味しているということになるのでしょうか。
"Zact"が上の式におけるZを意味しているならば、Wは99×263ではなく263×99になると考えているのですが、その点においても間違っているのでしょうか。
Yamashita Okito
先のメールの私の説明に誤りがありました。
その訂正を含めて正確に説明したいと思いますので、やや長文になりますことご容赦ください。
話を明確にするために、プログラム上の関数名や変数名、具体的な変数のサイズを用いて記述します。
まずおさらいから。
estimate current = vb_job_current.m の機能は、J-current
を計算するために必要な
空間平滑化フィルタ W と Z-current を計算することです (J=WZ)。
curr.mat 内の変数 "Zact" は z-current となります。平滑化フィルタ情報は
"Jinfo" に収めれています。
つぎに vb_job_current内での計算ですが、主に以下の3つを行っています。
1. bayes_parm.area_key で指定された領域において、空間平滑化フィルタ W
を計算する。 (今回は全脳領域です)
2. bayes_parm.area_key で指定された領域において z-current を計算する。
（今回の場合は全脳です）
3. 1,2で計算された W, Z を curr_parm.area_key
で指定された領域に制限する。(今回の場合 ROIです)
ここで、ステップ３の部分が、今回のややこしい結果を引き起こす原因となっています。
全脳の J-current, Z-current, 空間平滑化フィルタを、 J (20004xT), W
(20004x20004), Z(20004xT), ROI の index を ix (99x1) とすると、
まずステップ1，2では
J = WZ
を満たすようにW, Jが計算されます。
次に、上記ステップ3 では、『ROI内のJ-currentの計算』に必要な値を、 W, Z
から抽出します。
つまり、次の式が満たされるように、? のインデックスを選びます。
J(ix) = W(ix, ?) Z(?)
ここで ? は、 ix
より大きい集合になるでしょうか？ 小さい集合になるでしょうか？
W として空間平滑化なし(恒等行列)のケースを考えると ? は ix
と同じになります。
Wとしてある幅を持った空間平滑化フィルタのケースを考えると、 ？ は ix
より大きくなります。
メール最後に付記しました、1次元の単純なケースを考えてもらえれば理解してもらえるかと思います。
以下、回答です。
> ".curr.mat"の"Zact"は J = WZにおける J を意味しているということになるのでしょうか。
Z です。
> "Zact"が上の式におけるZを意味しているならば、Wは99×263ではなく263×99になると考えているのですが、その点においても間違っているのでしょうか。
ROI内のJが復元されるように、Z,W が選ばれるため、非直感的ですが、W は
99x263、Zは263xT で間違いありません。
=========================================================
１次元の頂点の並びを考え、
J(n)はZ(n-1),Z(n),Z(n+1)の平均値で与えられるとします。
この時、平滑化フィルタWは
W(n-1,n) = W(n,n) = W(n+1,n) = 1/3
他のWの要素は、0です。
この時、J(2)に影響を与えるZは、Z(1),Z(2),Z(3)となります。
JはZを平滑化して得られるのですが、Jに影響を与えるZの範囲は、
拡がる結果となります。
===========================================================
山下
--
-----------------------------------------------------
Department of Computational Brain Imaging,
Neural Information Analysis Laboratories, ATR
Department Head : Okito Yamashita
E-Mail : oyam...@atr.jp
Address : 2-2-2 Hikaridai,Seika-cho,Soraku-gun,Kyoto,Japan
Phone : +81-774-95-1073
--------------------------------------------------------
shibusawa
丁寧なご説明を頂き、誠に感謝致します。
ご説明頂いた内容で、全ての疑問が晴れました。ありがとうございました。
稚拙な質問にお時間をお取りしてしまい、申し訳ございませんでした。
今後ともお伺いすることがあるかもしれませんが、どうぞよろしくお願い申し上げます。