*fieldtripでのMRIの処理の流れ(理解は不十分ですが)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 概要 =自分の理解の範囲で ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
MEGの信号源を求めるのにMRIが必要ですが、役割としては
(参照 tutorialのConstruct a headmodel for MEG source analysisのIntroとかbackground)
1. volume conduction model (head model)を求める
bemとか言ってるのはこのことのようです。ちなみにbemは計算法(boundary element model)の名前のようです。
bemを使ってsingle shell modelを計算する、みたいな感じだと思います。
* 脳波で問題になりやすいです(と思います)が、信号源(大脳)からセンサー(or電極)に行くまでに信号がゆがむ程度が脳と骨と皮膚とかで違うので、それを分けて情報として保存する必要があるということ。
* MEGでは組織の差はほとんどない(透磁率が一定)ので大脳+小脳がどんな形してるかがわかればよい=single shellモデル
* 脳波では分ける必要がある(透電率が異なる)=3 compartment modelとか、そんな感じ
fildtripではvolとかhdmとかで変数を格納していることが多いです。=volとかhdmを求める作業がvolume conduction modelの作業
2. source modelを求める
1との区別はしばしばサラッと流されてますが、信号源をどこに置くか、ということだと思います。
beamformer(特にfieldtrip上)では1のvolの情報とsource modelは一致しますが、
electaのbeamformerではsourceから小脳、脳幹を外したり、mneでは大脳表面だけにsourceを置いたりするので、必ず一致するわけではないです。
3. forward solution
信号源(大脳)を細かく分割して1つの小さな信号源(beamformerの時は各ボクセル)を分けたときに、各信号源の電流が1,2を勘案したときにどのようにMEGセンサーに記録されるか、を意味するのだと思います。
leadfield行列(gridという名前で格納していることが多い)を求める、という作業と一致する(と思います)。
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ fieldtripでの処理の流れ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
参照 wiki.cimec.unitn.it/tiki-index.php?page=ForwardModels (CIMEC?のwiki)
fieldtrip.fcdonders.nl/tutorial/headmodel_meg (fieldtripのtutorial)
1. vol (hdm)を求める=head model (volume conduction model)を求める
やることは
(1) mriを読み込んで (ft_read_mri)
(2) nasion, LPA, RPAをポイントして座標系を指定する (ft_volumerealign)
(3) mriから脳とかを分離する(ft_volumesegment)
(4) vol(hdmも同じ)を求める(ft_prepare_singleshell)
(1) mriの読み込み
mri=ft_read_mri('orig.mgz')
* 経験上の注意点:dicomから読み込むのはうまくいかないことが多かった。freesurferで処理して出てくるorig.mgzとかT1.mgzあたりを読むのはトラブルが少ない(もちろんfreesurferが組み込まれている必要があるが)
(2) nasion, LPA, RPAのポイント(座標系の指定)
cfgRlgn=[];%cfgでももちろんok
cfgRlgn.method='interactive';%guiで指定するということ。(通常はないが)前もって座標がわかっていればそれを使えるらしい。
cfgRlgn.coordsys='neuromag';%座標系.指定しないと'ctf'になる
*経験上の注意点:fieldtripのパスにspmのパスが含まれているが、こちらを使った方がよかった
=spm8にパスを通してうまくいかないことがあった=spm8のパスを外して処理する必要があるかも知れない
mriRlgn=ft_volumerealign(cfgRlgn,mri);%ウィンドウがポップアップしてnasion, LPA, RPAを指定する
*経験上の注意点:matlabのウィンドウに作業が出てくるが、十字の線が出てくるのでそれをnasionに合わせて'n'をタイプ、LPAに合わせて'l'をタイプ・・みたいにする。nasion, LPA, RPAの順番は適当でよい。
最後作業が終わったらウインドウの上で'q'をタイプ。ウィンドウが消える。
! ウィンドウが消えなかったら作業をチョンボしたと思った方がよいです。
(3) segmentation
cfgSeg=[];%cfgでももちろんok
cfgSeg.coordsys='neuromag';
segmentedmri=ft_volumesegment(cfgSeg,mriRlgn);
* fieldtripのtutorialではcfg.output='brain'の指定があるが、なくても動いた
(4) head model (volume conduction model)の計算=single shell model (MEGの場合)の計算
cfgVol=[];%cfgでもok
vol=ft_prepare_singleshell(cfgVol, segmentedmri);
* 余談 singleshellでないsphere modelとかも計算はできて、その時はft_prepare_headmodelにオプション指定するとよいみたい
-------------------------- うまくvolができてるか確認する方法があります -----------------------------
確認したら脳が横を向いていたことがありましたorz。
原因はcoordsys='neuromag'が'ctf'になっていたため、でしたが、何事も確認が必要です。
脳が変な形に切り取られていることもありましたorz。LPAとRPAを逆に指定していました。
(1) MEGのデータ(たとえばhoge_raw.fifとします)にMEGセンサー、polhemus(スタイラスで登録した位置情報)が入っていますので、それを読み込みます。
hs=ft_read_headshape('hoge_raw.fif');% polhemus
sens=ft_read_sens('hoge_raw.fif');%MEGセンサー
(2) プロットします (順番は適当で大丈夫)
ft_plot_vol(vol);%脳が表示されます
ft_plot_sens(sens,'style','b*');%MEGセンサーが表示されます
ft_plot_headshape(hs);%polhemusのデータが表示されます
matlabのfigureの機能でぐるぐる回すことができるのでちゃんと向きとか位置関係が合ってるか確認しましょう。
2. source modelを求める
これに相当する作業は不要のようです。1のvolを再利用するということのようです。
3. forward solution = leadfield行列を求める
(1) 前準備:MEGのデータの情報を取り出す(たぶんMEGセンサーの位置とかpolhemusの位置情報とか必要なので)
cfgRaw=[];
cfgRaw.dataset='hoge_raw.fif';
data_raw=ft_preprocessing(cfgRaw);
(2) leadfield行列の計算 (ft_prepare_leadfield) fieldtripのtutorialそのまんまです
cfgGrid = [];
cfgGrid.grad = data_raw.grad;
cfgGrid.vol = vol;
cfgGrid.reducerank = 2; % ここは意味がわかりません
cfgGrid.channel = 'MEG';
cfgGrid.grid.resolution = 1; % use a 3-D grid with a 1 cm resolution
[grid] = ft_prepare_leadfield(cfgGrid);
gridとかgradとか出てきてわかりにくいですね。
-------------------------- 参考: CIMECのtutorialでは ---------------------------------------
雰囲気でだいたいわかると思いますが、こちらではreducerankの指定をしていないのが気になります。不要?
cfg=[]; cfg.vol=hdm; cfg.channel=gradlabels; cfg.grid.xgrid='auto'; cfg.grid.ygrid='auto'; cfg.grid.zgrid='auto'; cfg.grid.resolution=1; %1cm ... beware of units!! cfg.grad=data_raw.grad; lf=ft_prepare_leadfield(cfg);