Scripts to convert hdr/img format into nii format

Nifti file havs two kinds of extensions, hdr/img or nii.
Sometimes hdr/file nifti files can be confusing because the extension is the same as ANALYZE format.
There are some ways to convert hdr/img into nifti.
One is to use fslchfiletype included in FSL. This is the easiest, but it requires FSL.
The other is to use SPM. However, some scripting is needed.
So I wrote tiny scripts for converting hdr/img format into nii format.

If you prefer SPM, use img2nii_spm.m. Save the script into a folder under MATLAB path and type “img2nii_spm” from MATLAB command window.

If you prefer FSL, use img2nii_fsl.sh. Save the script into a folder, add executable attribute using “chmod 755 img2nii_fsl.sh”.
Wild card can be used for selecting files.

Download img2nii_spm.m (Right-click and save).
Download img2nii_fsl.sh (Right-click and save).

SPM12におけるspm_defaults.mの設定 【重要】

SPM12の正式版が出て、気づいたのですが、spm_defaults.mの設定に関して、重要な変更が行われました。

今後は、spm_defaults.mは直接いじってはいけないことになりました。
これは、非常に理にかなっていることで、これまで、spm_defaults.mを修正しても、アップデートによってこれが上書きされてしまうことがありました。

これを回避するために、新たな方法が導入されました。
spm_defaults.mの最初の記載のスクリーンショットは以下のようになります。

spm_defaults_screenshot

シンプルにいうと、以下のようになります。

  • spm_defaults.mは編集しない
  • MATLABのパスが通っているところにspm_my_defaults.mを作成する

私がSPMをインストールして行う作業は、defaults.stats.maxmemを編集することです。これは、統計の際に、Matlabにメモリをどの程度まで使わせてあげるかという設定です。デフォルトは2^26バイトで、64MBになります。せっかくメモリを多く積んでいる場合、これはもったいないです。そのため、私は2^31(=2GB)としています。

Matlabから以下のようにタイプします。

edit spm_my_defaults.m

そこに以下のように記載します。

global defaults
defaults.stats.maxmem = 2^31;

スクリーンショットを以下に示します。

spm_my_defaults_screenshot

これを保存するのですが、MATLABのパスが通っているところならばいいのですが、
無難なのは、ドキュメントフォルダの下のMATLABフォルダもしくは、SPM12のフォルダ内に保存したらよいと思います。

このうえで、SPMを起動すると、きちんと認識されると、Matlabのコマンドウィンドウ内に以下の表示がなされるはずです。

spm_screenshot_modified_defaults

これで設定が終了です。これは重要な設定ですので、チェックするようにしてください。

包括脳MRIチュートリアルが12月13日に開催されます

毎年開催している包括脳MRIチュートリアルですが、
今年は12月13日(土)に開催されることになりました。

詳しい案内はこちらをご覧ください。

VBMチュートリアルに関しては、席はあと20席程度残っているとのことです。

今回は、

・画像のQC(岩手医大山下先生)
・多重比較補正についての理解(京大川口先生)
・SPMの結果を論文に載せるまでの具体的手順

といった内容を計画しています。

関心のある方は、こちらにある申し込み方法をご覧のうえ、ご応募ください。

Linuxでfind, sedを使いながら文字列を再利用してファイル名の途中に文字列を挿入する方法

最近、以下のようなことをやる必要がありました。

ファイル名が

201105_ID.jpg
20120608_ID.jpg
20130704_ID.jpg
201409_ID.jpg

となっているファイル(年月が必ずあり、日が時についている)に対して、
年の後にアンダースコアを入れたい、つまり、

2011_05_ID.jpg
2012_0608_ID.jpg
2013_0704_ID.jpg
2014_09_ID.jpg

としたいわけです。

もし、最初の4文字が全く同一の文字列ならばただ置換すればよいので、全く問題ありません。
たとえば、すべて最初が2011からはじまるのであれば、
renameを使って、

find . -name '2011*' | rename 's/^2011/2011_/'

とするだけです。

しかし、今回は、4文字ということでは規則性がありますが、それぞれの中はバラバラの値です。

いろいろ調べる中で、いくつかポイントが見つかりました。

 

  • findのオプション、-execを組み合わせる
  • findで合致したファイル名を再利用するには、{}を使う
  • sedにおいて検索で使った文字列を再利用するには、&を使う
  • bashで文字列を引数にとるためには、-cオプションを使う
  • find -execオプションは必ず;で終わらないといけないが、bashではエスケープシーケンスが必要になるために、必ず\;となる。

最初に、今回うまくいったワンライナーをご紹介します。

find . -name '201*' -exec bash -c "mv {} \`echo {} | sed -e 's/201./&_/'\`" \;

これを(あとで自分が忘れないために)もう少し解説します。

find . -name '201*'

この部分は、findの基本的な使い方です。
find <検索したい場所> -name <ファイル名> です。
ワイルドカードを使いたい場合は、シングルクオートかダブルクオートでくくらなければいけません。
*の場合は、どちらでも大丈夫です。

-exec bash -c

findのオプションで-execを使うと、コマンドを実行することができます。

今は、複雑なことをやりたいので、シェルを呼び出します。そのために、bashとしています。
man bashとすると以下のように書かれています。

オプション
bash は以下のオプションを起動時に解釈します (組み込みコマンド set の説
明で述べられている 1 文字のシェルオプションも使えます):

-c string -c オプションが指定されると、コマンドが string から読み込まれ
ます。 string の後に引き数があれば、これらは 位置パラメータ
(positional parameter: $0 からはじまるパラメータ) に代入され

今は複数のコマンドを読み込みたいので、-cオプションを指定します。

"mv {} \`echo {} | sed -e 's/201./&_/'\`" \;

これが今回のキモです。

まず、bash -cの後は文字列ですので、それを明記するために全体をダブルクオーテーションでくくっています。
そして、find -execの最後は;で終わるとなっているのですが、エスケープシーケンスが必要のために\;とします。

残るは mv {} \`echo {} | sed -e ‘s/201./&_/’\` です。

mv の後の{} はfindで見つかったファイル名を意味します。ここで、findの結果を利用できるわけです。

そして、\`echo {} | sed -e ‘s/201./&_/’\`

がmvの後半部分になります。これは全体がバッククオート`でくくられています。エスケープさせる必要があるために、\`となっています。

バッククオートの中を見ると、

echo {} | sed -e ‘s/201./&_/’

となっていて、先ほどと同様に{}はfindの結果を利用するわけで、echo {}で文字列を表示したうえで、それをsedを使って文字列置換をするわけです。

sedの文字列置換は正規表現を用いることができますので、検索文字列は2011,2012,2013,2014のいずれかですから、201.としました。
そして、その後が大事なのですが、sedの検索文字列をそのまま使える方法が、&となります。
今の場合、見つけた文字列の後にアンダースコアをつけたいので、&_としました。

これで無事にファイル名を変更することができました。数千のファイルをたった1行のコマンドで一気に変更できましたので、かなりすっきりしました。といってもここにたどり着くまでに数週間かかりましたが…。

bashでできる技なので、Linuxに限らず、Macでも同様にできるかと思われます。
ただ、これは我流でたどり着いた方法なので、もっとエレガントな方法がありましたら、ぜひ教えていただきたいです。

SPM8とSPM12を同一のMatlabのライセンスで使用する方法

SPMのメーリングリストでSPMの中の人が明らかにしていますが、SPM12の正式版がそろそろ公開されるようです。SPM12が正式に公開されました。(2014年10月1日)

SPM12はいろいろよい点も加わっていますが、まだ、プラグインで対応していないものも多々あります。
このため、しばらくは、SPM8とSPM12が同時に動かせる環境にあった方がよいですよね。

いい方法がありますので、ご紹介します。

以下の2つのファイル(spm8.mとspm12.m)をMatlabのパスが通っているフォルダ(ディレクトリ)に保存します。
Windowsの場合、「ドキュメント」の下のMATLABフォルダがよいでしょう。
Macの場合、「書類」の下のMATLABフォルダがよいでしょう。
これらはパスが通っているフォルダです。

そして、最初の4行と最後の方のaddpathからはじまる行を自分の環境にあわせて調整します。

このスクリプトがしているのは3つです。

  1. SPMのパスとツールボックスのパスを変数として指定します。
  2. SPMのパスが既にある場合、spm_rmpathを用いてそのパスを取り除きます。
  3. 改めてSPM8もしくはSPM12のパスを設定します。

以下にspm8.mのスクリプトを例に説明していきます。

%%spm8.m
%%specify spm8 and other toolbox path
spm8_path = 'c:\spm\spm8';
pickatlas_path = 'c:\spm\spm8\toolbox\wfu_pickatlas';
%spm-related_path = 'path_to_the_application';

%からはじまる行はすべてコメントなので気にしなくていいです。
spm8_path = のところで、SPM8がインストールされているフォルダを設定します。
pickatlas_path = のところで、(例として)PickAtlasのフォルダを設定します。
同様に、他にも設定しなくてはいけないプラグインがあったら、ここに設定します。
(これまで、設定する必要があったのは、PickAtlasとAnatomyツールボックスぐらいです。
その他はあえて設定する必要はありません)
これを後に使います。

% remove spm path
while true
  try
    spm_rmpath;
  catch
    break
  end
end

ここはいじる必要がありません。
まず、spm_rmpathというコマンドをトライします。
もし、SPMにパスが通っていれば、このコマンドが実行され、SPMのパスが取り除かれます。
もし、パスが通っていなければ、実行されずにそのまま次に進みます。

% add spm8 and related path
addpath(spm8_path,pickatlas_path);

次に、addpathを使ってパスを設定します。
addpathの中に先程設定したパスを書きこんでいきます。順序が大事で、先にspm8_pathを書き、その次にツールボックスのパスを書きます。

% run spm8
spm;

最後にSPMを起動します。

関心のある方は、下記の2つをダウンロードして、ご自身の設定にあわせて修正してご利用ください。

spm8.mをダウンロード(右クリックで「名前をつけ保存」から保存してください)
spm12.mをダウンロード

心理のためのMatlabチュートリアル:解答例と解説

何人かの方々から、心理のためのMatlabチュートリアルの練習の解答がないか問い合わせを受けていました。
これまで、ちまちま作っていたのですが、練習Fまではできていますので、公開させていただきます。
練習GとHに関してはまだできていないので、できたら、追加させていただきます。

心理のためのMatlabチュートリアルの解答例をダウンロード

Update of Lin4Neuro (10 Aug 2014)

Lin4Neuro was updated. This time the change was limited to kernel update and security updates.

  • Kernel was updated to 3.13.0-32-generic.
  • All security updates were applied.

I’m working on Xubuntu 14.04, but it seems to have several problems with neuroimaging software packages. So for now I keep updating Lin4Neuro based on Xubuntu 12.04.

You can download the latest Lin4Neuro from here.

RでのROC解析:ROCRパッケージを使ったROC曲線とAUCの求め方

研究でROC解析を行う必要があり、Rでどうやったらできるのか調べてみました。
そうしたところ、ROCRというパッケージが公開されており、比較的簡単にROC解析を行い、グラフを作成できることがわかりました。

  • ROCRパッケージのインストール(Ubuntu)
  • 既にRはインストールされているとします。Ubuntuの場合、ROCRパッケージはapt経由で簡単に入手できます。

    $ sudo apt-get install r-cran-rocr

    これでインストール完了です。

  • ROCRを使うための準備
  • ROC解析に必要なものは、何らかの指標と、それが属するグループの一覧です。具体例を挙げると、以下のようになります。
    第1列に指標、第2列に属するグループ(0か1)が記載されています。
    これをroc_data.txtとという名前で保存することとします。保存したディレクトリをRのワーキングディレクトリとします。

    0.9706 1
    0.9572 1
    0.4854 1
    0.8003 1
    0.1419 1
    0.4218 1
    0.9157 1
    0.7922 1
    0.9595 1
    0.6557 1
    0.0357 1
    0.8491 1
    0.934 1
    0.6787 1
    0.7577 1
    0.7431 1
    0.3922 1
    0.6555 1
    0.1712 1
    0.706 1
    0.4797 0
    0.4551 0
    0.0374 0
    0.081 0
    0.2984 0
    0.7597 0
    0.1404 0
    0.3853 0
    0.0238 0
    0.5513 0
    0.0551 0
    0.306 0
    0.4991 0
    0.6909 0
    0.7593 0
    0.3472 0
    0.0614 0
    0.0507 0
    0.0575 0
    0.6407 0

  • ROCRの起動
  • ROCRはRを立ち上げた後に、library(ROCR)で起動できます。

    $ R
    > library(ROCR)

  • データの読み込み
  • 先ほどのroc_data.txtをrocdataという変数に読み込みます。変数名は何でもいいのですが、ここではそうします。
    read.tableという関数で表を読み込めるので、それを使います。

    > rocdata <- read.table("roc_data.txt")

  • ROCR ステップ1: prediction
  • ROCRのステップ1はpredictionで、値と属するグループを指定します。
    今、rocdataは20行2列の行列になっています。1列目はrocdata[,1]で、2列目はrocdata[,2]であらわすことができますので、以下のように記載します。

    pred <- prediction(rocdata[,1], rocdata[,2])

    ここで、変数predはpredictionの頭文字4文字です。もちろん、別の名前でもかまいません。
    私はよくカンマを忘れるので、カンマも忘れないようにしましょう。

  • ROCR ステップ2: performance
  • ROCRのステップ2はperformanceでここでプロットするための項目を決定します。今は、基本的なROC曲線、すなわち横軸がfalse positive rate, 縦軸がtrue positive rateで構成される曲線を描こうと考えます。以下のようにタイプします。

    perf <- performance(pred, "tpr", "fpr")

    ここでtprはtrue positive rateを、fprはfalse positive rateを意味します。

  • ROCR ステップ3: グラフの描画
  • 最後にROC曲線を描きます。非常に簡単です。

    plot(perf)

    そうすると、次のようなグラフが現れるはずです。

    roc-curve

  • グラフの保存
  • グラフをPNG形式で保存するには、次のように行うことで、roc-curve.pngという名前でワーキングディレクトリに保存されます。

    > png("roc-curve.png")
    > plot(perf)
    > dev.off()

  • AUCの算出
  • 先ほどのperformanceの際に"auc"と指定するとAUCも計算されます。ただ、1クッション入れる必要があります。具体的な方法はこちらのサイトに記載されていましたが、それを転載します。

    > auc.tmp <- performance(pred,"auc")
    > auc <- as.numeric(auc.tmp@y.values)

    performance(pred,"auc")の結果をauc.tmpという変数に代入し、
    auc.tmpの中からy.valuesの値を取り出して、その値を変数aucに代入します。

    最後に変数aucを表示させてみます。

    > auc
    [1] 0.8

    これでAUCが0.8だということがわかります。

  • おまけ:感度特異度曲線
  • 感度と特異度の曲線も簡単に書けます。predictionまで行った後に、次のようにします。

    > perf <- performance(pred, "sens", "spec")
    > png("sens-spec-curve.png")
    > plot(perf)
    > dev.off()

    これで、下図のような感度、特異度の曲線がsens-spec-curve.pngという名前で保存されます。

    sens-spec-curve

比較的簡単に求められるので便利です。

Workaround for Remastersys on Xubuntu 14.04

Remastersys is a very useful remastering tool though developing has ceased now.
Many reported that they could make Remastersys work on Ubuntu 14.04, but I encountered an error with Xubuntu 14.04.

I found workaround for this error.

Short answer: Comment out the line 245 beginning with [ “`grep lightdm /etc/X11/default-display-manager`” != “” ]

Below is the long answer.

  • Installation
  • Installation is simple.

    1. Import GPG key
    2. $ wget -O- http://www.remastersys.com/ubuntu/remastersys.gpg.key | sudo apt-key add -

    3. Add sources
    4. Add the line to /etc/apt/sources.list

      deb http://www.remastersys.com/ubuntu precise main

      Please note that you don’t have to replace “precise” with “trusty.”

    5. Installation
    6. $sudo apt-get update; sudo apt-get install remastersys

  • Error with “remastersys dist”
  • “sudo remastersys dist” terminated with error. The below is error message;

    $ sudo remastersys dist
    Distribution Mode Selected
    Enabling remastersys-firstboot
    Adding system startup for /etc/init.d/remastersys-firstboot …
    /etc/rc0.d/K20remastersys-firstboot -> ../init.d/remastersys-firstboot
    ……
    /etc/rc5.d/S20remastersys-firstboot -> ../init.d/remastersys-firstboot
    Checking filesystem type of the Working Folder
    /home/remastersys/remastersys is on a ext4 filesystem
    Making sure popularity contest is not installed
    Installing the Ubiquity GTK frontend
    Lightdm not setup properly. You must set your default desktop with lightdm prior to remastering

    I examined the script /usr/bin/remastersys with the keyword “Lightdm.” Then I found the description below in the line 245

    [ “`grep lightdm /etc/X11/default-display-manager`” != “” ] && [ ! -f /etc/lightdm/lightdm.conf ] && [ ! -f /usr/share/xsessions/ubuntu.desktop ] && echo “Lightdm not setup properly. You must set your default desktop with lightdm prior to remastering” && echo “Lightdm not setup properly. You must set your default desktop with lightdm prior to remastering” >> $WORKDIR/remastersys.log && exit 1

    Though this line checks if the file “/etc/lightdm/lightdm.conf” exists, The latest Xubuntu 14.04 doesn’t have that lightdm.conf file in /etc/lightdm.

    So I commented out this line, just adding “#” to the beginning of the line.

    Then I tried sudo remastersys dist again, with success!

    If you want to use Remastersys in Xubuntu 14.04, it’s worth trying it.

Update of Lin4Neuro (12 Jun 2014)

Lin4Neuro was updated. The changes of this version were as follows;

  • Kernel was updated to 3.13.0-29-generic.
  • ITK-SNAP was updated to 3.0.0.
  • MITK was updated to 2014.03.00.
  • MITK Diffusion was updated to 2014.03.
  • All security updates were applied.

Lin4Neuro is still based on Xubuntu 12.04. Now I’m working on Xubuntu 14.04, so hopefully the next version will be based on Xubuntu 14.04.

You can download the latest Lin4Neuro from here.

VBMミニ勉強会@秋葉原 6月28日15:00-18:00

5人の方の申し込みをいただいたので、開催決定です!
岩手医大の山下典生先生も助っ人で参加してくださることになりました!

※場所が確定しました。秋葉原駅徒歩5分のAP秋葉原です。
事前登録制です。参加希望の方はコメント欄からご連絡をお願いします。

すぐできるVBMが無事に出版に至りましたが、「本は買ったけれどもやる時間がない」という方もいらっしゃるのではないかと思います。また、「ここはどうなの?」と質問をしたいという方もいらっしゃると思います。

そこで、VBMのミニ勉強会を開催させていただきます。

日時:2014年6月28日(土)15:00-18:00
場所:AP秋葉原
参加費:3000円程度(会議室の実費代のみ。)
持参していただくもの:SPMがセットアップされている自分のパソコン
人数:最小5名〜最大20名
テキスト:「すぐできるVBM」(前もってご準備ください)
内容:実習形式
 ・VBMを実際に動かしてみる
 ・自分の解析でわからないことの質疑応答

参加希望の方は、下記コメントにお名前(ニックネーム可)とメールアドレスを書きこんでいただけますでしょうか。
VBMに関心のある方はご検討ください。

第1回 IMU脳画像解析セミナー「VBM/ROI解析を極める」@岩手

友人の山下典生先生が、2014年8月31日に『第1回 IMU脳画像解析セミナー「VBM/ROI解析を極める」』を開催します。

SPMでのVBMだけでなく、3D-SlicerやITK-SNAPも利用するとのこと。とても楽しみです。
早速申し込みしました。

関心のある方はこちらをごらんください。

VBMでの全脳容積(Total Brain Volume; TBV)の求め方

「すぐできるVBM」に、全脳容積の求め方を記載していませんでした…。

以下の方法で、求めることができます。

必要なものは以下のとおりです。

  • get_totals.m
  • 容積を簡便に求めることができるスクリプトは、Ged Ridgway氏によるget_totals.mです。右クリックで「名前をつけて保存」とし、SPMのディレクトリに保存してください。

  • 灰白質画像c1*.niiと白質画像c2*.nii
  • Segmentで一番最初に作られる画像であるc1*.niiとc2*.niiを準備します。ここで*はワイルドカードを意味し、「c1からはじまり、そのあとは何でもOKで、最後は.niiで終わる」という意味です。

それでは実際にやっていきましょう。

まずは、灰白質の容積を求めます。Matlabのコマンドウィンドウから以下のようにタイプします。

>> gmv = get_totals

(gmvは変数名なので、gray_matterでも、gmでもかまいません)
そうすると、ファイルを選択するダイアログがあらわれます。ここで、容積を求めたいc1画像を指定してください。ここでは、ためしに5つのファイルを選びました。そうすると、以下のような表示になると思います。(もちろん、数字は選択した画像で変わります)

gmv =
640.9410
664.4302
617.7647
874.1366
670.1864

次に、白質の容積を求めます。Matlabのコマンドウィンドウから以下のようにタイプします。

>> wmv = get_totals

同様に、ファイルを選択するダイアログがあらわれますので、容積を求めたいc2画像を指定します。重要なのは、c1画像と同じ順番で選択してください。以下のような表示になると思います。

wmv =
431.9188
463.1192
418.5571
572.9889
488.4050

そうしたら、2つを足して全脳容積を求めましょう。

>> tbv = gmv + wmv

tbv =
1.0e+03 *
1.0729
1.1275
1.0363
1.4471
1.1586

これは、1.0×10^3(すなわち1000)にそれ以下の数字をかけたものということです。つまり、

tbv =
1072.9
1127.5
1036.3
1447.1
1158.6

ということになります。

これを選択してコピーしても使えますが、ファイルに書き出しましょう。
コンマ区切りファイルというものに出力することができます。
出力ファイル名をtotal_brain_volume.csvとします。

>> csvwrite('total_brain_volume.csv',tbv)

csvwriteの書式は csvwrite(出力ファイル名, 変数名) となります。
変数tbvの内容をtotal_brain_volume.csvに出力したいので、上記のようになります。
そうすると、ワーキングディレクトリにtotal_barin_volume.csvが生成され、Excelのような表計算ソフトでも、テキストエディタでも開くことができます。

最後に少し応用しましょう。灰白質容積(GMV), 白質容積(WMV), 全脳容積(TBV)をすべてCSVファイルに出力したいとしましょう。

その場合は、まず、すべての情報をひとつの変数volumesにまとめます。
Matlabのコマンドウィンドウから次のようにタイプしてください。

>> volumes = [gmv wmv tbv]

そうすると下記のようになります。

volumes =
1.0e+03 *
0.6409 0.4319 1.0729
0.6644 0.4631 1.1275
0.6178 0.4186 1.0363
0.8741 0.5730 1.4471
0.6702 0.4884 1.1586

左からgmv, wmv, tbvという順番になっています。

これを、先程と同じようにしてCSVファイルに保存します。ファイル名をvolumes.csvとします。

>> csvwrite('volumes.csv',volumes)

このようにして、脳容積の情報を保存できます。

SPMでの「1対多」の解析

臨床では、1人の症例で健常者と比較して容積低下部位を求めたいことがあります。
SPMでどのようにするかというと、シンプルです。

  • DesignはTwo sample t-testを用います
  • Group 1に1例のみ選びます
  • Group 2に健常者の複数の画像を指定します
  • IndependenceはYesを選びます
  • (ここがキモ)VarianceはEqualを選びます。もし、Unequalにするとエラーとなります。1例では分散を計算することができないからです。

そうすると、Design matrixを作るための設定画面は下図のようになるはずです。

SingleCase

これを普通にEstimateして、Contrast managerで症例の容積が低下している領域を求めたかったら、[-1 1]とすれば大丈夫です。

Twitterで疑問を投げかけてくださった@silverjet_jpさん、ありがとうございました。

MRIConvertやdcm2niiを用いたDICOM→NIFTI変換の方法

画像解析のためには、DICOM画像をNIFTI画像に変換する必要があります。
主に使われるソフトウェアとして、MRI Convertとdcm2niiがあります。
開発者からスクリーンショットを使ってよいという許可をいただきましたので、
MRI Convertおよびdcm2niiの使用方法をPDFで公開します。

MRIConvertの使用方法(PDF)
dcm2niiの使用方法(PDF)