Hi Matt,
sorry to come back on this again but I am mainly trying to also
understand the CIFTI format, as I am also interested by formats to be
not only readable by closed source (even if freely distributed) binary
but could be read by others software which is reason why I spent hours
trying to reverse-engineer it (and finding bugs it seems but the
reason you mentionned is not the one causing problems in my analysis).
Now that I have each hemisphere timeseries stored in gifti separately
I was able to compare to what is stored in the cifti.
Loading the cifti in python with nibabel gave me an array of size
(1,1,1,1,91282,1200) and from the gifti file I got the data from left
cortex as time frame in darrays. Here is what I have done in ipython
to understand why my analysis results are wrong:
In [92]:
nii=nb.load('../REST/REST_fix/101915/MNINonLinear/Results/rfMRI_REST1_LR/rfMRI_REST1_LR_Atlas_hp2000_clean.dtseries.nii')
In [93]: left_ts=nb.gifti.read('101915.left_ts.gii')
In [94]: nii.shape
Out[94]: (1, 1, 1, 1, 91282, 1200)
In [95]: len(left_ts.darrays)
Out[95]: 1200
In [96]: left_ts.darrays[0].data.shape
Out[96]: (32492,)
In [97]: from lxml import etree
In [98]: cifti = etree.XML(nii.header.extensions[0].get_content())
In [99]: cifti[0][2][0].items()
Out[99]:
[('IndexOffset', '0'),
('IndexCount', '29696'),
('ModelType', 'CIFTI_MODEL_TYPE_SURFACE'),
('BrainStructure', 'CIFTI_STRUCTURE_CORTEX_LEFT'),
('SurfaceNumberOfNodes', '32492')]
In [100]: cifti[0][2][0][0].text[:100]
Out[100]: '0 1 2 3 4 5 6 8 9 10 11 12 13 14 15 16 17 18 19 20 68 69 70
71 72 73 74 75 76 77 78 79 80 81 82 83 8'
In [113]: idxs=np.array(cifti[0][2][0][0].text.split(' ')).astype(np.int)
# so ok until now it seems that the 29696 first columns of the cifti
file should be the one from the left cortex grayordinates.
In [102]: left_tss = np.c_[[d.data for d in left_ts.darrays]]
In [103]: left_tss.shape
Out[103]: (1200, 32492)
tss = np.squeeze(nii.get_data())
In [105]: tss.shape
Out[105]: (91282, 1200)
In [106]: np.allclose(tss[0],left_tss[:,0])
Out[106]: False
In [109]: np.allclose(tss[:1200,0],left_tss[:,0])
Out[109]: True
In [110]: np.allclose(tss[1200:2400,0],left_tss[:,1])
Out[110]: True
#now I do compare the all the data from left cortex (hard reshaping
the data) excluding those which are not stored in the cifti (how are
there magically restored when doing cifti-separate is a mystery??)
In [117]: np.allclose(tss.T.reshape(tss.shape)[:29696],left_tss[:,idxs].T)
Out[117]: True
so forcing the slicing in the other direction gave me the right data,
so it seems there is a problem either:
- in the way the cifti library is using the nifti2 format to store the
array dimensions
- or there is a mismatch between row/column-major order storage of the data.
Do you have any insight on this problem?
Sorry again but another question is whether the cifti format/nifti2
implementation used in HCP opened to public?
If not, do you plan to release it to the community or it is intended
to remain closed source?
Many thanks for your help.
Cheers.
Basile Pinsard
On Thu, Dec 5, 2013 at 4:35 PM, Glasser, Matthew
<[email protected]> wrote:
> So there is a bug in the .dlabel.nii files that left and right ROIs are not
> kept separate (because FreeSurfer gives them the same name in both
> hemispheres). I have just fixed this for the next data release, but it will
> have to be worked around here by processing left and right hemispheres
> separately. My suggestion for extracting timecourses from the ROIs is to
> use this approach:
>
> wb_command -cifti-create-label
> ${StudyFolder}/${Subject}/MNINonLinear/fsaverage_LR32k/${Subject}.L.aparc.a2009s.32k_fs_LR.dlabel.nii
> -left-label
> ${StudyFolder}/${Subject}/MNINonLinear/fsaverage_LR32k/${Subject}.L.aparc.a2009s.32k_fs_LR.label.gii
> -roi-left
> ${StudyFolder}/${Subject}/MNINonLinear/fsaverage_LR32k/${Subject}.L.atlasroi.32k_fs_LR.shape.gii
> wb_command -cifti-create-label
> ${StudyFolder}/${Subject}/MNINonLinear/fsaverage_LR32k/${Subject}.R.aparc.a2009s.32k_fs_LR.dlabel.nii
> -right-label
> ${StudyFolder}/${Subject}/MNINonLinear/fsaverage_LR32k/${Subject}.R.aparc.a2009s.32k_fs_LR.label.gii
> -roi-right
> ${StudyFolder}/${Subject}/MNINonLinear/fsaverage_LR32k/${Subject}.R.atlasroi.32k_fs_LR.shape.gii
>
> For left and right dlabel files, parcellate the dense timeseries to get
> average timecourses in each ROI:
> wb_command -cifti-parcellate <dtseries> <dlabel> COLUMN <ptseries>
>
> To get the average ROI timeseries out of the parcellated timeseries file:
> wb_command -nifti-information <ptseries> -print–matrix
>
> The alternative, if you were more comfortable with GIFTI files would be to
> split the CIFTI file into components:
>
> wb_command -cifti-separate <dtseries> -metric CORTEX_LEFT <func.gii> -metric
> CORTEX_RIGHT <func.gii>, would would produce GIFTI timeseries at which
> point, you could directly use the GIFTI label files.
>
> Peace,
>
> Matt.
>
> From: basile pinsard <[email protected]>
> Date: Thursday, December 5, 2013 3:03 AM
> To: Matt Glasser <[email protected]>
> Cc: "[email protected]" <[email protected]>
> Subject: Re: [HCP-Users] Q3 rfMRI
>
> Hi Matt,
>
> yes I mean homotopic connections.
> The individual dense connectivity map in the image you sent is not available
> from the HCP DB, am I wrong?
> Was this computed from dense timeseries ICA-FIX as provided online using the
> wb_command -cifti-correlation? Does it add global signal regression to get
> the negative parts?
>
> Maybe the way I extract the signal from cortical ROIS is wrong, but I
> heavily rely on the cifti header to determine which rows belong to each ROIs
> so I used the "brain models" to get the IndexOffset and IndexCount of the
> different brain models, then used the
> ??????.?.aparc.a2009s.32k_fs_LR.label.gii to get the ROIs in the surface
> brain models using only the vertices listed in the CIFTI_STRUCTURE_CORTEX_*
> brain model NodeIndices xml element (which certainly exclude the filled
> holes in surface and fs medial wall ROI).
> for instance in the 101915 first scan LR , there are the following brain
> models:
>
> model type structure
> index count
> CIFTI_MODEL_TYPE_SURFACE CIFTI_STRUCTURE_CORTEX_LEFT 0 29696
> CIFTI_MODEL_TYPE_SURFACE CIFTI_STRUCTURE_CORTEX_RIGHT 29696
> 29716
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_ACCUMBENS_LEFT 59412
> 135
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_ACCUMBENS_RIGHT 59547
> 140
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_AMYGDALA_LEFT 59687
> 315
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_AMYGDALA_RIGHT 60002
> 332
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_BRAIN_STEM 60334
> 3472
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_CAUDATE_LEFT 63806
> 728
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_CAUDATE_RIGHT 64534
> 755
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_CEREBELLUM_LEFT 65289
> 8709
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_CEREBELLUM_RIGHT 73998
> 9144
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_DIENCEPHALON_VENTRAL_L
> 83142 706
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_DIENCEPHALON_VENTRAL_R
> 83848 712
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_HIPPOCAMPUS_LEFT 84560
> 764
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_HIPPOCAMPUS_RIGHT 85324
> 795
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_PALLIDUM_LEFT 86119
> 297
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_PALLIDUM_RIGHT 86416
> 260
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_PUTAMEN_LEFT 86676
> 1060
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_PUTAMEN_RIGHT 87736
> 1010
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_THALAMUS_LEFT 88746
> 1288
> CIFTI_MODEL_TYPE_VOXELS CIFTI_STRUCTURE_THALAMUS_RIGHT 90034
> 1248
>
> but the problem might be in the mapping of ROIs from the aparc.a2009s+aseg
> gifti files.
>
> Sorry if it is not clear but I had to reverse engineer the stuff so I may
> have done it wrong and it might not be possible to do it without some prior
> hidden in the workbench.
>
> Thanks for you insight on this.
>
> Basile
>
>
> On Thu, Dec 5, 2013 at 12:52 AM, Glasser, Matthew <[email protected]>
> wrote:
>>
>> I'm not sure what you mean by homologous connectivity. Does that mean
>> connectivity between corresponding parts of the hemispheres (like is shown
>> in the attached png from subject 101915—seed vertex is white, corresponding
>> vertex is green)?
>>
>> These data are acquired with a very fast TR (0.720ms), which means that
>> the movement parameters may pick up physiological motion effects in addition
>> to actual subject movements. In slower TR dat sets these are aliased into
>> the timeseries so you cannot see them.
>>
>> Peace,
>>
>> Matt.
>>
>> From: basile pinsard <[email protected]>
>> Date: Wednesday, December 4, 2013 4:24 PM
>> To: "[email protected]" <[email protected]>
>> Subject: [HCP-Users] Q3 rfMRI
>>
>> Hi all,
>>
>> trying to process individual rfMRI from HCP, we found unexpected results
>> of connectivity:
>> - almost no homologous connectivity (see attachment) either by taking mean
>> timeseries in aparc.a2009s+aseg rois(as store in fsaverage32k gifti files
>> for surfaces) or by computing the 30gigs voxelwise dense matrix zscore and
>> averaging in the same rois.
>> We did it with and without regressing out mean gray timecourse from the
>> same dense timeseries because it seems to influence largely the
>> connectivity.
>> We also did in on the base preprocessed, and ica-fixed data to check,
>> results are quite similar.
>> Will the pipeline used for preprocessing be released publicly at some
>> point?
>>
>> - the motion estimated is quite shaky and unrealistic (for example the
>> scans of subject 101915 enclosed), do you know which method was used for
>> realignment of the data? If this estimated motion is used both for
>> interpolation of data and regression of nuisance signal, this has certainly
>> huge effect on the connectivity.
>>
>> Was the same data/preprocessing used for estimating the 40subject dense
>> correlation matrix ? This has expected connectivity structure.
>>
>> Thanks for your help.
>> Best.
>>
>> Basile Pinsard.
>>
>> _______________________________________________
>> HCP-Users mailing list
>> [email protected]
>> http://lists.humanconnectome.org/mailman/listinfo/hcp-users
>>
>>
>>
>> ________________________________
>>
>> The materials in this message are private and may contain Protected
>> Healthcare Information or other information of a sensitive nature. If you
>> are not the intended recipient, be advised that any unauthorized use,
>> disclosure, copying or the taking of any action in reliance on the contents
>> of this information is strictly prohibited. If you have received this email
>> in error, please immediately notify the sender via telephone or return mail.
>>
>> _______________________________________________
>> HCP-Users mailing list
>> [email protected]
>> http://lists.humanconnectome.org/mailman/listinfo/hcp-users
>
>
>
>
> ________________________________
>
> The materials in this message are private and may contain Protected
> Healthcare Information or other information of a sensitive nature. If you
> are not the intended recipient, be advised that any unauthorized use,
> disclosure, copying or the taking of any action in reliance on the contents
> of this information is strictly prohibited. If you have received this email
> in error, please immediately notify the sender via telephone or return mail.
>
> _______________________________________________
> HCP-Users mailing list
> [email protected]
> http://lists.humanconnectome.org/mailman/listinfo/hcp-users
_______________________________________________
HCP-Users mailing list
[email protected]
http://lists.humanconnectome.org/mailman/listinfo/hcp-users
