FSLを用いたシェルスクリプト演習：脳３次元T1強調MRIの水平断一覧画像の作成

あるプロジェクトにおいて、「MRI画像の一覧画像を作成してもらいたい」という依頼を受けました。

イメージとしては、下記のようなものです。

以前、FSLのシェルスクリプトチュートリアルを翻訳していた際に、演習に
Lightbox Viewerがあったのを思い出して、改めて、それを見ながら、さらに機能を追加したスクリプトを考えてみました。

作成しながら勉強になることがいくつもあったので、細かくわけてみたいと思います。

1. 剛体変換 (Rigid body transformation)

MRIは、人によって頭の位置が異なりますので、ある範囲内で画像を表示しようとすると、脳が変なところで切れてしまうことがあります。これを解決するために、個々人の脳を、標準脳に対して線形変換のうちの剛体変換 (rigid body transformation)を行おうと考えました。剛体変換は、脳を水平移動と回転だけするものです。つまり、脳の大まかな形は変わりません。水平移動と回転をx, y, z軸の３方向に対して行いますので、2(水平移動、回転)×3軸で、6つのパラメーターを使うことになります。

FSLでは、flirtというプログラムが線形変換を行うプログラムです。用法としては、以下になります。

flirt [オプション] -in 入力画像 -ref 参照画像 -out 出力画像

今、オプションに、剛体変換を指定したいと思います。flirtでは、-dof の後に数値を指定します。今は上記のように6パラメーターですので、-dof 6となります。

さらに、参照画像ですが、FSLの中に入っているMNI152_T1_1mm.nii.gzを使いたいと考えました。このファイルは、FSLの中のdata/standard/ディレクトリに入っています。FSLのパスは、${FSLDIR}で指定できます。これを使うことで、MacでもLinuxでも使えるスクリプトになります。

あとで、スクリプトを読みやすくするために、参照画像は変数 $ref で指定することにします。そして、出力画像は入力画像のファイル名の最初に”r”を付加することにします。

入力画像はforループを使って変数fにおさめられているとします。

すると、以下のようになります。

ref=${FSLDIR}/data/standard/MNI152_T1_1mm.nii.gz

for f in "$@"
do
  #Rigid body registration
  flirt -dof 6 -in $f -ref $ref -out r${f}
(以下続く)

これで、MNI152_T1_1mm.nii.gzと同じFOVの画像を得ることができます。

具体的に見てみましょう。

元々はこんな画像です。



これが、出力画像として、こんな画像になります。



首の部分がきれいに切れていることがわかります。また、AC-PCラインに沿って脳の位置がpitch方向に回転していることがわかります。

2. 水平断の抽出

次に線形変換後の画像から水平断の画像を抽出しようと考えました。線形変換後画像を見ると、z=26からz=145に脳がだいたい入っていることがわかりました。そして、5×6で表示したいと思っていましたので、30枚抽出できたらいいなと考えました。z=26からz=145だと120枚になりますので、4枚おきにすれば30枚選ぶことができます。

MRI画像から水平断を抽出し、PNG画像として保存するコマンドにFSLにはslicerというコマンドがあります。(3D-Slicerとは違うので注意が必要です)

水平断を抽出する用法としては、以下になります。

slicer 入力画像　-z -番号　出力画像

今、これを、z=26から4枚おきに画像を生成したいと思います。

4つごとに数字を出す方法に、UNIX系OSでは、seqというコマンドがあります。

seq 開始番号　間隔値　終了番号

なお、seqには-fオプションを使うことで、ゼロ埋めを行うことができます。つまり、26は026と表すことができます。これは、”%03g”となるそうです。

つまり、

seq -f "%03g" 26 4 145

となります。

今、出力画像は、ファイル名の後に連番をつけたいと思います。
ただ、ファイル名の拡張子はとりたいと思います。FSLには、remove_extというプログラムがあります。

これらをまとめて、以下のようにすることで、z=26から4枚ごとに30枚分の水平断画像を抽出することができます。バッククオートを使っていることに注意してください。

for i in `seq -f "%03g" 26 4 145`
  do
   file=`remove_ext $f`
   output=${file}_$i.png
   slicer r$f -z -$i $output
  done

これによって30枚のPNG画像を得ることができます。

3. PNG画像を配列し、1枚の画像とする

次に、この得られた30枚のPNG画像を5×6のタイル状に配置したいと思います。

FSLには、pngappendというプログラムがあります。

pngappendは次のような用法で使います。

pngappend 入力画像1 +/- 入力画像2 +/- 入力画像3 ... 出力画像

ここで、入力画像を”+”か”-“でつなぐのですが、”+”の時は、水平(横)方向に、”-“の時は垂直(縦)方向につなげることができます。

つまり、6つの画像(im1〜im6)を横3×縦2としたい場合、以下になります。

pngappend im1 + im2 + im3 - im4 + im5 + im6 出力画像

今、5×6の配置にしたいので、
im1 + im2 + im3 + im4 + im5 – im6 + im7 + …

としないといけません。
ただ、これを手で書くのは賢くないので、なんとか自動で書くように考えてみます。

いろいろ試行錯誤したあげく、今回は、以下の方法でいきたいと考えました。

まず、先ほど作成したPNG画像の一覧を ls で表示します。
その次に、5行目ごとに、”-” を行頭につけ、その他には、”+”を行頭につけます。
そして、その後、それらを改行を削除して1つに表示します。
これらを”inputarg”という変数に代入します。

ポイントは、「5行目ごとに」となります。

これにAWKを使うことにしました。
AWKでは、いくつか決められた変数があり、その中のひとつにNRというものがあります。これは、Number of Recordsの略で、行数を意味します。
今、6行目、11行目、16行目…にそれぞれ “-” をつけたいと思っています。
このルールがどこにあるかというと、「5で割ったときに余りが1」となります。

AWKにおいて、あまりは%で表示されます。今の場合だと、

NR%5==1

であらわすことができます。

また、AWKでは、1行1行処理していくわけですが、1行全部は”$0″で表示します。文字列を表示するには、print “文字列”となります。これを合体させると、もし、行数の前に”-“をつけたい場合には、

print " - " $0

で表示できます。

これらをまとめて、以下で表示できることがわかりました。

inputarg=`ls ${file}_???.png | \
awk '{ if (NR%5==1) {print " - " $0} else {print " + " $0} }'` 

こうすると無事に5行ごとに “-” が、それ以外は “+” がつきます。

その次に、改行をとりたいと思います。これは、UNIX系OSにある “tr” が便利です。trで改行を削除するには以下を使います。

tr -d '\n'

そして、最後に一番最初の “-” だけ取り除く必要があります。これは、sedを使うことで簡単にできます。

sed 's/-//'

これらをまとめると、以下のようになります。

inputarg=`ls ${file}_???.png | \
awk '{ if (NR%5==1) {print " - " $0} else {print " + " $0} }' | \
tr -d '\n' | sed 's/-//'`

出力画像は、ファイル名の前にLB (Lightboxの略とします)をつけることにしました。
従って、以下になります。

pngappend $inputarg LB_${file}.png

4. テンポラリファイルの削除

途中で作成した画像を削除します。

rm ${file}_???.png

以上で、完成です。
若干機能を加えてありますが、以下のようになりました。

#!/bin/bash
#Generate Lightbox View of 3D T1 image

FSLOUTPUTTYPE=NIFTI
ref=${FSLDIR}/data/standard/MNI152_T1_1mm.nii.gz

#Check if arguments are specified.
if [ $# -lt 1 ] ; then
 echo "Usage: $0 <filename> [options]" ;
 exit 1 ;
fi

for f in "$@"
do
  #Rigid body transformation
  if [ ! -e r${f} ]; then
    flirt -dof 6 -in $f -ref $ref -out r${f}
  fi

  #Obtain every 4 slices from z=26 to 145
  for i in `seq -f "%03g" 26 4 145`
  do
    file=`remove_ext $f`
    output=${file}_${i}.png
    slicer r$f -z -${i} $output
  done
  
  #Prepare argument for pngappend
  inputarg=`ls ${file}_???.png | \
    awk '{ if (NR%5==1) {print " - " $0} else {print " + " $0} }' | \
    tr -d '\n' | sed 's/-//'`

  #For debug
  #  echo $inputarg

  #Append png files 
  pngappend $inputarg LB_${file}.png
  
  #Delete temporary files
  rm ${file}_???.png
  #rm r$f
done

試してみたい方は、以下のスクリプトをダウンロードして試してみてください。

lightbox-3DT1.shをダウンロード（右クリックで名前をつけて保存で保存してください）