WISHLIST: Support CIVET-style surfaces and measures

[gdevenyi]

The MNI object format is a simple list of vertices (for this application)

http://www.bic.mni.mcgill.ca/users/mishkin/mni_obj_format.pdf

And the thickness/area measures are a flat text file in the same order as vertices.

[trislett]

Hi Gabriel,

Including support for CIVET is definitely a feature that I'm planning to add. Unfortunately, I'm having trouble downloading CIVET or even obtaining an MNI object file (can you supply one?). The TFCE algorithm is generic and it can run on any datatype providing that an adjacency set is defined.

If you are eager for CIVET support, I would suggest converting the MNI obj to ascii file, and then convert the ascii file to a freesurfer surface. Next, edit the tfce_mediation/tools/create_adjacency_list.py to read in the new freesurfer image, and project the distance by -.5 to the midthickness surface (assuming the mni obj is projected to the outer surface and not the white matter surface). TFCE should now be possible, and the I/O would just need to be slightly adjusted for reading/writing area or thickness text files.

Do you know how to read the MNI obj files into python/numpy? If so, I could directly set the vertices and faces information to create_adjacent_list and then upload the adjacency set numpy object for a 3mm geodesic distance.

Tris

[gdevenyi]

The MNI object format is very simple (see documentation link above). If you're only interested in the vertices, the first line is info on the number of vertices and then it's just a list of XYZ coordinates (everything past that is extra objects, colours, normals, etc)

It is supported by VTK, and as such pyVTK, which I've previously used to manipulate surfaces. http://www.vtk.org/doc/nightly/html/classvtkMNIObjectReader.html

Find at https://www.dropbox.com/sh/n2k0ivsedf99cnj/AABicFU_WGx9Vw8PmB0qpPWla?dl=0 left-right mid-surface for an OASIS subject, along with the the thickness measurement. Both are already resampled into the common MNI surface space.

[trislett]

EDIT: Outdated. See the final post.

Okay, I've made preliminary support for CIVET without massively over hauling my I/O.

1. The first step is to convert the MNI object to a freesurfer surface because I need to be working with triangles. Maybe Lea has a better solution. I had to use afni (I know...) to convert *obj file, then mris_convert to convert the ascii file to freesurfer formatted surface.

e.g.,

${AFNI_DIR}/ConvertSurface -i_mni  OASIS_OAS1_0001_MR1_mid_surface_rsl_left_81920.obj -make_consistent -o_fs OASIS_MIDTHICK_left
mris_convert OASIS_MIDTHICK_left.asc lh.OASIS_MIDTHICK
${AFNI_DIR}/ConvertSurface -i_mni OASIS_OAS1_0001_MR1_mid_surface_rsl_right_81920.obj -make_consistent -o_fs OASIS_MIDTHICK_right
mris_convert OASIS_MIDTHICK_right.asc rh.OASIS_MIDTHICK

2. The next step is to create the adjacency sets for the MNI surface using the script I have just push to misc_scripts. Note, I've removed the surface projection from the create_adjacency script.

e.g., for adjacency sets ranging from geodesic distances 1-4 mm:

$TFCE_MEDIATION_DIR/tfce_mediation/misc_scripts/create_adjacency_list_mni.py -d 1 4 --surfaces lh.OASIS_MIDTHICK rh.OASIS_MIDTHICK

3. Convert the cortical thickness text file to an mgh. This is necessary until a rewrite my surface I/O to accommodate text files. You can easily convert the mgh back to text file if you want to use brainview2 to see your results. Otherwise, you can visualize your results in freeview. For group analysis, you will need to paste each subjects surface values into a csv formatted file. TFCE_mediation assumes freesurfer output so you will rename the mgh files to the following convention ?h.all.{surface}.03B.mgh.

$TFCE_MEDIATION_DIR/tfce_mediation/tfce_mediation/misc_scripts/tm_mgh_txt_converter.py -i lh_cortical_thickness_for_all_subjects.csv
mv lh_cortical_thickness_for_all_subjects.csv.mgh lh.all.thickness.03B.mgh
$TFCE_MEDIATION_DIR/tfce_mediation/tfce_mediation/misc_scripts/tm_mgh_txt_converter.py -i rh_cortical_thickness_for_all_subjects.csv
mv rh_cortical_thickness_for_all_subjects.csv.mgh rh.all.thickness.03B.mgh

4. Run TFCE_mediation as normal; however, make sure use the -a option (Load custom adjacency set) for 'step1' analyses.

[trislett]

Also, one caveat. This is very experimental, and I haven't validated anything. That said, I would highly recommend zeroing out the surface connecting the hemispheres prior to any statistical analyses. tm_maths might be useful here...

[gdevenyi]

Oh wow, that was quick. Awesome.

I'll need to absorb the freesurfer-style and methods a bit, since I have zero experience there.

Is there a complete freesurfer example dataset so I can see how to input regressors and such?

[gdevenyi]

Random note: As of right now, you MNI commit has a commit time of "12 minutes from now" which I think means your local clock is off? Or maybe your timezone setting :)

[trislett]

Freesurfer comes with some examples, but they are limited. Also, the HCP extended structural data contains the recon-all outputs for each subject. In general, inputting regressors ect. is quite easy with TFCE_mediation. The readme page has an example. Basically, it is just loading a csv file or two csv files including a covariate file.

re: clock time. Yeah, it is something annoying that I haven't fixed yet. Every system that I install ubuntu 16.04 seems to have this issue. The system time keeps resetting, so I suspect it is something to do with the bios.

[trislett]

The previous instructions are now outdated. tm_tools now includes a script for converting surfaces.

1. Convert the MNI obj polygons to freesurfer triangular mesh

   tm_tools convert-surface -i_mni OASIS_OAS1_0001_MR1_mid_surface_rsl_left_81920.obj -o_fs OASIS_OAS1_0001_MR1_mid_surface_rsl_left_81920.srf
   tm_tools convert-surface -i_mni OASIS_OAS1_0001_MR1_mid_surface_rsl_right_81920.obj -o_fs OASIS_OAS1_0001_MR1_mid_surface_rsl_right_81920.srf

2. Create adjacency set for newly converted surfaces to be used later for applying surface based TFCE (~3-4 min)

$TFCE_MEDIATION_DIR/tfce_mediation/misc_scripts/create_adjacency_list_mni.py -d 1 4 --surfaces OASIS_OAS1_0001_MR1_mid_surface_rsl_left_81920.srf OASIS_OAS1_0001_MR1_mid_surface_rsl_right_81920.srf

3. Convert the cortical thickness text file to an mgh. For group analysis, you will need to paste each subjects surface values into a csv formatted file. TFCE_mediation assumes freesurfer output so you will rename the mgh files to the following convention ?h.all.{surface}.03B.mgh.

$TFCE_MEDIATION_DIR/tfce_mediation/misc_scripts/tm_mgh_txt_converter.py -i lh_cortical_thickness_for_all_subjects.csv
mv lh_cortical_thickness_for_all_subjects.csv.mgh lh.all.thickness.03B.mgh
$TFCE_MEDIATION_DIR/tfce_mediation/misc_scripts/tm_mgh_txt_converter.py -i rh_cortical_thickness_for_all_subjects.csv
mv rh_cortical_thickness_for_all_subjects.csv.mgh rh.all.thickness.03B.mgh

4. Run TFCE_mediation as normal; however, make sure use the -a option (Load custom adjacency set) for 'step1' analyses.

5. After statistical analysis, *mgh in the output folder can converted back to text file using 'tm_mgh_txt_converter.py -s'. Otherwise, the results can be visualized in freeview.

[vsteiger]

hey @trislett this is great, thank you very much.

Although, I'm running into an error...

tm_tools convert-surface -i_mni new_lefthc.obj -o_fs new_lefthc.srf

tm_tools: error: invalid choice: 'convert-surface' (choose from 'regressor-tools', 'plot-permutations', 'voxel-cluster-step1', 'voxel-cluster-step2', 'voxel-extract-mean-values', 'voxel-npy-to-nifti', 'voxel-make-subsample', 'vertex-freeview-quick', 'vertex-cluster-step1', 'vertex-cluster-step2', 'vertex-box-cox-transform', 'vertex-create-adjac-list', 'vertex-make-subsample', 'vertex-extract-mean-values', 'vertex-threshold-image', 'merge-images', 'remove-components', 'geodesic-fwhm')

[trislett]

Hi @vsteiger , It seems that you don't have the latest version installed. You need to pull the latest version from the repository and install it, or wait until I upload v1.1.1 to pip (in a few moments).

edit: the latest version is on pip now

[vsteiger]

@trislett , I pipped the "old" version...now, it works fine. Thx

[vsteiger]

@trislett for the group analysis... would then the "lh_cortical_thickness_for_all_subjects.csv" contain values of thickness for every subject in separate columns?

[trislett]

@vsteiger Yes. It has to be a comma separated file with each column representing a subject, and row representing the vertices. Make sure to convert the csv to an mgh using tfce_mediation/misc_scripts/tm_mgh_txt_converter.py.

[vsteiger]

I guess something with the .csv is wrong, since I'm running into an error:

$tfce_mediation step1-vertex-regress -r dudo_predictors.csv -s area -f 03B -a lh_adjacency_dist_3.0_mm.npy rh_adjacency_dist_3.0_mm.npy Loading prior adjacency set Traceback (most recent call last): File "/usr/local/bin/tfce_mediation", line 130, in <module> args.func(args) File "/usr/local/lib/python2.7/dist-packages/tfce_mediation/tmanalysis/STEP_1_vertex_tfce_multiple_regression.py", line 243, in run invXX = np.linalg.inv(np.dot(X.T, X)) File "/usr/local/lib/python2.7/dist-packages/numpy/linalg/linalg.py", line 526, in inv ainv = _umath_linalg.inv(a, signature=signature, extobj=extobj) File "/usr/local/lib/python2.7/dist-packages/numpy/linalg/linalg.py", line 90, in _raise_linalgerror_singular raise LinAlgError("Singular matrix") numpy.linalg.linalg.LinAlgError: Singular matrix

[trislett]

It seems like you either a column of ones, or one of your columns is equal to another is your regressors file. Could you please attach dudo_predictors.csv so I can check it?

[vsteiger]

column1: Group 1 or 2 column2: Group 1, dummies column3: Group 2, dummies

dudo_predictors.txt

[trislett]

Edit: Sorry, I misread your file. For dummy coding, you have one less column than the number of groups. Therefore, in your file you should have just one column.

i.e.,

1 1 1 1 0 0 0 0

[vsteiger]

hehe, I came across this solution but this resulted in:

Loading prior adjacency set
/usr/local/lib/python2.7/dist-packages/tfce_mediation/pyfunc.py:71: VisibleDeprecationWarning: boolean index did not match indexed array along dimension 0; dimension is 1152 but corresponding boolean dimension is 1215
  vertStat_out[bin_mask] = vertStat
Traceback (most recent call last):
  File "/usr/local/bin/tfce_mediation", line 130, in <module>
    args.func(args)
  File "/usr/local/lib/python2.7/dist-packages/tfce_mediation/tmanalysis/STEP_1_vertex_tfce_multiple_regression.py", line 254, in run
    write_vertStat_img('tstat_con%d' % tnum, tvals[tnum,num_vertex_lh:], outdata_mask_rh, affine_mask_rh, surface, 'rh', bin_mask_rh, calcTFCE_rh, all_vertex)
  File "/usr/local/lib/python2.7/dist-packages/tfce_mediation/pyfunc.py", line 71, in write_vertStat_img
    vertStat_out[bin_mask] = vertStat
IndexError: index 1152 is out of bounds for axis 0 with size 1152

The two surfaces do not contain the same number of vertices left 1152 and right 1215...that might be a problem....freesurfer's mri_info revealed that lh.all.thickness.03B.mgh has correct dimensions( 10x1x1x1152).

[gdevenyi]

Hehe @vsteiger is totally trying to do TFCE on MAGeTbrain surface outputs :)

[trislett]

@gdevenyi It seems like it. That's about the right number of vertices for the hippocampi.

@vsteiger it shouldn't matter if there are a different number of vertices. Could you please email me (tris.lett{at}gmail.com) the MNI object files that you trying TFCE_mediation on so I can bug test them? I don't need the SA or CT values.

Do you get mean values for the left and right hemisphere? Also, do you get the tstat output for the left hemisphere?

Last, lh.all.thickness.03B.mgh does not have the correct dimension. It should be of the shape [NumVertices,0,0,NumSubjects] or in your case 1152x1x1x10. Did you transpose your thickness values?

[trislett]

It was a small index issue caused by the hemispheres having a different number of vertices. It should now been fixed in v1.1.2.