FSLのバージョン情報がどこにあるか探したので備忘録として書いておきます。
cat /usr/local/fsl/etc/fslversion
で求められます。
FSLのバージョン情報がどこにあるか探したので備忘録として書いておきます。
cat /usr/local/fsl/etc/fslversion
で求められます。
bedpostx_gpu を走らせると、以下のエラーがでます。
/usr/local/fsl/bin/bedpostx_postproc_gpu.sh: 行 20: --cnonlinear/bin/merge_parts_gpu: そのようなファイルやディレクトリはありません
この解決法がFSLのMLで紹介されています。
https://www.jiscmail.ac.uk/cgi-bin/wa-jisc.exe?A2=FSL;ee0b1626.2112
具体的には、
${FSLDIR}/bin の中にある bedpostx_postproc_gpu.sh の
# last 2 parameters are subjdir and bindir parameters="" while [ ! -z "$2" ] do
を
# last 2 parameters are subjdir and bindir parameters="" while [ ! -z "${2+x}" ] do
に変更します。
while の後の test文 の中が、 $2 が ${2+x} になっています。
これで無事に動きます。
ご紹介まで。(金子貴久子先生、情報提供ありがとうございました)
FSL 6.0.5 ships eddy_cuda10.2 which literally uses CUDA 10.2.
I explored a way to use eddy_cuda10.2, PyTorch and Tensorflow concurrently. Below is How-To for Ubuntu 18.04/20.04.
1. 目的
2. リポジトリ
3. 必要なデータ
4. 注意点
5. 実行
5.1. run_FLICA.m
6. 結果
6.1. index.html
6.2. コンポーネントごとの結果(All_in_one_page)
6.3. モダリティごとの結果(FA, MD, AD, RD, etc.)
7. GBSSへの適応
7.1. index.html
7.2. コンポーネントごとの結果(All_in_one_page)
7.3. モダリティごとの結果(FA, MD, AD, RD, etc.)
1. 目的
2. fslstats
2.1. binary maskを使う場合
2.2. index maskを使う場合
3. 使用例
3.1. フォルダ構造
3.2. コード
3.3. 結果
1. 目的
2. JHU White-matter labels & tractography atlasでROI解析する方法
2.1. ディレクトリ構造
2.2. 実行方法
2.3. ソースコード
3. Desikan Killiany AtlasでROI解析する方法
3.1. 前提条件
3.2. ディレクトリ構造
3.3. 実行方法
3.4. ソースコード
4. FreeSurferのwm.seg.nii.gzをDWI空間に位置合わせする方法
4.1. 必要なファイル
4.2. ソースコード
5. 標準空間にあるROIを各被験者脳に位置合わせ
6. bertのwmparcをMNI空間に移動後、個人脳(Perfusion image)にレジストレーション
1. 目的
2. GBSSとは
2.1. データの準備
2.2. 解析パラメータの設定
2.3. gbss_1_T1wpreproc
:T1WIの前処理
2.4. gbss_2_DWIwpreproc
:拡散定量値の前処理
2.5. gbss_3_skelpreproc
:拡散定量値をスケルトンに投影
2.6. randomise
:スケルトンに投影された定量値画像を入力したGLMと並べ替え検定(permutation test)
3. おまけ
Gray matter-Based Spatial Statistics(GBSS)は、灰白質の統計解析をするための手法。TBSSの灰白質版。
灰白質の中心線(skeleton)に定量値を投影する。通常の脳画像の統計解析では、脳構造の個人差を除外するために空間的「平滑化」を用いる。しかし、平滑化の程度に原則がなく、平滑化をかけては情報があいまいになり、MRIの高空間分解能を生かせないという問題がある。一方、GBSSでは、灰白質の中心線と思われるところにskeletonを生成し、そこに個人ごとの定量値を投影するという手法をとる。これにより、平滑化せずに群間比較をすることができるため、平滑化による問題を回避できるという利点がある。
GBSS解析では、次のような処理をする。
gbss_1_T1wpreproc
:T1WIの前処理gbss_2_DWIwpreproc
:拡散定量値の前処理gbss_3_skelpreproc
:拡散定量値をスケルトンに投影randomise
:スケルトンに投影された定量値画像を入力したGLMと並べ替え検定(permutation test)まず、次のデータを準備する。
フォルダ構造は次のようにする。
ここでは、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)いることを想定する。各フォルダに入れるファイル名は同じにする必要がある(例:FA/Subj001.nii.gz, T1w/Subj001.nii.gzは同じファイル名)。
. ├── FA │ ├── Con0001.nii.gz │ ├── Con0002.nii.gz │ ├── ... │ └── Pat0010.nii.gz ├── FW │ ├── Con0001.nii.gz │ ├── Con0002.nii.gz │ ├── ... │ └── Pat0010.nii.gz ├── MD │ ├── Con0001.nii.gz │ ├── Con0002.nii.gz │ ├── ... │ └── Pat0010.nii.gz ├── T1w │ ├── Con0001.nii.gz │ ├── Con0002.nii.gz │ ├── ... │ └── Pat0010.nii.gz └── b0 ├── Con0001.nii.gz ├── Con0002.nii.gz ├── ... └── Pat0010.nii.gz
解析に用いるパラメータを設定する。
PEDIR
、ECHOSPACING
はepi_reg
コマンドのパラメータであり、それぞれ位相エンコード方向とecho spacing timeを指定する。SKELETON_THR
はスケルトン化する際のしきい値である。
MAP_LIST=$(ls | grep -v b0 | grep -v T1w | grep -v stats) SUBJ_LIST=$(ls T1w | cut -d . -f1) PEDIR='-y' # Phase encoding direction: "epi_reg" config ECHOSPACING='0.0380544' # Echo spacing time (s): "epi_reg" config SKELETON_THR=0.2 # Skeleton threshold
gbss_1_T1wpreproc
:T1WIの前処理T1WIの前処理内容は、次の通り。
以上の処理内容を実行するために、関数gbss_1_T1wpreproc
を定義。
function gbss_1_T1wpreproc() { echo "T1w preprocessing..." # T1w preprocessing for SUBJ in ${SUBJ_LIST}; do ## Segment GM ### Skull-stripping bet T1w/${SUBJ} T1w/${SUBJ}_skull_stripped -f 0.3 -R -S -B ### Segment T1w into CSF/GM/WM fast -t 1 -g -B -b -p -o T1w/${SUBJ}_fast T1w/${SUBJ}_skull_stripped ## Register GM to MNI ### Register T1w to MNI flirt -ref ${FSLDIR}/data/standard/MNI152_T1_2mm_brain \ -in T1w/${SUBJ}_fast_restore \ -omat T1w/${SUBJ}_indiv2MNI.mat fnirt --in=T1w/${SUBJ}_fast_restore \ --aff=T1w/${SUBJ}_indiv2MNI.mat \ --cout=T1w/${SUBJ}_indiv2MNI_warp \ --config=T1_2_MNI152_2mm ### applywarp to GM and move it into MNI applywarp --ref=${FSLDIR}/data/standard/MNI152_T1_2mm_brain \ --in=T1w/${SUBJ}_fast_pve_1 \ --warp=T1w/${SUBJ}_indiv2MNI_warp \ --out=T1w/${SUBJ}_fast_pve_1_inMNI & done ## Create 4D all_GM image ### GM_pve1 in MNI are merged into one volume mkdir stats fslmerge -t stats/all_GM $(imglob T1w/*_fast_pve_1_inMNI.nii.gz) ### Create merged GM_pve1 (all_GM) mask and apply it to all_GM fslmaths stats/all_GM -max 0 -Tmin -bin stats/mean_GM_mask -odt char fslmaths stats/all_GM -Tmean stats/mean_GM }
関数が定義できたら、実行する。
gbss_1_T1wpreproc
処理が完了すると、statsフォルダにすべての被験者の灰白質画像がまとめられた「all_GM.nii.gz」とそれを平均化した「mean_GM.nii.gz」が生成される。
gbss_2_DWIwpreproc
:拡散定量値の前処理拡散定量値の前処理の過程は、次の通り。
mean_GM_mask
で、3の拡散定量値マップをマスク処理function gbss_2_DWIwpreproc() { echo "Diffusion MAP preprocessing..." # Diffusion MAP preprocessing for SUBJ in ${SUBJ_LIST}; do ## Register b0 to MNI ### Register b0 to T1WI epi_reg --epi=b0/${SUBJ} \ --t1=T1w/${SUBJ} --t1brain=T1w/${SUBJ}_skull_stripped \ --echospacing=${ECHOSPACING} --pedir=${PEDIR} \ --wmseg=T1w//${SUBJ}_fast_seg_2.nii.gz \ --out=b0/${SUBJ}_BBR for MAP in ${MAP_LIST}; do ### applywarp to diffusion map and move it to T1w based on epi_reg flirt -in ${MAP}/${SUBJ} \ -out ${MAP}/${SUBJ}_inT1 \ -ref T1w/${SUBJ}_skull_stripped \ -init b0/${SUBJ}_BBR.mat -applyxfm ### applywarp to dMAP in T1w and move it to MNI applywarp --ref=${FSLDIR}/data/standard/MNI152_T1_2mm_brain \ --in=${MAP}/${SUBJ}_inT1 \ --warp=T1w/${SUBJ}_indiv2MNI_warp \ --out=${MAP}/${SUBJ}_inMNI done done ## 4D all_${dMAP} for MAP in ${MAP_LIST}; do ### Merge dMAP in MNI into one volume, all_${MAP} fslmerge -t stats/all_${MAP} $(imglob ${MAP}/*_inMNI.nii.gz) ### Masking dMAP in MNI using mean GM mask fslmaths stats/all_${MAP} -mas stats/mean_GM_mask stats/all_${MAP} done }
関数が定義できたら、実行する。
gbss_2_DWIwpreproc
処理が完了すると、拡散定量値の種類に応じてstatsフォルダに「all_
以下は、標準空間上における灰白質領域のFA画像である。
gbss_3_skelpreproc
:拡散定量値をスケルトンに投影拡散定量値をスケルトンに投影するまでの手順は、次の通り。
function gbss_3_skelpreproc() { echo "Projecting MAP into skeleton..." # Skeletonize dMAP ## Create GM skeleton tbss_skeleton -i stats/mean_GM -o stats/mean_GM_skeleton ## threshold GM skeleton fslmaths stats/mean_GM_skeleton -thr ${SKELETON_THR} -bin stats/mean_GM_skeleton_mask ## Create skeleton distancemap fslmaths stats/mean_GM_mask -mul -1 -add 1 -add stats/mean_GM_skeleton_mask stats/mean_GM_skeleton_mask_dst distancemap -i stats/mean_GM_skeleton_mask_dst -o stats/mean_GM_skeleton_mask_dst for MAP in ${MAP_LIST}; do ## Project masked dMAP in MNI into skeleton tbss_skeleton -i stats/mean_GM \ -p ${SKELETON_THR} \ stats/mean_GM_skeleton_mask_dst \ ${FSLDIR}/data/standard/LowerCingulum_1mm \ stats/all_${MAP} \ stats/all_${MAP}_skeletonised done }
処理が完了すると、statsフォルダに灰白質のスケルトンmean_GM_skeleton.nii.gz
、灰白質のスケルトンの距離マップmean_GM_skeleton_mask_dst.nii.gz
、拡散定量値が投影されたスケルトンall_<MAP>_skeletonised.nii.gz
が保存される。
randomise
:スケルトンに投影された定量値画像を入力したGLMと並べ替え検定(permutation test)まず、GLMのデザインマトリックス(計画行列)とコントラストを設定する。
今回は、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)のデータがある。これらのファイルを昇順に並び替えると、先に健常者10名の画像、次に患者10名の画像の順に並ぶ。
Con0001.nii.gz
Con0002.nii.gz
Con0003.nii.gz
...
Pat0010.nii.gz
次に、GLMのデザインマトリックス(計画行列)とコントラストを決める設定ファイルを生成する。design_ttest2 <出力ファイル> <健常者数> <患者数>
でコマンドを実行。
design_ttest2 stats/design 10 10
statsフォルダに、デザインマトリックス(design.mat)とコントラスト(design.con)が生成される。
デザインマトリックス(design.mat)の中身を確認。
/Matrixの一列目は健常者データであるかどうか、二列目は患者データであるかを0, 1で表している。行の順番は、被験者ファイルの順番(昇順)に対応する。したがって、これらは対応があるようにしておかなければならない。
cat stats/design.mat
/NumWaves 2
/NumPoints 20
/PPheights 1 1
/Matrix
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
コントラスト(design.con)の中身を確認してみる。
/Matrix一列目は健常者の偏回帰係数、二列目は患者の偏回帰係数に対するもので、行は別々のコントラストである。この場合、一行目は健常者>患者の検定、二行目は健常者<患者の検定に相当する。
cat stats/design.con
/NumWaves 2
/NumContrasts 2
/PPheights 1 1
/Matrix
1 -1
-1 1
デザインマトリックス(計画行列)とコントラストの確認ができたら、randomise
コマンド使ってGLMと並べ替え検定(permutation test)を実行する。
randomise
コマンドの各オプションは、次の通り。
for MAP in ${MAP_LIST}; do randomise -i stats/all_${MAP}_skeletonised \ -o stats/gbss_${MAP} \ -m stats/mean_GM_mask \ -d stats/design.mat \ -t stats/design.con \ -n 10000 --T2 -V -x --uncorrp -R done
次のようなファイルが、生成される。
TFCEを用いた”健常群>患者群”の検定で、FWE補正をされたP値マップ(gbss_FA_tfce_corrp_tstat1.nii.gz)を確認する。ここで、得られたP値マップは1-P値のマップであることに注意する。つまり、P<.05を有意とするのであれば、P値マップで0.95-1.00の値を見ればよい。
fsleyes ${FSLDIR}/data/standard/MNI152_T1_2mm_brain \ stats/mean_GM_skeleton -cm Green \ stats/gbss_FA_tfce_p_tstat1 -cm Red-Yellow -dr 0.95 1
スケルトンは細いため有意差が見づらい場合がある。そのような時、tbss_fill
コマンドが役に立つ。
tbss_fill
コマンドの基本的な使い方は、以下の通り。
tbss_fill <P値マップ> <しきい値> <平均FA画像> <出力画像>
TFCEを用いた”健常群>患者群”の検定で、FWE補正をされたP値マップ(gbss_FA_tfce_corrp_tstat1.nii.gz)の有意差があった領域のみを0.95でしきい値処理をして抽出し、その領域を膨張させる。
tbss_fill stats/gbss_FA_tfce_p_tstat1 0.95 stats/mean_GM stats/gbss_FA_tfce_p_tstat1_fill
赤く表示されている領域は、健常群が患者群よりも有意(FWE-corrected)にFA値が大きいことを示している。
大量に定量値があり、それらをすべて検定する場合、有意差があったかどうかをすべて確認するのは大変である。そこで、各定量値画像のP値マップが0.95以上の値を持つかどうかを判定し、有意差があった場合のみ、tbss_fill
コマンドを実行する。
for PMAP in $(ls stats/ | grep tfce_corrp); do PMAX=$(fslstats stats/${PMAP} -R | cut -d " " -f2) echo ${PMAP} >>stats/tmp1.txt echo ${PMAX} >>stats/tmp2.txt if [ $(echo "${PMAX} > 0.95" | bc) == 1 ]; then tbss_fill stats/${PMAP} 0.95 stats/mean_GM stats/${PMAP}_fill fi done paste stats/tmp* >stats/tmp_corrected_P_report.txt echo -e "$(cat stats/tmp_corrected_P_report.txt)\n\n\n$(cat stats/tmp_corrected_P_report.txt | sort -r -n -k 2)" \ >stats/corrected_P_report.txt rm stats/tmp*
上のコマンドを実行すると、statsフォルダに各検定とそのP値マップの最大値が記された「corrected_P_report.txt」が出力される。
検定結果を、ファイル名でソート(上段)したものと、P値でソートしたもの(下段)に分けて保存している。
cat stats/corrected_P_report.txt
gbss_FA_tfce_corrp_tstat1.nii.gz 0.992000
gbss_FA_tfce_corrp_tstat2.nii.gz 0.416000
gbss_FW_tfce_corrp_tstat1.nii.gz 0.361839
gbss_FW_tfce_corrp_tstat2.nii.gz 0.997261
gbss_MD_tfce_corrp_tstat1.nii.gz 0.389816
gbss_MD_tfce_corrp_tstat2.nii.gz 0.985748
gbss_FW_tfce_corrp_tstat2.nii.gz 0.997261
gbss_FA_tfce_corrp_tstat1.nii.gz 0.992000
gbss_MD_tfce_corrp_tstat2.nii.gz 0.985748
gbss_FA_tfce_corrp_tstat2.nii.gz 0.416000
gbss_MD_tfce_corrp_tstat1.nii.gz 0.389816
tbss_FW_tfce_corrp_tstat1.nii.gz 0.361839
1. 目的
2. VBAとは
2.1. ファイルの準備
2.2. 脳解剖学的標準化
2.3. 平滑化
2.4. GLMと並べ替え検定(permutation test)
3. おまけ
Voxel-Based Analysis(VBA)は、脳構造解析手法の一つであり、現在では脳科学の分野において幅広く用いられている。VBAは、観測者が関心領域を設定して解析する古典的なマニュアル計測とは異なり、自動処理によって全脳をボクセルごとに解析するため、評価の客観性が高い。
VBAの解析手順は、次の通り。
VBA解析をしたファイルを準備する。
ここでは、拡散MRIから計算したFA画像を準備する。健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)いることを想定する。
. ├── Con0001_FA.nii.gz ├── Con0002_FA.nii.gz ├── Con0003_FA.nii.gz ... └── Pat0010_FA.nii.gz
非線形変換を用いて、全ての被験者のFA画像を標準FA画像(FMRIB58_FA)に位置合わせする。
以下のコマンドを実行する。
mkdir preprocessing for SUBJ in $(ls | grep _FA.nii.gz | cut -d . -f1); do # 位置合わせ ## Affine transform flirt -in ${SUBJ}.nii.gz -ref ${FSLDIR}/data/standard/FMRIB58_FA_1mm.nii.gz \ -dof 12 -omat preprocessing/${SUBJ}_indiv2std.mat ## non-linear transform fnirt --in=${SUBJ}.nii.gz --aff=preprocessing/${SUBJ}_indiv2std.mat --config=FA_2_FMRIB58_1mm.cnf \ --cout=preprocessing/${SUBJ}_warp_indiv2std.nii.gz --iout=preprocessing/${SUBJ}_instd.nii.gz done
位置合わせ方法の詳細は、こちらの記事を参考に。
https://www.nemotos.net/?p=4513
非線形変換を用いた位置合わせで、脳解剖学的標準化をしたとしても、ボクセル単位でみると完全に標準化できているわけではない。そこで、そのようなエラーを抑えるために平滑化処理をする。
まず、標準FA画像に位置合わせをした、被験者ごとのFA画像をまとめて4D-FA画像(all_FA.nii.gz)とし、4D-FA画像の平均画像からマスク画像を生成する。
mkdir stats fslmerge -t stats/all_FA.nii.gz preprocessing/*_instd.nii.gz fslmaths stats/all_FA.nii.gz -Tmean -bin mask stats/all_FA_mask.nii.gz
次に、ガウシアンフィルタを用いた平滑化処理をする。
ここでは、よく用いられるσ=3 voxel(FWHM=約7mm)のガウシアンフィルタで平滑化する。
fslmaths stats/all_FA.nii.gz -fmean -s 3 stats/all_FA_smoothed_s3.nii.gz
ガウシアンフィルタで平滑化した灰白質画像を用いて、健常群と患者群の群間比較をする。
まず、GLMのデザインマトリックス(計画行列)とコントラストを設定する。
今回は、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)のデータがある。すべての被験者のFAをls
コマンドでみると、先に健常者10名のFA、次に患者10名のFAが並んでいることが分かる。
ls |grep nii
Con0001_FA.nii.gz
Con0002_FA.nii.gz
Con0003_FA.nii.gz
...
Pat0010_FA.nii.gz
次に、GLMのデザインマトリックス(計画行列)とコントラストを決める設定ファイルを生成する。design_ttest2 <出力ファイル> <健常者数> <患者数>
でコマンドを実行。
design_ttest2 stats/design 10 10
statsフォルダに、デザインマトリックス(design.mat)とコントラスト(design.con)が生成される。
デザインマトリックス(design.mat)の中身を確認。
/Matrixの一列目は健常者データであるかどうか、二列目は患者データであるかを0, 1で表している。行の順番は、被験者フォルダにあるファイルの順番(昇順)に対応する。したがって、これらは対応があるようにしておかなければならない。
cat stats/design.mat
/NumWaves 2
/NumPoints 20
/PPheights 1 1
/Matrix
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
コントラスト(design.con)の中身を確認してみる。
/Matrix一列目は健常者の偏回帰係数、二列目は患者の偏回帰係数に対するもので、行は別々のコントラストである。この場合、一行目は健常者>患者の検定、二行目は健常者<患者の検定に相当する。
cat stats/design.con
/NumWaves 2
/NumContrasts 2
/PPheights 1 1
/Matrix
1 -1
-1 1
デザインマトリックス(計画行列)とコントラストの確認ができたら、randomise
コマンド使ってGLMと並べ替え検定(permutation test)を実行する。
randomise
コマンドの各オプションは、次の通り。
design_ttest2 stats/design 10 10 randomise -i stats/all_FA_smoothed_s3.nii.gz \ -m stats/all_FA_mask.nii.gz \ -o stats/fslvba \ -d stats/design.mat \ -t stats/design.con \ -n 10000 --T2 -V -x --uncorrp -R
次のようなファイルが、生成される。
TFCEを用いた”健常群>患者群”の検定で、FWE補正をされたP値マップ(fslvba_tfce_corrp_tstat1.nii.gz)を確認する。ここで、得られたP値マップは1-P値のマップであることに注意する。つまり、P<.05を有意とするのであれば、P値マップで0.95-1.00の値を見ればよい。
fsleyes $FSLDIR/data/standard/MNI152_T1_2mm \ stats/fslvba_tfce_corrp_tstat1 \ -cm Red-Yellow -dr 0.95 1
検定数が多い場合、有意差があったかどうかをすべて確認するのは大変である。そこで、一目で有意差があるかどうかを判断できるように、各検定ごとのP値マップの最大値を自動計測し、テキスト(corrected_P_report.txt)としてまとめる。
for PMAP in $(ls stats/ | grep tfce_corrp); do PMAX=$(fslstats stats/${PMAP} -R | cut -d " " -f2) echo ${PMAP} >>stats/tmp1.txt echo ${PMAX} >>stats/tmp2.txt done paste stats/tmp* >stats/tmp_corrected_P_report.txt echo -e "$(cat stats/tmp_corrected_P_report.txt)\n\n\n$(cat stats/tmp_corrected_P_report.txt | sort -r -n -k 2)" \ >stats/corrected_P_report.txt rm stats/tmp*
上のコマンドを実行すると、statsフォルダに各検定とそのP値マップの最大値が記された「corrected_P_report.txt」が出力される。
検定結果を、ファイル名でソート(上段)したものと、P値でソートしたもの(下段)に分けて保存している。
cat stats/corrected_P_report.txt
fslvba_tfce_corrp_tstat1.nii.gz 0.982391
fslvba_tfce_corrp_tstat2.nii.gz 0.993181
fslvba_tfce_corrp_tstat2.nii.gz 0.993181
fslvba_tfce_corrp_tstat1.nii.gz 0.982391
1. 目的
2. TBSSとは
2.1. 拡散MRIからFA画像を計算
2.2. tbss_1_preproc
:TBSSの準備
2.3. tbss_2_reg
:FAを標準空間に非線形位置合わせ
2.4. tbss_3_postreg
:平均FA画像を生成し、FAスケルトンを生成
2.5. tbss_4_prestats
:被験者ごとのFA画像を平均スケルトンに投影
2.6. tbss_non_FA
:FA画像以外の定量値をスケルトンに投影
2.7. randomise
:スケルトンに投影された定量値画像を入力したGLMと並べ替え検定(permutation test)
3. おまけ
Tract-Based Spatial Statistics (TBSS) は、白質路の統計解析をするための手法。
神経線維束の中心線(skeleton)に定量値を投影する。通常の脳画像の統計解析では、脳構造の個人差を除外するために空間的「平滑化」を用いる。しかし、平滑化の程度に原則がなく、平滑化をかけては情報があいまいになり、MRIの高空間分解能を生かせないという問題がある。一方、TBSSでは、神経線維束の中心線と思われるところにskeletonを生成し、そこに個人ごとの定量値を投影するという手法をとる。これにより、平滑化せずに群間比較をすることができるため、平滑化による問題を回避できるという利点がある。
TBSS解析では、次のような処理をする。
tbss_1_preproc
:TBSSの準備tbss_2_reg
:FAを標準空間に非線形位置合わせtbss_3_postreg
:平均FA画像を生成し、FAスケルトンを生成tbss_4_prestats
:被験者ごとのFA画像を平均スケルトンに投影tbss_non_FA
:FA画像以外の定量値をスケルトンに投影randomise
:スケルトンに投影された定量値画像を入力したGLMと並べ替え検定(permutation test)拡散MRIからFA画像を計算し、被験者ごとのFA画像を準備する。
ここでは、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)いることを想定する。
. ├── Con0001_FA.nii.gz ├── Con0002_FA.nii.gz ├── Con0003_FA.nii.gz ... └── Pat0010_FA.nii.gz
その他の定量値画像も解析に含めたい場合は、FA画像と合わせて次のように用意する。
ここでは、MDとFWという定量値画像を用意した。この時、定量値画像ごとのフォルダを作成し、各被験者の定量値画像を格納するが、ファイル名はFA画像と全く同じにすること(例. FA画像: Con0001_FA.nii.gz, MD画像: MD/Con0001_FA.nii.gz)。
. ├── Con0001_FA.nii.gz ├── Con0002_FA.nii.gz ├── Con0003_FA.nii.gz ... └── Pat0010_FA.nii.gz ├── FW │ ├── Con0001_FA.nii.gz │ ├── Con0002_FA.nii.gz │ ├── Con0003_FA.nii.gz │ ... │ └── Pat0010_FA.nii.gz └── MD ├── Con0001_FA.nii.gz ├── Con0002_FA.nii.gz ├── Con0003_FA.nii.gz ... └── Pat0010_FA.nii.gz
tbss_1_preproc
:TBSSの準備ファイルの準備ができたら、tbss_1_preproc
コマンドでTBSSの準備をする。
tbss_1_preproc *_FA.nii.gz
処理が終わると、「FAフォルダ」と「origdataフォルダ」を生成し、それぞれにFA画像が格納される。FAフォルダには、さらにFAのマスク画像が生成される。
. ├── FA │ ├── Con0001_FA.nii.gz │ ├── Con0001_FA_FA_mask.nii.gz │ ├── Con0002_FA.nii.gz │ ├── Con0002_FA_FA_mask.nii.gz │ ├── Con0003_FA.nii.gz │ ├── Con0003_FA_FA_mask.nii.gz │ ... │ └── Pat0010_FA.nii.gz │ ├── Pat0010_FA._FA_mask.nii.gz ├── FW │ ├── Con0001_FA.nii.gz │ ├── Con0002_FA.nii.gz │ ├── Con0003_FA.nii.gz │ ... │ └── Pat0010_FA.nii.gz ├── MD │ ├── Con0001_FA.nii.gz │ ├── Con0002_FA.nii.gz │ ├── Con0003_FA.nii.gz │ ... │ └── Pat0010_FA.nii.gz └── origdata ├── Con0001_FA.nii.gz ├── Con0002_FA.nii.gz ├── Con0003_FA.nii.gz ... └── Pat0010_FA.nii.gz
「FA/slicesdir/index.html」をブラウザ(e.g., Chrome)で開くことで、各被験者のFA画像一覧をみることができる。
tbss_2_reg
:FAを標準空間に非線形位置合わせtbss_1_preproc
コマンドで、全ての被験者のFA画像を1x1x1mmの標準空間に非線形的な位置合わせする。通常は、-T
オプションで標準空間にある標準FA画像(FMRIB58_FA)に位置合わせをするが、-t
オプションを用いて任意の画像に位置合わせすることもできる(推奨)。また、-n
オプションでは、被験者の中で最も位置合わせ先としてふさわしいFA画像を見つけ出し、そのFA画像にすべての被験者のFA画像を位置合わせすることができる。
ここでは、TBSSで推奨されている-T
オプションを指定しすべての被験者FA画像を標準FA画像(FMRIB58_FA)に位置合わせをする。
tbss_2_reg -T
処理が完了すると、FAフォルダに結果が保存される。
「FA/*_to_target_warp.nii.gz」が、標準FA画像に位置合わせするためのwarp fieldである。
tbss_3_postreg
:平均FA画像を生成し、FAスケルトンを生成先程生成した標準FA画像に位置合わせするためのwarp fieldを用いて、各被験者のFA画像を標準空間(MNI152)に移動させる。その後、平均FA画像を生成し、その平均FA画像からFAスケルトンを生成する。tbss_3_postreg
では、これらの処理を-S
オプションで実行することができが、代わりに標準FA画像(FMRIB58_FA mean)とそのスケルトンを用いたい場合は-T
を指定する。
ここでは、TBSSで推奨されている-S
オプションを指定する。
tbss_3_postreg -S
処理後の画像は、「statsフォルダ」に格納される。
stats/ ├── all_FA.nii.gz # 標準空間上における各被験者のFA画像 ├── mean_FA.nii.gz # 平均FA画像 ├── mean_FA_mask.nii.gz # 平均FA画像のマスク └── mean_FA_skeleton.nii.gz # 平均FA画像から生成したスケルトン画像
tbss_4_prestats
:被験者ごとのFA画像を平均スケルトンに投影tbss_4_prestats
コマンドでは、まず平均FA画像から生成したスケルトン画像(mean_FA_skeleton.nii.gz)をしきい値処理(通常 0.2)をし、スケルトンのバイナリーマスク画像を生成する。次に、このスケルトンマスクからの距離マップ(distance map)が計算され、この距離マップを参考に、被験者ごとのFA画像をスケルトン画像に格納(投影)する。
tbss_4_prestats 0.2
statsフォルダに、新たに次のファイルが生成される。
stats/ ├── all_FA_skeletonised.nii.gz # スケルトンに投影されたすべての被験者のFA画像 ├── mean_FA_skeleton_mask.nii.gz # スケルトンマスク ├── mean_FA_skeleton_mask_dst.nii.gz # スケルトンマスクからの距離マップ(distance map) └── thresh.txt # バイナリースケルトンマスク画像を作る際のしきい値
tbss_non_FA
:FA画像以外の定量値をスケルトンに投影tbss_non_FA
コマンドで、FA画像以外の定量値をスケルトンに投影する。このとき、標準空間への移動やスケルトンを生成するためのパラメータは、FA画像で使ったものが適用される。
NONFA_LIST=$(ls -F | grep / | cut -d / -f 1 | grep -v stats| grep -v origdata) for MAP in ${NONFA_LIST}; do tbss_non_FA ${MAP} done
randomise
:スケルトンに投影された定量値画像を入力したGLMと並べ替え検定(permutation test)まず、GLMのデザインマトリックス(計画行列)とコントラストを設定する。
今回は、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)のデータがある。「origdataフォルダ」をみると、先に健常者10名のFA画像、次に患者10名のFA画像が並んでいることが分かる。
ls -1 origdata
Con0001_FA.nii.gz
Con0002_FA.nii.gz
Con0003_FA.nii.gz
...
Pat0010_FA.nii.gz
次に、GLMのデザインマトリックス(計画行列)とコントラストを決める設定ファイルを生成する。design_ttest2 <出力ファイル> <健常者数> <患者数>
でコマンドを実行。
design_ttest2 stats/design 10 10
statsフォルダに、デザインマトリックス(design.mat)とコントラスト(design.con)が生成される。
デザインマトリックス(design.mat)の中身を確認。
/Matrixの一列目は健常者データであるかどうか、二列目は患者データであるかを0, 1で表している。行の順番は、origdataフォルダにあるファイルの順番(昇順)に対応する。したがって、これらは対応があるようにしておかなければならない。
cat stats/design.mat
/NumWaves 2
/NumPoints 20
/PPheights 1 1
/Matrix
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
コントラスト(design.con)の中身を確認してみる。
/Matrix一列目は健常者の偏回帰係数、二列目は患者の偏回帰係数に対するもので、行は別々のコントラストである。この場合、一行目は健常者>患者の検定、二行目は健常者<患者の検定に相当する。
cat stats/design.con
/NumWaves 2
/NumContrasts 2
/PPheights 1 1
/Matrix
1 -1
-1 1
デザインマトリックス(計画行列)とコントラストの確認ができたら、randomise
コマンド使ってGLMと並べ替え検定(permutation test)を実行する。
randomise
コマンドの各オプションは、次の通り。
for MAP in FA ${NONFA_LIST}; do randomise -i stats/all_${MAP}_skeletonised \ -o stats/tbss_${MAP} \ -m stats/mean_FA_skeleton_mask \ -d stats/design.mat \ -t stats/design.con \ -n 10000 --T2 -V -x --uncorrp -R done
次のようなファイルが、生成される。
TFCEを用いた”健常群>患者群”の検定で、FWE補正をされたP値マップ(tbss_FA_tfce_corrp_tstat1.nii.gz)を確認する。ここで、得られたP値マップは1-P値のマップであることに注意する。つまり、P<.05を有意とするのであれば、P値マップで0.95-1.00の値を見ればよい。
fsleyes ${FSLDIR}/data/standard/FMRIB58_FA_1mm.nii.gz \ ${FSLDIR}/data/standard/FMRIB58_FA-skeleton_1mm.nii.gz -cm Green \ stats/tbss_FA_tfce_p_tstat1.nii.gz -cm Red-Yellow -dr 0.95 1
スケルトンは細いため有意差が見づらい場合がある。そのような時、tbss_fill
コマンドが役に立つ。
tbss_fill
コマンドの基本的な使い方は、以下の通り。
tbss_fill <P値マップ> <しきい値> <平均FA画像> <出力画像>
TFCEを用いた”健常群>患者群”の検定で、FWE補正をされたP値マップ(tbss_FA_tfce_corrp_tstat1.nii.gz)の有意差があった領域のみを0.95でしきい値処理をして抽出し、その領域を膨張させる。
tbss_fill stats/tbss_FA_tfce_p_tstat1.nii.gz 0.95 stats/mean_FA stats/tbss_FA_tfce_p_tstat1_fill.nii.gz
赤く表示されている領域は、健常群が患者群よりも有意(FWE-corrected)にFA値が大きいことを示している。
大量に定量値があり、それらをすべて検定する場合、有意差があったかどうかをすべて確認するのは大変である。そこで、各定量値画像のP値マップが0.95以上の値を持つかどうかを判定し、有意差があった場合のみ、tbss_fill
コマンドを実行する。
for PMAP in $(ls stats/ | grep tfce_corrp); do PMAX=$(fslstats stats/${PMAP} -R | cut -d " " -f2) echo ${PMAP} >>stats/tmp1.txt echo ${PMAX} >>stats/tmp2.txt if [ $(echo "${PMAX} > 0.95" | bc) == 1 ]; then tbss_fill stats/${PMAP} 0.95 stats/mean_FA stats/${PMAP}_fill fi done paste stats/tmp* >stats/tmp_corrected_P_report.txt echo -e "$(cat stats/tmp_corrected_P_report.txt)\n\n\n$(cat stats/tmp_corrected_P_report.txt | sort -r -n -k 2)" \ >stats/corrected_P_report.txt rm stats/tmp*
上のコマンドを実行すると、statsフォルダに各検定とそのP値マップの最大値が記された「corrected_P_report.txt」が出力される。
検定結果を、ファイル名でソート(上段)したものと、P値でソートしたもの(下段)に分けて保存している。
cat stats/corrected_P_report.txt
tbss_FA_tfce_corrp_tstat1.nii.gz 0.992000
tbss_FA_tfce_corrp_tstat2.nii.gz 0.416000
tbss_FW_tfce_corrp_tstat1.nii.gz 0.361839
tbss_FW_tfce_corrp_tstat2.nii.gz 0.997261
tbss_MD_tfce_corrp_tstat1.nii.gz 0.389816
tbss_MD_tfce_corrp_tstat2.nii.gz 0.985748
tbss_FW_tfce_corrp_tstat2.nii.gz 0.997261
tbss_FA_tfce_corrp_tstat1.nii.gz 0.992000
tbss_MD_tfce_corrp_tstat2.nii.gz 0.985748
tbss_FA_tfce_corrp_tstat2.nii.gz 0.416000
tbss_MD_tfce_corrp_tstat1.nii.gz 0.389816
tbss_FW_tfce_corrp_tstat1.nii.gz 0.361839
1. 目的
2. VBMとは
2.1. 3D-T1WIの準備
2.2. fslvbm_1_bet
: 3D-T1WIの脳頭蓋除去
2.3. fslvbm_2_template
: 被験者の脳から灰白質テンプレート画像を生成
2.4. fslvbm_3_proc
: すべての被験者の灰白質を灰白質テンプレート画像に位置合わせし、ボクセルをモジュレーション後、平滑化
2.5. randomise
:GLMと並べ替え検定(permutation test)
3. おまけ
Voxel-Based Morphometry(VBM)は、脳構造解析手法の一つであり、特に脳容積を対象に解析する。
VBMは、古典的なマニュアルの脳容積計測とは異なり、自動処理によって全脳を客観的に評価することができ、現在では脳科学の分野において幅広く用いられている。
VBM解析では、次のような処理をする。
fslvbm_1_bet
: 3D-T1WIの脳頭蓋除去fslvbm_2_template
: 被験者の脳から灰白質テンプレート画像を生成fslvbm_3_proc
: すべての被験者の灰白質を灰白質テンプレート画像に位置合わせし、ボクセルをモジュレーション後、平滑化randomise
:GLMと並べ替え検定(permutation test)各被験者の3D-T1WIを準備する。
ここでは、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)いることを想定する。
. ├── Con0001_T1w.nii.gz ├── Con0002_T1w.nii.gz ├── Con0003_T1w.nii.gz ... └── Pat0010_T1w.nii.gz
fslvbm_1_bet
: 3D-T1WIの脳頭蓋除去ファイルの準備ができたら、fslvbm_1_bet
コマンドを実行して脳頭蓋除去をする。3D-T1WIに首を含まない場合、-b
オプションを指定し、首を含む場合は-N
オプションを指定する。
ここでは、3D-T1WIに首を含んでいたとして、-N
オプションを指定する。
fslvbm_1_bet -N # 首を含む場合 # fslvbm_1_bet -b # 首を含まない場合
処理が完了すると「strucフォルダ」が生成され、その中に頭蓋除去後の画像が保存される。頭蓋除去後の画像は「*_brain.nii.gz」というファイル名で保存される。頭蓋除去前後の画像一覧をみるには、「struc/slicesdir/index.html」を開く。下地が頭蓋除去前の画像であり、赤い線が頭蓋除去後の画像である。
fslvbm_2_template
: 被験者の脳から灰白質テンプレート画像を生成次に、脳の解剖学的標準化で用いるターゲット画像(標準脳)を、頭蓋除去済みの被験者脳から生成する。
まず最初に、標準脳を生成するために用いる被験者の脳画像リスト(template_list)を作成する。このとき、バイアスがかからないように健常者と患者の人数が同じになるようにリストを作る。
ここでは、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)の合計20名を「template_list」に記載した。
cat template_list
Con0001_T1w.nii.gz
Con0001_T1w.nii.gz
Con0001_T1w.nii.gz
...
Pat0010_T1w.nii.gz
template_listは、3D-T1WI(*_T1w.nii.gz)と同じ階層に配置する。
. ├── template_list # ここに配置 ├── Con0001_T1w.nii.gz ├── Con0002_T1w.nii.gz ├── Con0003_T1w.nii.gz ... └── Pat0010_T1w.nii.gz
標準脳を生成するために用いる被験者の脳画像リスト(template_list)が準備できたら、頭蓋除去済みの脳画像を、灰白質・白質・脳脊髄液にセグメントする。次に、全ての被験者の灰白質(*_struc_GM.nii.gz)を元々ある灰白質標準脳(ICBM-152)にアフィン変換で位置合わせし、平均化する。このようにして得られた平均灰白質画像を、x軸に対して反転し再度平均化して、左右対称の平均灰白質画像(template_GM_init.nii.gz)を生成する。既に作成した「template_list」に記載されている灰白質画像を、先ほど作成した平均灰白質画像に位置合わせをし、平均化、さらに左右反転後に再度平均化する。最後に、元々ある灰白質標準脳(ICBM-152)に非線形位置合わせをし、2x2x2 mm^3の灰白質テンプレートを標準空間に生成する。
fslvbm_2_template
コマンドにはオプションがあり、各被験者の灰白質画像を元々ある灰白質標準脳(ICBM-152)に位置合わせする際に、アフィン変換を用いる場合は-a
オプションを、非線形変換を用いる場合は-n
オプションを指定する。
ここでは、fslvbm_2_template
に-a
オプションを指定して実行する。
fslvbm_2_template -a # Affine registration # fslvbm_2_template -n # non-linear registration
処理が完了するとstrucフォルダに灰白質テンプレート(template_GM_4D.nii.gz)が生成される。
fslvbm_3_proc
: すべての被験者の灰白質を灰白質テンプレート画像に位置合わせし、ボクセルをモジュレーション後、平滑化非線形変換を用いた位置合わせで、各被験者の灰白質画像を作成した灰白質テンプレート(template_GM_4D.nii.gz)に合わせ、statsフォルダに全被験者の灰白質4D画像を作る(GM_merg.nii.gz)。この時、被験者ごとの灰白質容積の違いを反映できるよう、非線形変換の際に収縮・拡大された度合いに応じて、灰白質画像のボクセル値を補正(モジュレーション)する。実際には、warp fieldのヤコビアンを灰白質画像の各ボクセルにかけ合わせる。その後、モジュレーションされた画像(GM_mod_merg.nii.gz)はσ=2, 3, 4mm (およそFWHM=2×2.3=4.6mmからFWHM=9mm)のガウシアンフィルタで平滑化する(GM_mod_merg_s[2,3,4].nii.gz)。
fslvbm_3_proc
以下は、モジュレーションした灰白質画像(GM_mod_merg.nii.gz)とσ=3mmのガウシアンフィルタで平滑化した灰白質画像(GM_mod_merg_s3.nii.gz)である。
randomise
:GLMと並べ替え検定(permutation test)ガウシアンフィルタで平滑化した灰白質画像を用いて、健常群と患者群の群間比較をする。
まず、GLMのデザインマトリックス(計画行列)とコントラストを設定する。
今回は、健常者10名(ID: Con0001~Con0010)と患者10名(ID: Pat0001~Pat0010)のデータがある。すべての被験者の3D-T1WIをls
コマンドでみると、先に健常者10名の3D-T1WI、次に患者10名の3D-T1WIが並んでいることが分かる。
ls |grep nii
Con0001_T1w.nii.gz
Con0002_T1w.nii.gz
Con0003_T1w.nii.gz
...
Pat0010_T1w.nii.gz
次に、GLMのデザインマトリックス(計画行列)とコントラストを決める設定ファイルを生成する。design_ttest2 <出力ファイル> <健常者数> <患者数>
でコマンドを実行。
design_ttest2 stats/design 10 10
statsフォルダに、デザインマトリックス(design.mat)とコントラスト(design.con)が生成される。
デザインマトリックス(design.mat)の中身を確認。
/Matrixの一列目は健常者データであるかどうか、二列目は患者データであるかを0, 1で表している。行の順番は、被験者ファイルの順番(昇順)に対応する。したがって、これらは対応があるようにしておかなければならない。
cat stats/design.mat
/NumWaves 2
/NumPoints 20
/PPheights 1 1
/Matrix
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
1 0
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
0 1
コントラスト(design.con)の中身を確認してみる。
/Matrix一列目は健常者の偏回帰係数、二列目は患者の偏回帰係数に対するもので、行は別々のコントラストである。この場合、一行目は健常者>患者の検定、二行目は健常者<患者の検定に相当する。
cat stats/design.con
/NumWaves 2
/NumContrasts 2
/PPheights 1 1
/Matrix
1 -1
-1 1
デザインマトリックス(計画行列)とコントラストの確認ができたら、randomise
コマンド使ってGLMと並べ替え検定(permutation test)を実行する。
randomise
コマンドの各オプションは、次の通り。
fsleyes $FSLDIR/data/standard/MNI152_T1_2mm \ stats/fslvbm_tfce_corrp_tstat1 \ -cm Red-Yellow -dr 0.95,1
次のようなファイルが、生成される。
TFCEを用いた”健常群>患者群”の検定で、FWE補正をされたP値マップ(fslvbm_tfce_corrp_tstat1.nii.gz)を確認する。ここで、得られたP値マップは1-P値のマップであることに注意する。つまり、P<.05を有意とするのであれば、P値マップで0.95-1.00の値を見ればよい。
fsleyes $FSLDIR/data/standard/MNI152_T1_2mm \ stats/fslvbm_tfce_corrp_tstat1 -cm Red-Yellow
検定数が多い場合、有意差があったかどうかをすべて確認するのは大変である。そこで、一目で有意差があるかどうかを判断できるように、各検定ごとのP値マップの最大値を自動計測し、テキスト(corrected_P_report.txt)としてまとめる。
for PMAP in $(ls stats/ | grep tfce_corrp); do PMAX=$(fslstats stats/${PMAP} -R | cut -d " " -f2) echo ${PMAP} >>stats/tmp1.txt echo ${PMAX} >>stats/tmp2.txt done paste stats/tmp* >stats/tmp_corrected_P_report.txt echo -e "$(cat stats/tmp_corrected_P_report.txt)\n\n\n$(cat stats/tmp_corrected_P_report.txt | sort -r -n -k 2)" \ >stats/corrected_P_report.txt rm stats/tmp*
上のコマンドを実行すると、statsフォルダに各検定とそのP値マップの最大値が記された「corrected_P_report.txt」が出力される。
検定結果を、ファイル名でソート(上段)したものと、P値でソートしたもの(下段)に分けて保存している。
cat stats/corrected_P_report.txt
fslvbm_tfce_corrp_tstat1.nii.gz 0.982391
fslvbm_tfce_corrp_tstat2.nii.gz 0.993181
fslvbm_tfce_corrp_tstat2.nii.gz 0.993181
fslvbm_tfce_corrp_tstat1.nii.gz 0.982391
FSLを用いて、拡散テンソルイメージング(DTI)をするには、dtifit
コマンドを用いる。
dtifit
のヘルプは、次の通り。
Usage: dtifit -k <filename> dtifit --verbose Compulsory arguments (You MUST set one or more of): -k,--data dti data file -o,--out Output basename -m,--mask Bet binary mask file -r,--bvecs b vectors file -b,--bvals b values file Optional arguments (You may optionally specify one or more of): -V,--verbose switch on diagnostic messages -h,--help display this message --cni Input confound regressors --sse Output sum of squared errors -w,--wls Fit the tensor with weighted least squares --kurt Output mean kurtosis map (for multi-shell data) --kurtdir Output parallel/perpendicular kurtosis maps (for multi-shell data) --littlebit Only process small area of brain --save_tensor Save the elements of the tensor -z,--zmin min z -Z,--zmax max z -y,--ymin min y -Y,--ymax max y -x,--xmin min x -X,--xmax max x --gradnonlin Gradient Nonlinearity Tensor file
基本的な使い方は、以下の通り。
dtifit -k <DWI images> -o <OUTPUT> -m <MASK> -r <b-vectors file> -b <b-values file>
まず、次のファイルを用意する。
. ├── DWI.nii.gz # 拡散MRI ├── DWI_mask.nii.gz # 拡散MRIのマスク画像 ├── bvals # b-values └── bvecs # b-vectors
ファイルの用意ができたら、dtifit
を次のように実行する
dtifit -k DWI.nii.gz -o output -m DWI_mask.nii.gz -r bvecs -b bvals
処理が終わると次のファイルが生成される。
FA, MD, L1の画像は、以下。
1. 目的
2. 脳画像解析における拡散MRIの問題点
2.1. 磁化率効果に起因する画像の幾何学的歪み(位相エンコードをAP方向として撮像した場合)
2.2. 撮像中の頭の動き
3. 使用するソフトウェア
3.1. TOPUP
3.1.1. TOPUP
3.2. –indexについて
3.2.1. TOPUP
3.3. EDDY
3.4. 頭の動き補正の効果
4. サンプルコード
1. 目的
2. コマンド
3. 使用例
3.1. FSLアルゴリズムを用いる場合
3.2. FreeSurferアルゴリズムを用いる場合
4. 結果
MRtrixの5ttgen
を用いて、次の5つの組織(five-tissue-type: 5TT)画像を生成する。
5ttgen
のヘルプは、次の通り。
SYNOPSIS Generate a 5TT image suitable for ACT USAGE 5ttgen [ options ] algorithm ... algorithm Select the algorithm to be used to complete the script operation; additional details and options become available once an algorithm is nominated. Options are: freesurfer, fsl, gif, hsvs DESCRIPTION 5ttgen acts as a 'master' script for generating a five-tissue-type (5TT) segmented tissue image suitable for use in Anatomically-Constrained Tractography (ACT). A range of different algorithms are available for completing this task. When using this script, the name of the algorithm to be used must appear as the first argument on the command-line after '5ttgen'. The subsequent compulsory arguments and options available depend on the particular algorithm being invoked. Each algorithm available also has its own help page, including necessary references; e.g. to see the help page of the 'fsl' algorithm, type '5ttgen fsl'. Options common to all 5ttgen algorithms -nocrop Do NOT crop the resulting 5TT image to reduce its size (keep the same dimensions as the input image) -sgm_amyg_hipp Represent the amygdalae and hippocampi as sub-cortical grey matter in the 5TT image Additional standard options for Python scripts -nocleanup do not delete intermediate files during script execution, and do not delete scratch directory at script completion. -scratch /path/to/scratch/ manually specify the path in which to generate the scratch directory. -continue <ScratchDir> <LastFile> continue the script from a previous execution; must provide the scratch directory path, and the name of the last successfully-generated file. Standard options -info display information messages. -quiet do not display information messages or progress status. Alternatively, this can be achieved by setting the MRTRIX_QUIET environment variable to a non- empty string. -debug display debugging messages. -force force overwrite of output files. -nthreads number use this number of threads in multi-threaded applications (set to 0 to disable multi-threading). -config key value (multiple uses permitted) temporarily set the value of an MRtrix config file entry. -help display this information page and exit. -version display version information and exit.
基本的な使い方は、以下の通り。5ttgen
のアルゴリズムは、freesurfer, fsl, gif, hsvsがあるが、ここではfreesurferとfslのアルゴリズムについて使い方を解説する。
5ttgen [アルゴリズム] <入力画像> <出力画像>
FSLアルゴリズムを用いる場合、3D-T1WI(T1w.nii.gz)が必要となる。また、オプションとして3D-T2WIも入力することができる。
5ttgen fsl T1w.nii.gz 5tt.nii.gz
FreeSurferアルゴリズムを用いる場合、Freesurferの生成ファイルであるaparc+aseg.mgz(asegとついたファイル)が必要となる。
FreeSurferの使い方は、こちらの記事を参考にするとよい。
aparc+aseg.mgzが準備できたら、以下のコマンドを実行する。
5ttgen freesurfer aparc+aseg.mgz 5tt.nii.gz
5ttgen
で生成された画像は、5ボリュームデータであり、各ボリュームと対応する組織は次の通り。
以下に、FSLとFreeSurferのアルゴリズムを用いて5ttgen
した結果(緑)を示す。
1. 目的
2. Boundary-Based Registration (BBR)とは
3. コマンド
4. 目的使用例
4.1. 目的前準備
4.2. 目的Boundary-Based Registration: BBB
Doug Greve氏によって開発されたEPI画像(拡散MRIや機能的MRI)用の位置合わせツールで、EPI画像と構造MRI画像(例:T1WI)との位置合わせで、灰白質・白質境界(Boundary)を頼りに位置合わせをする。以下のウェブサイトに、詳細な説明と分かりやすい図がある。
Boundary-Based Registration (BBR) に基づいた構造MRIと拡散MRIの位置合わせをFSLで実装するには、epi_reg
を用いる。
epi_reg
のヘルプは、以下の通り。
Usage: epi_reg [options] --epi=<EPI image> --t1=<wholehead T1 image> --t1brain=<brain extracted T1 image> --out=<output name> Optional arguments --fmap=<image> : fieldmap image (in rad/s) --fmapmag=<image> : fieldmap magnitude image - wholehead extracted --fmapmagbrain=<image> : fieldmap magnitude image - brain extracted --gdc=<image> : Gradient-distortion corection warpfield --wmseg=<image> : white matter segmentation of T1 image --echospacing=<val> : Effective EPI echo spacing (sometimes called dwell time) - in seconds --pedir=<dir> : phase encoding direction, dir = x/y/z/-x/-y/-z --weight=<image> : weighting image (in T1 space) --nofmapreg : do not perform registration of fmap to T1 (use if fmap already registered) --noclean : do not clean up intermediate files -v : verbose output -h : display this help message e.g.: epi_reg --epi=example_func --t1=struct --t1brain=struct_brain --out=epi2struct --fmap=fmap_rads --fmapmag=fmap_mag --fmapmagbrain=fmap_mag_brain --echospacing=0.0005 --pedir=-y Note that if parallel acceleration is used in the EPI acquisition then the *effective* echo spacing is the actual echo spacing between acquired lines in k-space divided by the acceleration factor.
基本的な使い方は、以下の通り。
epi_reg --epi=<b=0画像(spin-echo EPI)> --t1=<頭蓋除去前の3D-T1WI> --t1brain=<頭蓋除去後の3D-T1WI> --echospacing=<echo spacing時間(sec)> --pedir=<位相エンコード方向> --out=<出力画像>
次のファイルを準備する。
. ├── CSF_GM_WM_seg.nii.gz:FASTで作成したCSF, GM, WMのセグメント ├── DWI.nii.gz:拡散MRI画像(b≠0) ├── DWI_b0.nii.gz:拡散MRIのb=0画像(spin-echo EPI) ├── T1_skull_stripped.nii.gz:頭蓋除去後の3D-T1WI └── T1w.nii.gz:頭蓋除去前の3D-T1WI
頭蓋除去とFASTを用いたCSF, GM, WMのセグメントは、以下の記事を参考にするとよい。
また、拡散MRIからb値ごとに画像を抽出する方法は、以下の記事を参考にするとよい。
ここで、構造MRI(3D-T1WI)空間にあるCSF, GM, WMのセグメント(CSF_GM_WM_seg.nii.gz)を、拡散MRIの空間に移動させることを考える。
そのために、まずBBBを使って拡散MRIのb=0画像(spin-echo EPI)を、構造MRI(3D-T1WI)に位置合わせする。
この時、echo spacing時間( sec)--echospacing
と位相エンコード方向--pedir
を前もって調べておく必要がある。これらの情報は、dcm2niix
を用いてDICOM形式からNIfTI形式に変換する際に出力されるJSONファイルに記載されている。dcm2niix
の使い方は、以下の記事を参考にするとよい。
準備ができたら、以下のようにepi_reg
コマンドを実行して、拡散MRIから構造MRIへ位置合わせする変換行列を生成する。
epi_reg --epi=DWI_b0.nii.gz --t1=T1w.nii.gz --t1brain=T1_skull_stripped.nii.gz --echospacing=0.0380544 --pedir=-y --out=DWI2T1w
次に、先ほどの変換行列(拡散MRI空間→構造MRI空間)の逆変換行列(構造MRI空間→拡散MRI空間)を生成する。
convert_xfm -omat T1w2DWI.mat -inverse DWI2T1w.mat
逆変換行列(構造MRI空間→拡散MRI空間)を構造MRI空間にいるCSF, GM, WMのセグメント(CSF_GM_WM_seg.nii.gz)に適用することで、セグメントを拡散MRI空間に移動させる。このとき、-applyxfm
と-interp nearestneighbour
としていることに注意する。FLIRTの使い方の詳細はこちら。
flirt -in CSF_GM_WM_seg.nii.gz -ref DWI_b0.nii.gz -out CSF_GM_WM_seg_DWIspace.nii.gz -init T1w2DWI.mat -applyxfm -interp nearestneighbour
処理前後のセグメントの様子をみてみる。
fsleyes DWI.nii.gz \ CSF_GM_WM_seg.nii.gz -cm random \ CSF_GM_WM_seg_DWIspace.nii.gz -cm random
軸位断で後頭葉付近を拡大してみる。
私のメインマシンは Lin4Neuro 18.04 ですが、そろそろ Lin4Neuro 20.04 への移行を考えています。
今、実験機には NVIDIA GeForce RTX 2070 が備え付けられています。
これを使って、FSL 6.0.5 の eddy をGPUが使えるように設定し、なおかつ、Tensorflow, Pytorch といった Deep Learning のフレームワークも使えるようにしたいと思います。
FSL 6.0.5 にはデフォルトで CUDA 10.2 に対応した eddy_cuda10.2 が配布されています。なので、CUDA 10.2を入れることにします。
なお、これは Ubuntu 18.04 でも全く問題なくできることがわかりましたので、タイトルを変更しました。