基于ANTs and FSL 的預(yù)處理腳本

最近在寫一個toolbox箱玷，所以嘗試了各種主流的預(yù)處理方法和算法歧蕉，通過不斷的嘗試夫嗓，的出來了下面的這個預(yù)處理的pipeline比較好用裹唆，這里記錄一下推薦給大家究珊，主要的代碼是用python實現(xiàn)，但是其實還是調(diào)用的ANTs跟fsl,核心代碼如下所示

import os
from multiprocessing import Pool as Pool
Standard_image_T1w = '../Template/MNI152_T1_1mm_brain.nii.gz'
Standard_t1w='../Template/MNI152_T1_1mm.nii.gz'
Standard_mask = '../Template/MNI152_T1_1mm_brain_mask.nii.gz'
#Not BIDS format
def T1w_prep(input_dir,output_dir,sub):
    #T1_directory = os.path.join(input_dir,'anat','')
    os.system('mkdir '+output_dir+'tmp/')
    os.system('mkdir '+output_dir+'tmp/T1/')
    os.system('mkdir '+output_dir+'T1')
    os.system('mkdir '+output_dir+'Seg')
    T1_files = os.listdir(input_dir)
    for i in range(len(T1_files)):
        if T1_files[i].find('T1w.nii') != -1 :
            T1_name = T1_files[i]
            T1_file = os.path.join(input_dir,T1_files[i])
    #preprocess T1 image竿裂，核心代碼如下
    #Crop
        if os.path.exists(output_dir+'tmp/T1/crop_'+T1_name) == 0:
            os.system('robustfov -i '+T1_file+' -r '+output_dir+'tmp/T1/crop_'+T1_name)
        #Skullstrip
        if os.path.exists(output_dir+'tmp/T1/BEBrainExtractionBrain.nii.gz') == 0:
            if T1_name.find('.gz'):
                os.system('antsBrainExtraction.sh -d 3 -a '+output_dir+'tmp/T1/crop_'+T1_name+' -e '+Standard_t1w+' -m '+Standard_mask+' -o '+output_dir+'tmp/T1/BE')
            else:
                os.system('antsBrainExtraction.sh -d 3 -a '+output_dir+'tmp/T1/crop_'+T1_name+'.gz'+' -e '+Standard_t1w+' -m '+Standard_mask+' -o '+output_dir+'tmp/T1/BE')

        #Registration
        if os.path.exists(output_dir+'T1/T1w2std_'+T1_name) == 0:
            os.system('antsRegistrationSyN.sh -d 3 -f '+Standard_image_T1w+' -m '+output_dir+'tmp/T1/BEBrainExtractionBrain.nii.gz  -o ' \
                +output_dir+'tmp/T1/rega2t -n 20')
            os.system('antsApplyTransforms -d 3 -i '+output_dir+'tmp/T1/crop_'+T1_name+' -o '+output_dir+'T1/T1w2std_'+T1_name+' -r '+Standard_t1w+ \
                '-t '+output_dir+'tmp/T1/rega2t1Warp.nii.gz -t '+output_dir+'tmp/T1/rega2t0GenericAffine.mat')

        #Segmentation
        if os.path.exists(output_dir+'Seg/'+sub+'_CSF_mask.nii.gz') == 0 or os.path.exists(output_dir+'Seg/'+sub+'_WM_mask.nii.gz') == 0:
            os.system('fast -S 1 -t 1 -o '+output_dir+'Seg/seg -g -n 3 -b -I 10 '+output_dir+'tmp/T1/rega2tWarped.nii.gz')
            os.system('mv '+output_dir+'Seg/seg_seg_0.nii.gz '+output_dir+'Seg/'+sub+'_CSF_mask.nii.gz')
            os.system('mv '+output_dir+'Seg/seg_seg_2.nii.gz '+output_dir+'Seg/'+sub+'_WM_mask.nii.gz')
            #os.system('mv '+output_dir+'Seg/seg_seg_1.nii.gz '+output_dir+'Seg/'+sub+'_GM_mask.nii.gz') # mask 灰質(zhì)的一般不需要

唯一的缺點是玉吁，在進(jìn)行配準(zhǔn)的時候是基于迭代優(yōu)化求解的，所以需要花費(fèi)比較多的時間腻异，但是效果還是比直接用cat12做出來的好很多的进副。當(dāng)然，這個用python調(diào)用也是支持并行的，并行代碼可以參考下面

pool = Pool(processes = n)
for i in range(num,num+n):
    sub = sublist[i]
    pool.apply_async(prep_main,(sub,))# prep_main 是封裝的一個預(yù)處理的函數(shù)影斑，參考一下寫法就可以
pool.close()
pool.join()