FreeSurfer
脳MRI構造画像の自動解析手法として広く用いられているものは、大きく分けてvoxel-basedなもの(VBM)とsurface-basedなもの(SBM)に分けられる。
VBMでは構造画像に線形・非線形の変形を施してテンプレート脳と合致させて、その位置合わせ後の各画像についてvoxelごとに比較を行う。VBMのためのツールにはSPM12やFSLなどがある。
一方、SBMでは、構造画像のsegmentation結果から脳表面および皮質・白質の境界面を3次元的に再構成し、この3次元モデルをテンプレート脳と位置合わせし、位置合わせ後の3次元モデルの頂点 ごとに比較を行う。SPMのためのツールにはFreeSurferがある。ここでは、FreeSurferによる構造画像の解析についてメモしておく。
概要
FreeSurferはLinuxとmacOSに対応する脳画像解析パッケージで、MGH関連の研究機関the Athinoula A. Martinos Center for Biomedical ImagingのLaboratory for Computational Neuroimagingで開発されており、無料で公開されているオープンソースソフトウェアである。
小さなコマンド群とそれらを組み合わせるスクリプトやビューワーで構成されている。
構造画像に対しては
- 頭蓋除去、B1バイアスフィールド補正、皮質・白質のセグメンテーション
- 皮質の三次元モデル再構成
- 皮質および皮質下の領域ラベリング
- 皮質モデルのアトラスへの非線形registration
- 形態の違いについての統計的な解析
を行うことができる。
セットアップ
FreeSurferでは環境変数$SUBJECTS_DIR以下に被験者IDをディレクトリ名とするディレクトリを作成して、その中に前処理後の結果を保存するようになっている。
デフォルトでは$SUBJECTS_DIRは$FREESURFER_HOME/subjectsになっているが、プロジェクトごとに変更するのが良いだろう。
前処理
解析を始めるには、まずは構造画像に対して前処理を行う必要がある。この過程は全自動だが、1例あたり5時間くらいかかる。
前処理は、かなり複雑なパイプラインから構成されているが、recon-allスクリプトを使って簡単に実行できる。
recon-all \
-all \ # すべての前処理過程を行う
-subjid SUBJ \ # 被験者ID
-i /path/to/your/dataset/sub-SUBJ/anat/sub-SUBJ_T1w.nii \ # 入力する構造画像
-sd /path/to/your/subject/directory \ # 前処理後の結果を保存するSUBJECTS_DIR
-parallel -openmp NUM_CORES
CPUがマルチコアの場合は-parallelオプションを使うことで一部の処理は並列で実行されるようになるが、手元の環境では8スレッドで実行して前処理にかかる時間は半分程度にしかならなかった。そのため、このオプションは単一の画像を処理する際には有効であるものの、多数の画像をまとめて前処理する場合はparallelコマンドなどを使ってrecon-allそのものを並列で実行した方が効率的である。
パイプラインの中身
パイプラインは体積画像を扱うvolume-based streamと皮質の3次元モデルを扱うsurface-based streamとに大別される。
大雑把に言うと、volume-based streamで皮質下領域のセグメンテーションと白質表面の抽出までを行う。この白質の体積データの構成が終わるとsurface-based streamに移行し、左右の大脳半球それぞれに対して独立して処理がなされる。まず白質の体積表現から表面の3次元モデルを構成し、それをT1画像をもとに脳表に達するまで膨張させて皮質表面のモデルも構成する。さらにその表面を引き伸ばして、球体まで変形し、球面上で曲率のマップからテンプレートとの位置合わせを行う。
recon-allではこれらの前処理過程をいくつかのstepに分けており、それらを独立して実行することができる。多少細かくパイプラインの各課程について述べておく：
- Step 1
- mgzに変換
- （ 入力画像が複数ある場合は）平均化および動き補正
- 1mmのisometric voxelにconform
- MNI305に合わせるaffine transformの計算（※変形を画像にapplyするわけではない）
- 信号強度の不均一性の補正・標準化
- 頭蓋除去
- Step 2
- Gaussian Classifier Atlasへの位置合わせ
- 頸部除去
- 皮質下構造のセグメンテーション・ラベリング
- 白質の抽出、内部の充填
- ーーー（ここでvolume-based streamからsurface-based streamに移行）ーーー
- 白質表面のテッセレーション、トポロジーの修正
- 皮質表面モデルの構成
- 膨張モデルの構成
- Step 3
- 球面モデルへの変形
- 曲率を用いたatlasへの位置合わせ
- 皮質のパーセレーション・ラベリング
- ラベルの体積画像への出力
なお、recon-all -helpでより細かく各ステップを確認できる。
前処理の結果できるファイルたち
前述した通りFreeSurferは小さなexecutableたちをshell scriptでorchestrateするような構成になっており、各コマンドの処理の結果作成された中間ファイルがすべて残される。そのため、25 MB程度(だいたい16 bits = 2 bytes / voxel x 256 x 256 x 200)の構造画像から派生した前処理後のファイルたちは全部で300 MB程度にまで膨れ上がる。
処理後のディレクトリ構成は以下の通り：
- /mri : 体積画像
- /orig
- 001.mgz : 入力された構造画像（をmgzに変換したもの）
- rawavg.mgz : 複数の入力画像がある場合、それらを動きを補正し平均したもの
- orig.mgz : rawavg.mgzを1mm等方ボクセルにresampleしたもの（これ以降の体積画像は1mm等方ボクセルになっている）
- nu.mgz : N3により不均一性を補正したもの
- T1.mgz : nu.mgzにさらに信号値の標準化を施したもの
- brainmask.mgz : 頭蓋除去後の画像
- aseg.mgz : 皮質下領域のセグメンテーション結果
- wm.mgz : 白質を抽出したもの
- filled.mgz : wm.mgzの中を埋めたもの
- aparc.(a2009s|DKTatlas)+aseg.mgz : aseg.mgzに皮質パーセレーションの情報を加えたもの
- /orig
- /label
- /scripts : 実行時のログ
- recon-all.log
- recon-all.done
- recon-all.cmd
- recon-all-status.log
- /stats : 統計情報
- aseg.stats
- brainvol.stats
- wmparc.stats
- ?h.curv.stats
- ?h.aparc.stats
- ?h.aparc.pial.stats
- /surf : 皮質の3次元モデルおよび頂点ごとの形態データ（頂点は各メッシュですべて対応している）
- ?h.pial : 脳表のメッシュ
- ?h.white : 白質表面のメッシュ
- ?h.inflated : 脳溝の内側も見えるように引き伸ばされたメッシュ
- ?h.sphere : 球面に変形されたメッシュ
- ?h.sphere.reg : 曲率による位置合わせ後のメッシュ
- ?h.area : 各頂点ごとの面積
- ?h.area.pial
- ?h.curv : 各頂点ごとの曲率
- ?h.thickness : 各頂点ごとの厚さ
- /tmp
- /touch
- /trash
画像ビューワ
FreeViewというコマンドで、FreeSurferの独自形式でフォーマットされた体積画像やメッシュを表示することができる。
Quality Control
ROI-based analysis
FreeSurfer wikiのAnatomical ROIページも参照のこと。
皮質下領域の自動セグメンテーションや皮質の自動パーセレーションの結果に基づいて、各領域ごとに値を集計したものがstatsディレクトリ以下に保存される。これはテキストファイルで、例えばaseg.statsの中身は以下のようである：
$ cat /path/to/your/dataset/derivatives/freesurfer/SUBJID/stats/aseg.stats
# Title Segmentation Statistics
#
# generating_program mri_segstats
# cvs_version 7.1.1
# cmdline mri_segstats --seed 1234 --seg mri/aseg.mgz --sum stats/aseg.stats --pv mri/norm.mgz --empty --brainmask mri/brainmask.mgz --brain-vol-from-seg --excludeid 0 --excl-ctxgmwm --supratent --subcortgray --in mri/norm.mgz --in-intensity-name norm --in-intensity-units MR --etiv --surf-wm-vol --surf-ctx-vol --totalgray --euler --ctab /home/hiro/freesurfer/ASegStatsLUT.txt --subject SUBJID
# sysname Linux
# hostname HP-Z8-G4
# machine x86_64
# user hiro
# anatomy_type volume
#
# SUBJECTS_DIR /path/to/your/dataset/derivatives/freesurfer
# subjectname SUBJID
# BrainVolStatsFixed-NotNeeded because voxelvolume=1mm3
# Measure BrainSeg, BrainSegVol, Brain Segmentation Volume, 1044146.000000, mm^3
# Measure BrainSegNotVent, BrainSegVolNotVent, Brain Segmentation Volume Without Ventricles, 1026548.000000, mm^3
# Measure VentricleChoroidVol, VentricleChoroidVol, Volume of ventricles and choroid plexus, 14437.000000, mm^3
# Measure lhCortex, lhCortexVol, Left hemisphere cortical gray matter volume, 240124.600036, mm^3
# Measure rhCortex, rhCortexVol, Right hemisphere cortical gray matter volume, 241704.181495, mm^3
# Measure Cortex, CortexVol, Total cortical gray matter volume, 481828.781531, mm^3
# Measure lhCerebralWhiteMatter, lhCerebralWhiteMatterVol, Left hemisphere cerebral white matter volume, 181579.500000, mm^3
# Measure rhCerebralWhiteMatter, rhCerebralWhiteMatterVol, Right hemisphere cerebral white matter volume, 182983.500000, mm^3
# Measure CerebralWhiteMatter, CerebralWhiteMatterVol, Total cerebral white matter volume, 364563.000000, mm^3
# Measure SubCortGray, SubCortGrayVol, Subcortical gray matter volume, 52019.000000, mm^3
# Measure TotalGray, TotalGrayVol, Total gray matter volume, 634930.781531, mm^3
# Measure SupraTentorial, SupraTentorialVol, Supratentorial volume, 914981.000000, mm^3
# Measure SupraTentorialNotVent, SupraTentorialVolNotVent, Supratentorial volume, 897383.000000, mm^3
# Measure Mask, MaskVol, Mask Volume, 1438526.000000, mm^3
# Measure BrainSegVol-to-eTIV, BrainSegVol-to-eTIV, Ratio of BrainSegVol to eTIV, 0.849547, unitless
# Measure MaskVol-to-eTIV, MaskVol-to-eTIV, Ratio of MaskVol to eTIV, 1.170426, unitless
# Measure lhSurfaceHoles, lhSurfaceHoles, Number of defect holes in lh surfaces prior to fixing, 21, unitless
# Measure rhSurfaceHoles, rhSurfaceHoles, Number of defect holes in rh surfaces prior to fixing, 12, unitless
# Measure SurfaceHoles, SurfaceHoles, Total number of defect holes in surfaces prior to fixing, 33, unitless
# Measure EstimatedTotalIntraCranialVol, eTIV, Estimated Total Intracranial Volume, 1229062.307009, mm^3
# SegVolFile mri/aseg.mgz
# SegVolFileTimeStamp 2020/10/30 02:11:11
# ColorTable /home/hiro/freesurfer/ASegStatsLUT.txt
# ColorTableTimeStamp 2020/07/24 00:51:51
# InVolFile mri/norm.mgz
# InVolFileTimeStamp 2020/10/29 20:04:17
# InVolFrame 0
# PVVolFile mri/norm.mgz
# PVVolFileTimeStamp 2020/10/29 20:04:17
# Excluding Cortical Gray and White Matter
# ExcludeSegId 0 2 3 41 42
# VoxelVolume_mm3 1
# TableCol 1 ColHeader Index
# TableCol 1 FieldName Index
# TableCol 1 Units NA
# TableCol 2 ColHeader SegId
# TableCol 2 FieldName Segmentation Id
# TableCol 2 Units NA
# TableCol 3 ColHeader NVoxels
# TableCol 3 FieldName Number of Voxels
# TableCol 3 Units unitless
# TableCol 4 ColHeader Volume_mm3
# TableCol 4 FieldName Volume
# TableCol 4 Units mm^3
# TableCol 5 ColHeader StructName
# TableCol 5 FieldName Structure Name
# TableCol 5 Units NA
# TableCol 6 ColHeader normMean
# TableCol 6 FieldName Intensity normMean
# TableCol 6 Units MR
# TableCol 7 ColHeader normStdDev
# TableCol 7 FieldName Itensity normStdDev
# TableCol 7 Units MR
# TableCol 8 ColHeader normMin
# TableCol 8 FieldName Intensity normMin
# TableCol 8 Units MR
# TableCol 9 ColHeader normMax
# TableCol 9 FieldName Intensity normMax
# TableCol 9 Units MR
# TableCol 10 ColHeader normRange
# TableCol 10 FieldName Intensity normRange
# TableCol 10 Units MR
# NRows 45
# NTableCols 10
# ColHeaders Index SegId NVoxels Volume_mm3 StructName normMean normStdDev normMin normMax normRange
1 4 7085 7398.2 Left-Lateral-Ventricle 18.8217 14.1739 0.0000 83.0000 83.0000
2 5 690 717.9 Left-Inf-Lat-Vent 33.4015 17.0720 1.0000 81.0000 80.0000
3 7 15413 15943.1 Left-Cerebellum-White-Matter 78.1067 8.3401 14.0000 100.0000 86.0000
4 8 51109 50979.1 Left-Cerebellum-Cortex 50.3740 10.4914 6.0000 89.0000 83.0000
...
40 85 136 148.6 Optic-Chiasm 80.4118 20.1430 35.0000 118.0000 83.0000
41 251 940 913.3 CC_Posterior 99.5362 10.8682 40.0000 110.0000 70.0000
42 252 526 468.3 CC_Mid_Posterior 85.4468 17.7985 12.0000 107.0000 95.0000
43 253 422 396.3 CC_Central 90.8957 16.3831 40.0000 108.0000 68.0000
44 254 406 385.3 CC_Mid_Anterior 97.4286 15.2189 40.0000 115.0000 75.0000
45 255 781 760.5 CC_Anterior 103.5237 14.4599 31.0000 121.0000 90.0000
これをもとに、関心領域の体積や信号値の平均などの情報を解析することができる。
.statsファイルを直接読み込んで、自身で解析してももちろんよいが、複数の被験者についてstatsファイルを集計して表を生成するpythonスクリプトasegstats2table、aparcstats2tableがFreeSurferに 標準で添付されている。asegstats2tableは体積画像に基づく情報(aseg.stats, wmparc.statsなど)、aparcstats2tableは皮質に関する情報の集計に用いる。
SUBJECTS_DIR="/path/to/your/dataset/derivatives/freesurfer" asegstats2table --subjects `cat subid_list` --tablefile aseg.stats.tsv
皮質のvertex-wise analysis
頂点ごとのデータ(thickness, sulc, area, curvなど)を一般線形モデル(general linear model; GLM)により解析することができる。
前処理
QCache を作成する。これにより、被験者の頂点に関連付けられた各測定値を、fsaverageの表面に対応付け、さらにGaussianカーネルによるsmoothingも行われる。
$ recon-all -s <subjid> -qcache
QDECによるGUI-baseでの解析
mri_glmfitによるCUI-baseでの解析
https://surfer.nmr.mgh.harvard.edu/fswiki/FsTutorial/GroupAnalysis