【FSL】FLICAを用いたマルチモーダル独立成分分析


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. 目的

  • FLICAを用いた画像のマルチモーダル独立成分分析(ICA)

2. リポジトリ

FLICAに必要なスクリプトを次のコマンドで取得。今のところプライベートリポジトリにしています。ご興味のある方はコメントしてください。

git clone https://github.com/Hexans/FLICA_corrected.git

3. 必要なデータ

必要なデータは、モダリティごとに全ての被験者の画像がまとめられた4D画像である。例えば、次のようなファイルを用いる。

  • VBM: GM_mod_merg_s?.nii.gz
  • TBSS: all_??_skeletonised.nii.gz
  • FreeSurfer: ?h.thick.fsaverage.mgh and ?h.pial.area.fsaverage.mgh
  • etc.

さらに、共変量(年齢・疾患グループ・性別)との相関もみることができる。共変量ごとに1列に数値をまとめテキストファイル(.txt)として保存する。

  • 年齢: age.txt
  • 疾患グループ: group.txt
  • 性別: sex.txt

4. 注意点

NANボクセルがある場合、FLICA処理でエラーが発生する。

例えば、次のようなファイルがあったとする。

.
├── all_median_FA_skeletonised.nii.gz
├── all_median_ICVF_skeletonised.nii.gz
├── all_median_ISO_skeletonised.nii.gz
├── all_median_MD_skeletonised.nii.gz
├── all_median_MVF_skeletonised.nii.gz
└── all_median_OD_skeletonised.nii.gz

NANボクセルがあるかは、FSLのfslstatsで検出できる。-Mオプションは、non zeroボクセルにおける平均値を算出するオプションである。

Method:

fslstats <INPUT> -M

次のコマンドを実行することで、すべての画像に対して、NANボクセルがありなしをチェックする。

In:

for image in *nii.gz; do
    echo "${image} : $(fslstats $image -M)"
done

以上のコマンドを実行した例が、以下である。この場合、「all_median_ICVF_skeletonised.nii.gz」と「all_median_MVF_skeletonised.nii.gz」にNANボクセルが含まれていることが分かる。

Out:

all_median_FA_skeletonised.nii.gz : 0.157853 
all_median_ICVF_skeletonised.nii.gz : -nan 
all_median_ISO_skeletonised.nii.gz : 0.246872 
all_median_MD_skeletonised.nii.gz : 1.141417 
all_median_MVF_skeletonised.nii.gz : -nan 
all_median_OD_skeletonised.nii.gz : 0.516476

NANボクセルをゼロに置換するには、FSLのfslmathsコマンドを用いる。

Method:

fslmaths <INPUT> -nan <OUTPUT>

すべての画像に対して、NANボクセルがある場合、0に置換する処理をするには、以下のコマンドを実行する。以下の処理が完了すると、「corrected_nan」フォルダに、NANボクセルが含まれない画像が保存される。

mkdir corrected_nan
for image in *nii.gz; do
    fslmaths $image -nan corrected_nan/${image}
done

NANボクセル補正後のデータに、NANが含まれていないかチェックするには、以下のコマンドを実行する。

In:

for image in $(ls corrected_nan); do
    echo "${image} : $(fslstats corrected_nan/$image -M)"
done

「all_median_ICVF_skeletonised.nii.gz」と「all_median_MVF_skeletonised.nii.gz」にあったNANボクセルが、解消されていることが分かる。

Out:

all_median_FA_skeletonised.nii.gz : 0.157853 
all_median_ICVF_skeletonised.nii.gz : 0.422611 
all_median_ISO_skeletonised.nii.gz : 0.246872 
all_median_MD_skeletonised.nii.gz : 1.141417 
all_median_MVF_skeletonised.nii.gz : 0.086699 
all_median_OD_skeletonised.nii.gz : 0.516476

5. 実行

FLICAを実行するには、MATLABから次のコマンドを実行。

ICAのComponent数を被験者数の4分の1以下(Components < subjects/4)に設定し、Overfittingを避ける。また、被験者ごとのノイズ推定を省く(opts.lambda_dims = ”;)ことでもOverfittingを避けることができる。

5.1. run_FLICA.m

%% Set up:
addpath([getenv('FSLDIR') '/etc/matlab/'])
addpath('~/Documents/software/flica')

%% Load data
Yfiles = {'all_FA_skeletonised.nii.gz',
    'all_MD_skeletonised.nii.gz',
    'all_AD_skeletonised.nii.gz',
    'all_RD_skeletonised.nii.gz'};
% NOTE that these should be downsampled to around 20k voxels in the mask,
% per modality... hoping to increase this memory/cpu-related limitation.
[Y, fileinfo] = flica_load(Yfiles);
fileinfo.shortNames = {'FA', 'MD', 'AD', 'RD'};

%% Run FLICA
clear opts
opts.num_components = 5;
opts.maxits = 100;
% opts.lambda_dims = '';  # Reduce subject-noise estimation to avoid overfitting

Morig = flica(Y, opts);
[M, weights] = flica_reorder(Morig);

%% Save results
outdir = [pwd '/result/'];
mkdir(outdir)
flica_save_everything(outdir, M, fileinfo);

%% Produce correlation plots
clear des
des.Group = load('group.txt');
des.Age = load('age.txt');
des.Gender = load('gender.txt');
flica_posthoc_correlations(outdir, des);

実行が成功すると次のように画像を読み込み、計算が始まる。

Loading "all_FA_skeletonised.nii.gz"... 
Converting back to single... Generating mask... Flattening data matrix... Freeing memory... 
Loading "all_MD_skeletonised.nii.gz"... 
Converting back to single... Generating mask... Flattening data matrix... Freeing memory... 
Loading "all_AD_skeletonised.nii.gz"... 
Converting back to single... Generating mask... Flattening data matrix... Freeing memory... 
Loading "all_RD_skeletonised.nii.gz"... 
Converting back to single... Generating mask... Flattening data matrix... Freeing memory... 
Setup... doing PCA... more setup... done (0.72 sec)

Iteration 1... eta H W X.... lambda.... done updates (0.2 sec) F (0.061 sec)No previous lastF
Iteration 2... eta H W X.... lambda.... done updates (0.17 sec) F (0.029 sec), dF = 6.38536e+07
Iteration 3... eta H W X.... lambda.... done updates (0.16 sec) F (0.027 sec), dF = 21559.9
Iteration 4... eta H W X.... lambda.... done updates (0.13 sec) F (0.019 sec), dF = 253.743
...
Iteration 100... eta H W X.... lambda.... done updates (0.13 sec) F (0.019 sec), dF = 253.743

計算が終わったら、bashでƒƒflica_report.shを実行することで、HTML形式で結果をまとめることができる。まとめられた結果(HTML)を見るには「index.html」を開く。

flica_report.sh result  # 引数としてresultフォルダを指定

6. 結果

6.1. index.html

index.htmlの開くと次のようにFLICAの、まとめの結果を見ることができる。

6.2. コンポーネントごとの結果(All_in_one_page)

ICAの各Componentごとの結果は、index.htmlの「All_in_one_page」から確認することができる。

6.3. モダリティごとの結果(FA, MD, AD, RD, etc.)

モダリティごとの結果は、index.htmlの「モダリティ名」をクリックすることで確認できる。以下の図は、FAの結果である。

7. GBSSへの適応

GBSSの結果に対してもFLICAが使えるか、不明であったためできるか試してみた。

特に問題なく、使えそうな感じ。

7.1. index.html

7.2. コンポーネントごとの結果(All_in_one_page)

7.3. モダリティごとの結果(FA, MD, AD, RD, etc.)

【QSM】QSM解析 ~MEDIの使い方~


1. 目的
2. 準備
2.1. MEDI toolboxにパスを通す
2.2. QSM解析実行のコード準備
2.3. Write_DICOM.mを修正
2.4. QSM_3Dを集める
3. 実行
3.1. 作業フォルダ(analyzeフォルダ)へ移動
3.2. runme.mを実行
3.3. Outputファイルの確認
4. 複数の被験者データをまとめて処理したい場合
4.1. run_sequential.m


1. 目的

Cornell大学が開発したQSM software (MEDI)を使って解析

(Cornell MRI Research Lab. QSMより、URL: http://pre.weill.cornell.edu/mri/pages/qsm.html)

2. 準備

####MEDI toolboxのダウンロード
こちらからMEDI_toolboxをダウンロード
MEDI_toolboxの中身はこのような感じ。

2.1. MEDI toolboxにパスを通す

MATLABを起動し、パスの設定からサブフォルダも追加をクリックしfunctionsフォルダを選択した後、保存。

2.2. QSM解析実行のコード準備

QSM解析を実行するためのコードは以下。
以下をrunme.m(※任意の名前でOK)で保存し、作業フォルダとなるanalyzeフォルダ(※任意の名前でOK)に保存する。

clear all;clc;close all;

% Set path
% MEDI_set_path

% STEP 1: Import data
[iField,voxel_size,matrix_size,CF,delta_TE,TE,B0_dir,files]=Read_DICOM('AXL_QSM');

%%%% Generate the Magnitude image %%%%
iMag = sqrt(sum(abs(iField).^2,4));

%%%%% provide a Mask here if possible %%%%%%
if (~exist('Mask','var'))                     
%     Mask = genMask(iField, voxel_size);
    Mask = BET(iMag,matrix_size,voxel_size);
end

%%%%% provide a noise_level here if possible %%%%%%
if (~exist('noise_level','var'))
    noise_level = calfieldnoise(iField, Mask);
end

%%%%% normalize signal intensity by noise to get SNR %%%
iField = iField/noise_level;

% STEP 2a: Field Map Estimation
%%%%%Estimate the frequency offset in each of the voxel using a 
%%%%%complex fitting %%%%
[iFreq_raw N_std] = Fit_ppm_complex(iField);

% STEP 2b: Spatial phase unwrapping %%%%
iFreq = unwrapPhase(iMag, iFreq_raw, matrix_size);

% STEP 2c: Background Field Removal
%%%% Background field removal %%%%
[RDF shim] = PDF(iFreq, N_std, Mask,matrix_size,voxel_size, B0_dir);

% CSF Mask for zero referencing
R2s=arlo(TE,abs(iField));
Mask_CSF = extract_CSF(R2s, Mask, voxel_size);

% STEP 3: Dipole Inversion
save RDF.mat RDF iFreq iFreq_raw iMag N_std Mask matrix_size...
     voxel_size delta_TE CF B0_dir Mask_CSF;

%%%% run MEDI %%%%%
QSM = MEDI_L1('lambda',1000, 'smv', 5, 'merit');

%%% save to DICOM
% ignore warnings...
Write_DICOM(QSM,files,'QSM')

2.3. Write_DICOM.mを修正

解析後のQSMの保存には、Write_DICOM.mを使用する。(runme.mを実行すれば自動で実行される。)

修正前のファイルでは、dicomwriteでエラーが発生する。
そのため、MEDI_toolbox/functions/Write_DICOM.mを以下のコードにそっくりそのまま書き換える。
やっていることは、dicomwriteのCreateModeをCopyからCreateへチェンジ。

function Write_DICOM(M,files,outdir,opts)
%WRITE_DICOM Summary of this function goes here
%   Detailed explanation goes here


defopts.SeriesDescription = 'QSM';
defopts.SeriesInstanceUID = [];
defopts.SeriesNumber = [];
defopts.Convert = @(x) convert(x);
defopts.Window = 0.500;
defopts.Level = 0;
defopts.flag3D = 0;
% defopts.EchoNumber = [];
% defopts.EchoTime = 0.0;
        
if nargin<4; opts=struct; end
deffields=fieldnames(defopts);
for i=1:length(deffields)
    if ~isfield(opts, deffields{i})
        opts.(deffields{i})=defopts.(deffields{i});
    end
end

if isfield(files,'M')
    filenames=files.M;
elseif isfield(files,'R')
    filenames=files.R;
elseif isfield(files, 'info3D')
    opts.flag3D=1;
else
    error('No filenames (M nor R) found.');
end
    
flag_signed=min(M(:))<0;

if ~opts.flag3D
    if size(M,3) ~= size(filenames,1)
        error([num2str(size(filenames,1)) ' filenames given for ' num2str(size(M,3)) ' slices.']);
    end
end

if isempty(opts.SeriesInstanceUID)
   opts.SeriesInstanceUID=dicomuid;
end
progress='';

mkdir(outdir);

warning('off','images:dicom_add_attr:invalidAttribChar');
load_flag=1;
insert_flag=~opts.flag3D;
for slice=1:size(M,3)
    if load_flag
        if opts.flag3D
            filename=files.info3D.filename;
        else
            filename=filenames{slice,end};
        end
        info = dicominfo(filename);
        info.SeriesDescription = opts.SeriesDescription;
        info.SeriesInstanceUID = opts.SeriesInstanceUID;
        if isempty(opts.SeriesNumber)
            opts.SeriesNumber=info.SeriesNumber*100;
        end
        info.SeriesNumber = opts.SeriesNumber;
        info.SOPInstanceUID = dicomuid;
        info.InstanceNumber = slice;
        if opts.flag3D
            load_flag=0;
        end
    end
    if opts.flag3D
        item=files.info3D.slice2index{slice};
%         info.PerFrameFunctionalGroupsSequence.(item).PlanePositionSequence.Item_1.ImagePositionPatient;
        info.PerFrameFunctionalGroupsSequence.(item).Private_2005_140f.Item_1.SOPInstanceUID = dicomuid;
    end
    im = M(:,:,slice);
    if (isfield(info, 'SmallestImagePixelValue'))
        info.SmallestImagePixelValue=opts.Convert(min(im(:)));
    end
    if (isfield(info, 'LargestImagePixelValue'))
        info.LargestImagePixelValue=opts.Convert(max(im(:)));
    end
    if (isfield(info, 'RescaleIntercept'))
        info.RescaleIntercept=0;
    end
    if (isfield(info, 'RescaleSlope'))
        info.RescaleSlope=1;
    end
    info.WindowCenter=opts.Convert(opts.Level);
    info.WindowWidth=opts.Convert(opts.Window);
%     if opts.flag3D
%         info.PerFrameFunctionalGroupsSequence.Item_1.PixelValueTransformationSequence.Item_1.RescaleIntercept=0;
%         info.PerFrameFunctionalGroupsSequence.Item_1.Private_2005_140f.Item_1.RescaleIntercept=0;
%         info.PerFrameFunctionalGroupsSequence.Item_1.PixelValueTransformationSequence.Item_1.RescaleSlope=1;
%         info.PerFrameFunctionalGroupsSequence.Item_1.Private_2005_140f.Item_1.RescaleSlope=1;
%         info.PerFrameFunctionalGroupsSequence.Item_1.FrameVOILUTSequence.Item_1.WindowCenter=opts.Convert(opts.Level);
%         info.PerFrameFunctionalGroupsSequence.Item_1.FrameVOILUTSequence.Item_1.WindowWidth=opts.Convert(opts.Window);
%     end
	info.SamplesPerPixel=1;
    info.BitsAllocated=16;
    info.BitsStored=16;
    info.HighBit=15;
    info.PixelRepresentation=flag_signed;
    if size(M,3)==slice
        insert_flag=1;
    end
    if insert_flag
        outfile=fullfile(outdir,sprintf('IM%05d.dcm', slice));
        print_progress(outfile);
        if opts.flag3D
           f=fieldnames(info.PerFrameFunctionalGroupsSequence);
           f=setdiff(f,files.info3D.slice2index,'stable');
           for i=1:length(f)
               info.PerFrameFunctionalGroupsSequence=rmfield(info.PerFrameFunctionalGroupsSequence, f{i});
           end
           for i=1:length(files.info3D.slice2index)
               item=files.info3D.slice2index{i};
               info.PerFrameFunctionalGroupsSequence.(item).Private_2005_140f.Item_1.InstanceNumber=1;
               info.PerFrameFunctionalGroupsSequence.(item).PixelValueTransformationSequence.Item_1.RescaleIntercept=0;
               info.PerFrameFunctionalGroupsSequence.(item).Private_2005_140f.Item_1.RescaleIntercept=0;
               info.PerFrameFunctionalGroupsSequence.(item).PixelValueTransformationSequence.Item_1.RescaleSlope=1;
               info.PerFrameFunctionalGroupsSequence.(item).Private_2005_140f.Item_1.RescaleSlope=1;
               info.PerFrameFunctionalGroupsSequence.(item).FrameVOILUTSequence.Item_1.WindowCenter=opts.Convert(opts.Level);
               info.PerFrameFunctionalGroupsSequence.(item).FrameVOILUTSequence.Item_1.WindowWidth=opts.Convert(opts.Window);   
           end
           info.NumberOfFrames=length(files.info3D.slice2index);
           sz=size(M);
           M=reshape(opts.Convert(M), sz(1), sz(2), 1, []);
           M=permute(M, [2 1 3 4]);
           if isfield(files, 'slices_added')
               if files.slices_added
                   warning('Removing empty slice at bottom of volume');
                   M=M(:,:,1:end-1);
               end
           end
           %dicomwrite(M,outfile, ...
           %    'CreateMode', 'copy', 'WritePrivate', true, info);
           dicomwrite(M,outfile, ...
               'WritePrivate', true, info);
        else
            if isfield(files, 'slices_added')
                if files.slices_added
                    warning('Removing empty slice at bottom of volume');
                    M=M(:,:,1:end-1);
                end
            end
            %dicomwrite(opts.Convert(M(:,:,slice))',outfile, ...
            %    'CreateMode', 'copy', 'WritePrivate', true, info);
            %
            dicomwrite(opts.Convert(M(:,:,slice))',outfile, ...
                'WritePrivate', true, info);
        end
    end
end
fprintf('\n');


    function print_progress(arg)
        num=length(progress);
        num=num-numel(regexp(progress, '\\\\'));
        for ii=1:num; fprintf('\b'); end
        progress=['Writing file ' arg];
        progress=regexprep(progress,'\','\\\\');
        fprintf(progress);
    end

    function y=convert(x)
        if flag_signed
            y=int16(x*1000);
        else
            y=uint16(x*1000);
        end
    end
end

2.4. QSM_3Dを集める

各フォルダには、強度画像と位相画像が入っている。
TEを7回変えて撮像している。
国際脳QSMの撮像は、1回撮像に128 slicesなので強度画像と位相画像はそれぞれ、896枚(=128 slices ×7 phase)

(Cornell MRI Research Lab. QSMより、URL: http://pre.weill.cornell.edu/mri/pages/qsm.html)

QSM_3Dにある強度画像と位相画像をanalyzeフォルダのrawdata(※必ずrawdata)にまとめて保存。

3. 実行

3.1. 作業フォルダ(analyzeフォルダ)へ移動

MATLABを起動し、analyzeフォルダへ移動。

3.2. runme.mを実行

runme.mをMATLABへDrag & Dropし、QSM解析を実行。

3.3. Outputファイルの確認

runme.mを実行後、

  • QSM
  • results

フォルダが生成される。

QSMフォルダにQSM画像がDICOM形式で保存される。

(Cornell MRI Research Lab. QSMより、URL: http://pre.weill.cornell.edu/mri/pages/qsm.html)

4. 複数の被験者データをまとめて処理したい場合

被験者ごとにフォルダを作成し一つのフォルダにまとめます。

さらに、後で紹介するまとめてMEDIを実行するスクリプトも入れておきます。

各被験者フォルダには強度画像と位相画像のDICOMがまとめて入ったrawdataフォルダがあります。

以下のスクリプトを実行すると、すべての被験者のQSM mapが計算できます。

4.1. run_sequential.m

clear all;clc;close all;

% change direct to study folder
selpath = uigetdir('Select the folder including all subject');
cd(selpath);
% list folder
folderlist = dir(selpath);
folderlist = folderlist(~ismember({folderlist.name}, {'.', '..'}));
folderlist = folderlist([folderlist.isdir]);

% run MEDI
for i = 1:length(folderlist)
    disp(['processing '  folderlist(i).name ' ...'])
    main(folderlist(i).name, selpath);
end

% define function of MEDI processings
function main(folder, basepath)
    % initialize
    clearvars -except selpath folderlist folder basepath i
    % move to subject folder
    cd(folder)

    % STEP 1: Import data
    [iField,voxel_size,matrix_size,CF,delta_TE,TE,B0_dir,files]=Read_DICOM('rawdata');

    %%%% Generate the Magnitude image %%%%
    iMag = sqrt(sum(abs(iField).^2,4));

    %%%%% provide a Mask here if possible %%%%%%
    if (~exist('Mask','var'))                     
    %     Mask = genMask(iField, voxel_size);
        Mask = BET(iMag,matrix_size,voxel_size);
    end

    %%%%% provide a noise_level here if possible %%%%%%
    if (~exist('noise_level','var'))
        noise_level = calfieldnoise(iField, Mask);
    end

    %%%%% normalize signal intensity by noise to get SNR %%%
    iField = iField/noise_level;

    % STEP 2a: Field Map Estimation
    %%%%%Estimate the frequency offset in each of the voxel using a 
    %%%%%complex fitting %%%%
    [iFreq_raw N_std] = Fit_ppm_complex(iField);

    % STEP 2b: Spatial phase unwrapping %%%%
    iFreq = unwrapPhase(iMag, iFreq_raw, matrix_size);

    % STEP 2c: Background Field Removal
    %%%% Background field removal %%%%
    [RDF shim] = PDF(iFreq, N_std, Mask,matrix_size,voxel_size, B0_dir);

    % CSF Mask for zero referencing
    R2s=arlo(TE,abs(iField));
    Mask_CSF = extract_CSF(R2s, Mask, voxel_size);

    % STEP 3: Dipole Inversion
    save RDF.mat RDF iFreq iFreq_raw iMag N_std Mask matrix_size...
        voxel_size delta_TE CF B0_dir Mask_CSF;

    %%%% run MEDI %%%%%
    QSM = MEDI_L1('lambda',1000, 'smv', 5, 'merit');

    %%% save to DICOM
    % ignore warnings...
    Write_DICOM(QSM,files,'QSM')

    % back to study folder
    cd(basepath)
end

macOSでのSPM12のコンパイル方法

ある方から、Apple M1のmacでSPMを起動しようとするとspm_check_installation(‘basic’)でエラーが出て起動しないという相談を受けました。

コンパイルしたら問題は解決しました。コンパイル方法を共有します。

ただし、その後、SPMのMLでこのディスカッションに乗ってみたところ、コンパイルは不要だよということも教えていただきました。なので、コンパイルに挑戦してみたい人向けと思ってください。(普通は不要です)

続きを読む

AC-PC自動設定スクリプト(SPM12対応版)

かつてauto_reorient.mというスクリプトを配布していましたが、これは現在のSPM12で動かなくなってしまいました。理由は単純でSPM12からspm_affregという機能がなくなってしまったからです。もし、過去のSPMからspm_affreg.mを持ってきたら問題なく使えるのですが、この機会に別の手法を考えてみました。

シンプルな方法は、MNI標準脳にCo-registrationすることです。これだけでかなりあいます。
しかし、画像の原点があまりにも違うところに設定されているとエラーが出ることがあります。
そこで、以前、山下先生に教わった方法を採用し、まず、originを画像の中心に設定し、そのうえで、SPM12に搭載されているicbm152.niiにco-registartionするスクリプトを書いてみました。

多くの画像で試してみましたが、それなりにうまくいきますし、処理速度も速いです。

よかったら試してみてください。

acpc_coreg.mをダウンロード(右クリックで保存してください)

Matlabのパスが通っているフォルダにこのファイルを置いていただき、

Matlabから

acpc_coreg

とタイプするだけです。

BashからMatlabスクリプトを実行する方法

先日、ある方と「BashからMatlabを呼び出せないだろうか」という話をしていました。もし、これができたら、シェルスクリプトから、Matlabを呼び出せるので、シェルとMatlabを完全に連携できるわけです。

結論としては、以下でできました。

  • Short answer
  • Matlabのスクリプト名を sample_code.m とすると、以下でできます。

    $ matlab -nodesktop -nosplash -r 'sample_code; exit'
    

    コツは2つです。

  • スクリプト名ではなく、コマンドとして指示するため、.mは外す
  • Matlabから抜けるために exit を追加する

続きを読む

Matlab R2016a日本語版では、SPM12を起動するときにエラーが出る

タイトル通りなのですが、

Matlab R2016a日本語版では、SPM12を起動しようとするとエラーが出ます。

詳しい現象は、こちらに説明されていますので、下記をご覧ください。

SPM12がMATLAB R2016aで動かない

こちらに解決法が書いてあり、「デスクトップの言語を英語にする」と書いてあります。

Linuxではこれは簡単です。

ターミナルを起動し、

LANG=C matlab &

で英語環境でMatlabが起動します。

Macで同じことができるかと思ったらダメでした。
Macでは、ターミナルでLANG変数を変更しても、Matlab起動時に無視されてしまうとのこと。Mathworks本家が言っているのでどうしようもないです。

Mac プラットフォームのロケール設定

決して賢くはないですが、なんとかする方法がわかりましたので書いておきます。

Matlabを起動する前に、

システム環境設定→言語と地域

で、優先する言語を”English”にします。

そして、そのウィンドウを閉じようとすると、再起動しますかと聞かれますが、再起動しなくてOKです。

そこで、Matlabを起動すると、英語モードで起動します。

そうすると、SPMは問題なく起動します。

終わったあとに、もう一度 System Preferences… で同様にJapaneseを優先言語にもってきます。

これで終了すれば、日本語に戻ります。

もっとスマートな方法があるかと思いましたが、ありませんでした…。

今後、時間があるときに、SPMのMailing listにバグ報告を出しておこうと思いますが、
一番の対策は、Matlab R2015b以降にはしばらくアップデートしないということになります。

ご参考まで。

A Matlab script to generate ROI masks using an Atlas in SPM12

SPM12 introduces some useful functions such as spm_atlas or new atlas “labels_Neuromorphometrics.” We find the description about labels_Neuromorphometrics in SPM12 Release note.

Maximum probability tissue labels derived from the “MICCAI 2012 Grand Challenge and Workshop on Multi-Atlas Labeling” are available in files tpm/labels Neuromorphometrics.{nii,xml}. These
data were released under the Creative Commons Attribution-NonCommercial (CC BY-NC) with no end date. Users should credit the MRI scans as originating from the OASIS project and the labeled
data as “provided by Neuromorphometrics, Inc. under academic subscription”. These references should be included in all workshop and final publications. See spm templates.man for more details about the generation of this file.

I wanted to generate masks of some regions using this labels_Neuromorphometrics.

Below is the tiny script which generates masks from your preferred atlas.
Running script brings up a file selector. You can choose any atlas you want.
Then it brings up another dialog which lists the region within the atlas. You can choose as many regions as you want, and the scripts generates masks whose file name is the name of the regions.


%generate_masks_from_atlas.m
%This script generate mask files from any atlases you prefer.
%K. Nemoto 25 April 2015

xA=spm_atlas('load');
S=spm_atlas('select',xA);

for i = 1:size(S,2)
    fname=strcat(S{i},'.nii');
    VM=spm_atlas('mask',xA,S{i});
    VM.fname=fname;
    spm_write_vol(VM,spm_read_vols(VM));
end

Download generate_masks_from_atlas.m (right click and save as)

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

2019.10.11: すべての回答ができましたので、アップデートしました。

何人かの方々から、「心理のためのMatlabチュートリアル」の練習の解答がないか問い合わせを受けていました。
回答の完全版ができましたので、公開させていただきます。

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