[pymvpa] [pynifti] NiftiDatatset bug?

[Michael Hanke]

Hi Scott,

[taking the discussion to the new list. Please drop CC on reply ]


On Thu, Jan 08, 2009 at 04:54:34PM -0500, Scott Gorlin wrote:
<snip>
> I agree with you except that the dimension index of an ND array should 
> be consistent no matter how many dimensions they have.  In other words, 
> shape[0] should have the same meaning for the same type of data 
> regardless of how many samples (or u,v,w dims) there are.  It seems self 
> evident this should be the case from the pymvpa point of view at least, 
> since shape[0] is always treated as the number of samples (which is 1 
> for a 3d volume).  With proper mapping handled in the NiftiDataset this 
> is less of an issue, but it should also be easy to manipulate this 
> manually when not using pymvpa - as I'm sure you intend pynifti for 
> general usage.

Right. From the PyMVPA point it is that easy: shape[0] is the #samples.
That works well in non-4d cases, i.e. each sample is a 4d volume on its
own, hence your input file is 5D.

> Either that, or the data should be handled privately (in NiftiImage) in 
> an intelligent way that can translate the correct dimensions into the 
> actual dimensions in the data array.  Perhaps a function like 
> getFullData(t=True, u=False, v=False, w=False) which extracts a 4d 
> (T,X,Y,Z) even if there's only one sample?  Otherwise you get bugs like 
> this.  For instance, I don't know ahead of time if i'll be loading a 3D, 
> 4D, or 4D with 1 timepoint volume (which apparently is treated 
> differently, as my cludgy hack works that way) and every time i need to 
> use the array I have to check for where the dimensions are.

Why? The dimensions are always the reverse of x,y,z,t,u,v,w, no?

I wonder what your '4d with 1' image header looks like. See this:

  In [3]: N.arange(24).reshape(1,2,3,4).shape
  Out[3]: (1, 2, 3, 4)

  In [4]: nifti.NiftiImage(N.arange(24).reshape(1,2,3,4))
  Out[4]: <nifti.niftiimage.NiftiImage object at 0xa529a6c>

  In [5]: nim=nifti.NiftiImage(N.arange(24).reshape(1,2,3,4))

  In [6]: nim.data.shape
  Out[6]: (1, 2, 3, 4)

  In [7]: nim.save('test.nii.gz')

  In [8]: nim2=nifti.NiftiImage('test.nii.gz')

  In [9]: nim2.data.shape
  Out[9]: (1, 2, 3, 4)

  In [10]: nim.header['dim']
  Out[10]: [4, 4, 3, 2, 1, 1, 1, 1]


This is exactly what I would expect -- and if I got you right this is
what you want too, right?

> To further question this - what would happen to the Nifti data if the 
> UVW dims are not 1?  Presumably this is rare... but would they keep 
> adding on to the beginning of the array?
Yes.

>  So the number of samples now 
> falls in the 2nd, 3rd, or 5th position? 
No, it is always the first.

> Now perhaps I extract a subset 
> of U,V,W 4-D volumes and squeeze away the preceding singleton dims - now 
> there's no guarantee to the number of dims I'm left with.  Or are U,V,W 
> mapped after Z? And what about single-slice, or other 1d and 2d images 
> encoded in an nii?  It seems like wanting to squeeze out dims stands at 
> ends with dim compatibility in this case.
Hmm, I did not get this.

> While I'm proposing crazy ideas, it would also be simpler if the index 
> just happened to be reversed - (Z,Y,X,T,W,V,U).  That way the proper 
> volume is always consistent and extra dimensions can be safely dropped.  
> I think in Matlab, nii's are typically loaded like this (I usually load 
> using Freesurfer tools, not sure about SPM).  It would also be easy to 
> concatenate samples or uvw coords without having to permute the array first.

I don't see the advantage -- apart from the fact that this would change
the order of information in the file and increase IO operation time.


> I understand these ideas are probably over the top complicated in the 
> face of backwards compatibility - so perhaps they can be my 'modest 
> proposal'.  But I hope they illustrate my point - IMHO it's crucial to 
> easily extract and concatenate 3D volumes from a NiftiImage without 
> having to figure out what dimension 'samples' is stored in  every time 
> it's necessary.

We should not confuse the needs of PyNIfTI with the ones of PyMVPA.

But in general it is quite simple: The last three dimensions/axes are
always the spatial ones, the one in front of them is the time and the
rest is something else.

> In other words - if you are heading to a new release with pynifti, yes, 
> i think this issue is worth addressing there :)

Alright, but I still have no clear concept of what the 'issue' is in the
first place -- let discuss some more ;-)

> sorry if this seems like a rant - just my two drops of oil...

Don't worry -- I am German. Start ranting ;-)


Michael


-- 
GPG key:  1024D/3144BE0F Michael Hanke
http://apsy.gse.uni-magdeburg.de/hanke
ICQ: 48230050