这个例子展示如何用FAST 分割三种组织，计算出三种组织的体积，最后使用R语言导入到一个表里并且计算Brain volume
FAST(FMRIB's Automated Segmentation Tool)主要分割全脑灰质GM、白质WM和脑脊液CSF
在运行FAST之前，应该使用BET剥头皮
处理上百个数据会耗费很长时间，为了提高效率，可以利用qsub和PBS将任务提交到cluster上
结构像文件结构如下图所示(设当前路径为pwd)，每一个人的结构像都放在相应文件夹里并且命名为t1.nii
结构说明
如果结构像同级放在一个文件夹下，可以参考Shell script - 处理同一文件夹下的一批数据来修改。
1. bet和fast的使用
分割结构像的脚本tissueSegmentOnSingle.sh存放于pwd路径
segmentOnSingle.sh如下
#!/bin/bash
# PBS -l nodes=1:ppn=8
# PBS -q short
bet t1.nii t1_brain.nii.gz # 剥头皮后生成t1_brain.nii.gz，这一步还算比较快
fast -t 1 -n 3 -H 0.1 -I 4 -l 20.0 -b -o t1_brain.nii t1_brain.nii.gz # fast分割，这一步不会很快
# fast运行结束以后生成
# t1_brain_bias.nii.gz
# t1_brain_mixeltype.nii.gz
# t1_brain_pve_0.nii.gz: CSF （intensity = 属于对应组织的概率）
# t1_brain_pve_1.nii.gz: GM（intensity = 属于对应组织的概率）
# t1_brain_pve_2.nii.gz: WM（intensity = 属于对应组织的概率）
# t1_brain_pveseg.nii.gz
# t1_brain_seg.nii.gz
重点关注的结果(Inspection by freeview)
t1_brain_pve_0.nii.gz: CSF
t1_brain_pve_0.nii.gz: CSF
t1_brain_pve_1.nii.gz: GM
t1_brain_pve_1.nii.gz: GM
t1_brain_pve_2.nii.gz: WM
t1_brain_pve_2.nii.gz: WM
tissueSegmentOnSingle.sh用于处理单个结构像
对于多个数据，应该在循环里多次将该脚本作为作业提交到服务器，如存放在pwd的脚本qsubTissueSegment.sh所示
#!/bin/bash
for n in {1..2..1} ###### 根据你的样本量修改，下面的编号也是如此
do
number=`echo $n|awk '{printf("%04d",$0)}'` # 自动补0
patient=`echo "S"$number` # Modify here (Now patient: S0001)
echo $patient
qsub -d `pwd`/$patient/ `pwd`/qsubTissueSegment.sh
# -d 指明single脚本的工作路径
done
900人左右的数据进入short队列后，5个半小时左右运行结束
2. 输出三种脑组织的体积
现在每个人的文件夹里都有那几个文件了，下面开始利用fslstats将所有人的数据写到一个txt里
fslstats简单例子 参考帖子：
$ fslstats t1_brain_pve_0.nii.gz -M -V
0.738002 394821 394821.000000
# -m : output mean
# -M : output mean (for nonzero voxels)
# -v : output <voxels> <volume>
# -V : output <voxels> <volume> (for nonzero voxels)
忽略所有的零体素（就是脑袋外的体素不算在内)
第一个数字表示整个图像上相应组织的平均体素，一个百分比
第二个是整张图体素数目
第三个数字是图像的总体积(单位：mm³)
将第一和第三个数字相乘可以得到相应组织的总体积。
$ fslstats t1_brain_pve_1.nii.gz -M -V | awk '{print $1 * $3}'
490514
现在可以写出批处理的脚本getTissueVolume.sh了(根据实际情况适当修改)
#!/bin/bash
for n in {1..2..1} ###### Modify here
do
number=`echo $n|awk '{printf("%04d",$0)}'` # 自动补0
patient=`echo "S"$number` # Modify here (Now patient: S0001)
echo $patient
csf=`fslstats $patient/t1_brain_pve_0.nii.gz -M -V | awk '{print $1 * $3}'`
gm=`fslstats $patient/t1_brain_pve_1.nii.gz -M -V | awk '{print $1 * $3}'`
wm=`fslstats $patient/t1_brain_pve_2.nii.gz -M -V | awk '{print $1 * $3}'`
echo -e $csf" "$gm" "$wm >> res # 自动追加换行符
done
fslstats运行速度比较快，平均1.5s就能处理完一个人的数据
运行结果
3.使用R语言整理出全脑体积并导出表格
全脑体积=三种体积之和
original = read.table("~/Share/temp/res.txt")
colnames(original) = c("CSF", "GM", "WM")
L = length(original$CSF)
serial = rep(0, L)
for (i in 1:L) {
serial[i] = paste("S", sprintf("%04d", i), sep = "") # 自动补0
}
result = data.frame(
serial = serial,
CSF = original$CSF,
GM = original$GM,
WM = original$WM,
brainVolume = original$CSF + original$GM + original$WM
)
write.csv(result, "~/Share/result.csv", row.names = F)
全脑体积在1.4×10^6mm³左右上是正常的
网友评论