tbss_non_FA discrepancy

[psharp1289]

I believe I may have found a bug in the create_tbss_non_fa() dMRI workflow, or I was doing something completely wrong :)

When I ran this, and compared it with the FSL tbss_non_fa command, I kept getting SIGNIFICANT discrepancies. So, I went into the source code and compared the "project radial diffusivity onto mean FA map (tbss_skeleton)" step in FSL with the projectfa node supposedly wrapping this same command in Nipype. However, when I compared the command.txt file for this step in Nipype with the source code in FSL's tbss_non_fa, there was a mismatch. Essentially, in FSL, it always uses the "'a" option for alternative skeleton, and inputs the 4D all_RD file, and under the -p option, it has all_FA as the data_file.

So, when I went into the source code for nipype, and added an "all_FA_file" input into the infosource node, and switched the workflow.connect piece to have the all_RD (i believe its maskedgroup) output from tbss_3 as the 'alt_data_file' input to projectfa (the wrapper for tbss_skeleton), the discrepancies disappeared!

Please let me know anyone out there if a) this makes sense and b) if this is worth changing on nipype. Also, I should say that before I began tampering with the tbss.py sourcecode, the tbss_all() on my FA files worked completely fine. In other words, when I compared the tbss_all FA output to the regular FSL tbss Pipeline, the all_FA_skeletonised files were identical. Such was not the case for all_RD_skeletonised files, and the values generated by FSLs command compared to nipype's made more sense. (nipype output significantly higher RD values).

Paul

[psharp1289]

Here is the FSL command from their source code: "${FSLDIR}/bin/tbss_skeleton -i mean_FA -p $thresh mean_FA_skeleton_mask_dst ${FSLDIR}/data/standard/LowerCingulum_1mm allFA all${ALTIM}skeletonised -a all$ALTIM"

where all_${ALTIM} in my example is all_RD. notice all_FA is in the placeholder for the 4d Data file under the -p option, and the -a all_RD is there too.

here is what I changed in tbss.py to allow for this:

1. see new inputnode, all_FA_file and
   1. see (maskgroup, projectfa, [('out_file', 'alt_data_file')]) step in tbss_non_fa.connect command. This is to get the -a option with the all_RD file.

def create_tbss_non_FA(name='tbss_non_FA'): """ A pipeline that implement tbss_non_FA in FSL

Example
-------

>>> from nipype.workflows.dmri.fsl import tbss
>>> tbss_MD = tbss.create_tbss_non_FA()
>>> tbss_MD.inputs.inputnode.file_list = []
>>> tbss_MD.inputs.inputnode.field_list = []
>>> tbss_MD.inputs.inputnode.skeleton_thresh = 0.2
>>> tbss_MD.inputs.inputnode.groupmask = './xxx'
>>> tbss_MD.inputs.inputnode.meanfa_file = './xxx'
>>> tbss_MD.inputs.inputnode.distance_map = []
>>> tbss_MD.inputs.inputnode.all_FA_file = './xxx'

Inputs::

    inputnode.file_list
    inputnode.field_list
    inputnode.skeleton_thresh
    inputnode.groupmask
    inputnode.meanfa_file
    inputnode.distance_map
    inputnode.all_FA_file

Outputs::

    outputnode.projected_nonFA_file

"""

# Define the inputnode
inputnode = pe.Node(interface=util.IdentityInterface(fields=['file_list',
                                                             'field_list',
                                                             'skeleton_thresh',
                                                             'groupmask',
                                                             'meanfa_file',
                                                             'distance_map',
                                                             'all_FA_file']),
                    name='inputnode')

# Apply the warpfield to the non FA image
applywarp = pe.MapNode(interface=fsl.ApplyWarp(),
                       iterfield=['in_file', 'field_file'],
                       name="applywarp")
if fsl.no_fsl():
    warn('NO FSL found')
else:
    applywarp.inputs.ref_file = fsl.Info.standard_image("FMRIB58_FA_1mm.nii.gz")
# Merge the non FA files into a 4D file
merge = pe.Node(fsl.Merge(dimension="t"), name="merge")
#merged_file="all_FA.nii.gz"
maskgroup = pe.Node(fsl.ImageMaths(op_string="-mas",
                                   suffix="_masked"),
                    name="maskgroup")
projectfa = pe.Node(fsl.TractSkeleton(project_data=True,
                                    #projected_data = 'test.nii.gz',
                                    use_cingulum_mask=True
                                  ),
                    name="projectfa")

tbss_non_FA = pe.Workflow(name=name)
tbss_non_FA.connect([
                (inputnode, applywarp, [('file_list', 'in_file'),
                                        ('field_list', 'field_file'),
                                        ]),
                (applywarp, merge, [("out_file", "in_files")]),

                (merge, maskgroup, [("merged_file", "in_file")]),

                (inputnode, maskgroup, [('groupmask', 'in_file2')]),

                (maskgroup, projectfa, [('out_file', 'alt_data_file')]),
                (inputnode, projectfa, [('skeleton_thresh', 'threshold'),
                                        ("meanfa_file", "in_file"),
                                        ("distance_map", "distance_map"),
                                        ("all_FA_file", 'data_file')
                                        ]),
            ])

# Define the outputnode
outputnode = pe.Node(interface=util.IdentityInterface(
                                        fields=['projected_nonFA_file']),
                     name='outputnode')
tbss_non_FA.connect([
        (projectfa, outputnode, [('projected_data', 'projected_nonFA_file'),
                                ]),
        ])
return tbss_non_FA

Paul

[satra]

@psharp1289 - a pull request with the changes would be much appreciated.

(please see https://github.com/nipy/nipype/blob/master/CONTRIBUTING.md for more info)

[psharp1289]

ok, great. I've only worked with github minimally in a previous python course. I will try it!