MRtrix3 预处理和纤维追踪全过程

step 1 预处理

1.1转格式

将nifti原始图像，bvec,bval等文件打包成mif格式

#convert nifti to native format mif
	mrconvert -fslgrad dwi.bvec dwi.bval dwi.nii \
		dwi.mif -force
#creat a mask
	dwi2mask dwi.mif - | maskfilter - dilate preproce_mask.mif -npass 3 -force

1.2 denoise

dwidenoise dwi.mif dwi_denoise.mif -noise noiselevel.mif \
		   -mask preproce_mask.mif -force

1.3 Gibbs Ring artifact

mrdegibbs dwi_denoise.mif dwi_denoise_degibbs.mif -force

1.4 头动和涡流矫正

使用FSL的涡流工具进行扩散图像预处理；如果可能，包括使用 FSL 的 topup 工具进行失真校正

dwifslpreproc dwi_denoise_degibbs.mif dwi_preprocessed.mif \
	-rpe_none -pe_dir AP \
	-eddy_options "--data_is_shelled --slm=linear --niter=5" -force

1.5 场矫正

# bias filed correction
	dwibiascorrect ants dwi_preprocessed.mif biascorr.mif -bias biasfield.mif

1.6 对齐到T1

#aligment to T1w
	flirt.fsl -dof 12 -cost normmi -ref T1w.nii.gz -in b0.nii -omat T_fsl.txt
	transformconvert T_fsl.txt b0.nii T1w.nii.gz flirt_import T_DWItoT1.txt	
	mrtransform -linear T_DWItoT1.txt biascorr.mif dwi_align.mif -force

1.7 使用5ttgen进行脑组织分割

也就是把脑组织分割为皮层灰质，次皮层灰质，白质，CSF和病理组织
5ttgen可以使用的算法包括freesurfer, fsl, gif, hsvs

1. fsl 算法 要求输入的图像必须是 T1-weighted image； 算法会默认进行大脑组织提取（BET），如果已经BET过，需添加-premask 参数
2. freesurfer 算法 要求输入freesurfer分割的图像，任何包括aseg的都行
3. 如果不想软件对影像进行修建可使用 -nocrop 参数
4. 可以使用-nocleanup参数保存过程文件

示例：

# 5tt segmentation
mrconvert aparc.a2009s aseg.mgz aparc.a2009s aseg.nii.gz
5ttgen freesurfer aparc.a2009s aseg.nii.gz 5ttseg.mif -force
5tt2gmwmi 5ttseg.mif 5tt_gmwmi.mif -force

step 2 拟合张量

不需要可以跳过

#creat a mask
dwi2mask align.mif - | maskfilter - dilate align_mask.mif -force
#creat diffusion tensor
dwi2tensor -mask align_mask.mif dwi_align.mif dwi_tensor.mif -force
tensor2metric dwi_tensor.mif -fa dwi_fa.mif -ad dwi_ad.mif \
              	-rd dwi_rd.mif -adc dwi_adc.mif \
              	-cl dwi_cl.mif -force

已经预处理过可以跳过

step 3 CSD方法进行纤维追踪

3.1 预步骤

3.1.1 估计响应函数

首先使用dwi2response命令通过使用指定算法将模型拟合到 DWI 数据来估计响应函数。它以 DWI 图像作为输入，并生成包含不同组织类型的估计响应函数的文本文件（可以后续输入到dwi2fod或tckgen命令中）。
可以选择的算法有6种：

tournier：根据DWI数据的球形去卷积估计单纤维体素的反应函数。
dhollander：基于DWI数据的迭代去噪和去卷积，估计单纤维体素的反应函数。
msmt_5tt：多壳多组织算法基于DWI和解剖图像的组合，估计多种组织类型的响应函数。
fa 用旧的FA-阈值启发法估计单纤维体素的反应函数
manual 仅仅使用输入的mask图像（即预先选择的体素）来推导响应函数。
tax 使用Tax等人（2014）的递归校准算法估计单纤维体素的反应函数

示例：

dwi2response tournier dwi_align.mif wm_response.txt

3.1.2 估计纤维方向分布图（FOD）

dwi2fod通过对DWI数据进行模型拟合，使用指定的算法来估计FOD。它将DWI图像、估计的响应函数作为输入，并输出一个包含估计FOD图像（可以后续输入到tckgen或tcksift2命令中）。

可以选择的算法有三种：

• csd：这种算法是为single-shell数据设计的（只有一个b值)。提供的响应函数文本文件应该只包含一行，对应于用于CSD的b值。

• msmt_csd：multi-shell数据可以使用这个算法，用于估计白质 FOD、灰质和 CSF（至少三个唯一的 b 值来估计三种组织）。如果只有两个唯一的 b 值可用，则也可以仅估计两种组织，例如白质和 CSF。

#estimate fiber orientation distribution
#csd
dwi2fod csd dwi_align.mif wm_response.txt wm_FOD.mif -force
#msmt_csd
dwi2fod msmt_csd dwi_align.mif wm_response.txt wm_fod.mif gm_response.txt gm.mif csf_response.txt csf.mif

3.2 追踪

示例1：

tckgen wm_FOD.mif tck.tck -seed_random_per_voxel 5tt_gmwmi.mif 10 -act 5ttseg.mif -backtrack -crop_at_gmwmi -maxlength 250 -info

• wm_FOD.mif 这里为dwi2fod步骤生成的FOD图
• tck.tck 输出的tractogram文件
• -seed_random_per_voxel 在mask图像中为每个体素植入固定数量的流线
• -act 在追踪过程中，使用解剖约束束缚技术框架；提供的图像必须是5TT（五种组织类型）格式的。
• -backtrack 如果遇到不良的结构终止，允许截断和重新跟踪tracks
• -crop_at_gmwmi 在穿过GM-WM交界时，更精确地裁剪流线的端点
• -maxlength 设置任何tracks的最大长度，单位是毫米（默认：100）。
• -info 呈现过程信息

示例2：

tckgen wm_FOD.mif tck.tck -seed_image dwi_mask.mif -mask dwi_mask.mif -select 2000000

3.3 计算流线权重

追踪之后，对tractogram的流线进行过滤。
MRtrix3中的tcksift2命令是一个用于过滤和减少tractogram中流线数量的工具，同时保留其基本特征。它将tractogram和一个预先计算好的FOD图像作为输入 ，然后通过SIFT2 spherical-deconvolution informed filtering of tractograms）方法来实现这一目的。 SIFT2方法估计大脑每个体素中的纤维方向分布（FOD），然后使用这一信息来权衡tractogram中每条流线的贡献。这样做的结果是减少了流线的数量，同时保留了tractogram的重要特征。这样所得到的tractogram更紧凑，更容易观察到，同时仍然保留了底层白质通路的重要特征，能显著提高tractogram的准确性和可解释性。

相比tcksift，tcksift2不拟合线性模型，而是使用图形套索算法估计沿流线的扩散信号的稀疏反协方差矩阵。这使得对流线权重的估计更加准确，并能产生更高质量的tractogram。然而，tcksift2通常比tcksift慢，因为该算法的计算量更大。

# filtering and reducing the number of streamlines in a tractogram
tcksift2 -force tck.tck wm_FOD.mif sift2.txt -act 5ttseg.mif

输出文件为包含tractogram中每条流线的贡献信息的txt文件（可在tck2connectome作为输入）。

然后使用tckedit得到过滤后的tractogram (这一步可选)

tckedit tck.tck -tck_weights_in sift2.txt sift2_tck.tck

或者你也可以使用tcksift命令对tractogram进行命令，然后基于过滤后的tck文件构建SC。
tcksift命令用于过滤全脑纤维跟踪数据，使流线密度与 FOD lobe integrals相匹配，输出的是过滤后的tractogram文件。它通过对每条流线的扩散信号进行线性模型拟合来估计每条流线的权重。然后，这些权重被用来过滤图谱，从而减少了流线的数量和内存占用。tcksift通常比tcksift2快，但由于线性模型的假设，可能会导致精度的损失。

tcksift tck.tck wm_FOD.mif tck_sift.tck -term_number 100000   

step 4 connectome

使用tck2connectome 构建纤维连接矩阵，至少需要输入纤维追踪的结果文件.tck，以及配准到个体空间的atlas文件

示例：

 tck2connectome -force -symmetric -zero_diagonal -assignment_end_voxels -stat_edge sum -tck_weights_in sift2.txt tck.tck warp_atlas.nii.gz connectome.txt

• -symmetric 这个参数就是生成对称的矩阵，如果不选上，输出的矩阵只有上三角，下三角都是0
• -zero_diagonal 定义对角线为0
• -assignment_end_voxels 每条流线的端点都被分配到包含该端点的体素和网络中最近的节点,这意味着每条流线对体素和节点的连接组都有贡献
• -stat_edge sum 这里是以流线数总和为边，还可以计算mean,min,max
• -tck_weights_in sift2.txt 其中sift2.txt 通过tcksift2这个步骤生成，存储了每条流线对tractogram的贡献。通过这个命令对每条流线的贡献进行加权，减少假阳性流线的数量。
• -tck.tck 输入tractogram
• -warp_atlas.nii.gz 配准到个体空间的atlas/roi
• -connectome.txt 生成的连接矩阵

其它例子：

#Default usage:
tck2connectome tracks.tck nodes.mif connectome.csv -tck_weights_in weights.csv -out_assignments assignments.txt

#Generate a matrix consisting of the mean streamline length between each node pair
tck2connectome tracks.tck nodes.mif distances.csv -scale_length -stat_edge mean

#Generate a connectome matrix where the value of connectivity is the "mean FA"
tcksample tracks.tck FA.mif mean_FA_per_streamline.csv -stat_tck mean
tck2connectome tracks.tck nodes.mif mean_FA_connectome.csv -scale_file mean_FA_per_streamline.csv -stat_edge mean

最后可以用MRview查看生成的connectom文件。

补充

tckgen 支持的追踪算法如下：
（大概翻译了一下，看不懂可以去终端看原文）

- FACT: Fiber Assigned by Continuous Tracking. 确定性纤维追踪算法的一种，将 4D扩散加权图像 (DWI)作为输入。它基于沿扩散主方向的连续追踪。这种方法计算效率高，可以快速产生结果，但它没有考虑到纤维路径方向的不确定性或交叉纤维的存在。
- iFOD1: First-order Integration over Fiber Orientation Distributions. 一种概率追踪算法，它把以球面谐波（Spherical Harmonic SH）为基础的纤维方向分布（Fiber Orientation Distribution，FOD）图像作为输入。在每个流线步骤中，从个体(trilinear interpolated) 的FOD图像中随机取样。一条流线更有可能遵循 FOD幅值大的方向；但它也可能穿过FOD幅值小的方向（很少 ）。
- iFOD2 (default): Second-order Integration over Fiber Orientation Distributions. Mrtrix的默认算法，它iFOD2是 iFOD1的改进版，它使用了更先进的 FOD 估计算法。IFOD2采用了CSD 方法（reduced spherical deconvolution），通过求解球谐域中的linear inverse problem来估计 FOD。这种方法允许更准确地估计交叉纤维区域的 FOD。
- NullDist1 / NullDist2: Null Distribution tracking algorithms. 也是概率算法，这些算法根据随机方向样本生成流线；也就是说，不使用与纤维方向相关的图像信息，流线轨迹完全由随机抽样确定。 NullDist2 算法要与 iFOD2 结合使用； NullDist1 应与任何一阶算法结合使用。
- SD_STREAM: Streamlines tractography based on Spherical Deconvolution (SD). 一种确定性追踪算法，它将以球谐 (SH) 基表示的纤维取向分布 (FOD) 图像作为输入。在每个流线步骤中，通过局部trilinear interpolation 对FOD 进行采样，并根据当前流线切线方向对球体执行Newton优化，以便定位最近的 FOD 幅值峰值的方向。
- SeedTest: A dummy streamlines algorithm used for testing streamline seeding mechanisms. 任何图像都可以用作输入（该图像不会以任何方式使用）。对于由 seeding mechanism生成的每个种子点，包含与该种子位置对应的单个点的流线将被写入输出的track文件。
- Tensor_Det: 一种确定性算法，将 4D 扩散加权图像 (DWI) 系列作为输入。它沿着每一点的扩散张量的主要特征向量的方向追踪纤维。这种方法可以处理交叉的纤维，但当dti拟合得不好的时候，容易出现错误。
- Tensor_Prob: 一种概率算法，将 4D 扩散加权图像 (DWI) 系列作为输入。在每个图像体素内，执行residual bootstrap以获得每个流线的该体素中 DWI 数据的unique realisation。然后在每个流线步骤通过trilinear interpolation对这些数据进行采样，拟合扩散张量模型，并且流线遵循该张量的主要特征向量的方向。

比如若想使用DTI方法进行纤维追踪，可参考：

tckgen dwi_denoised_preproc_bias.mif track.tck -algorithm Tensor_Det -select 20000 -mask dwi_mask.mif -seed_image seed_mask.mif

单纯使用iFOD2可参考

tckgen input.mif output.tck -algorithm iFOD2 -select 20000 -cutoff 0.2 -seed_image seed_mask.mif -mask dwi_mask.mif

• input.mif 是MIF格式的输入弥散MRI数据。
• output.tck是输出的tractogram文件。
• -algorithm iFOD2 指定了使用iFOD2 跟踪算法。
• -select 100000 指定在tractogram中生成的流线的数量。
• -cutoff 0.2 设置FA（分数各向异性）阈值，低于该阈值将停止追踪。
• -seed_image seed_mask.mif 将用于追踪的初始种子点设置为mask文件。
• -mask mask.mif 设置生成的tractogram的最终mask，比如可以以整个大脑为mask，也可以以白质为mask。

这篇好文章是转载于：学新通技术网