[Insight-users] Do RecursiveHessian or RecursiveGaussian requireisotropic voxels?

Luis Ibanez luis.ibanez at kitware.com
Tue Nov 22 07:18:43 EST 2005


Hi Jerome,

Thanks for the additional details about your processing.

Your approach of computing the Hessian only in the points
that are needed is certainly a good way of improving the
performance of the feature extraction.

It will be great if you write the class that computes the
Hessian locally, and post it as an article to the Insight
Journal.

  http://www.insightsoftwareconsortium.org/InsightJournal/


In this way we will be able to include it in the following
release of ITK (2.6 for February 2006).

You may want to follow the model of the GaussianFunction:

http://www.itk.org/Insight/Doxygen/html/classitk_1_1GaussianDerivativeImageFunction.html


Note that papers to the Insight Journal are more on the style
of a technical report. They should describe your motivation
for introducing this class (which is quite straight forward
in your case), and should provide tests for it, and examples
of its usage. The ultimate goal is to enforce REPRODUCIBILITY.

(a fundamental feature of scientific research that is totally
  neglected by the "Vanity Journals" in our field, because
  they are designed just for allowing the authors to claim
  publications in their CVs and annual reports).

REPRODUCIBILITY is what separates Science from Superstition.

---

Note also that reports of negative results are also welcome
and encouraged in the Insight Journal. This is because we
are interested in actually solving research problems, and not
in the hedonistic and pedantic exercise of reading more authors
claiming that their ideas are better than anything published
so far.

Publication of Negative results, in a format that permits its
reproducibility is a fundamental aspect of scientific research.


When you report that you attempted to use the DuctExtractor
on your image and

                   "it didn't work for you".

    That claim is useless in its current form in an email.

     It is useless for you and it is useless for us.

     It doesn't help you to solve your problem

         It doesn't help us to help you


What we need is that you report a full and detailed description
of how you used that program. We need the images that you used
as input, and we need the exact numerical values of the parameters
that you put in the GUI.


     With that information we *can* reproduce what you did,

                    and then

     we can figure out why it didn't work,

                    and then

     we can let you know how to make it work,

     by sending back to you a similar detailed description of
     the numerical parameters that may be appropriate for you
     image and for your application.


   That, will be useful for you because you solve your problem.
   That, will be useful for us because we extend the toolkit.
   That, will be useful to the community because they will learn.


      Our generation has been so long in the Darkness of the
      "Publish and Perish" that even the most obvious and
      basic behaviors of research have been rusted and forgotten.

      Too many bright minds have been lost so far,...

      Hopefully the next generation can still be spared
      from the purposeless lives of "Vanity Publishing".



Please post your report to the Insight Journal.

Instructions on how to submit can be fount at

http://www.insightsoftwareconsortium.org/wiki/index.php/Insight-Journal


Let us know if you find any difficulty during the
submission of your report.



Best Regards,



    Luis



----------------------
SCHMID, Jerome wrote:
> Luis Ibanez wrote:
> 
>  >
>  > Hi Jerome,
>  >
>  > You are right !
>  >
>  >    I'm going to tell you that acquiring datasets
>  >    with 5mm interslice spacing is just *bad imaging*.
>  >
>  >
>  > We hope that one day MRI radiologist will start looking
>  > at the orthogonal views of the dataset they acquire and
>  > perceive how much they are violating the Shannon sampling
>  > theorem.
>  >
>  > We have heard the standard replies:
>  >
>  >  *  "It would take too long..."
>  >  *  "The information in the slice is enough..."
>  >
>  > The sad reality is that algorithm developers are not
>  > working close enough to those who design the image
>  > acquisition protocols, and radiologist still look
>  > only at slices, not at 3D view of the the datasets
>  > they acquire. So they rarely realize, or care about
>  > the fact that those are not really 3D datasets, they
>  > are just collections of isolated and uncorrelated
>  > slices.
>  >
>  >
>  > ....
> 
> Yes I read the beautiful and quite explicit long comment in the software 
> guide /  src  code and  you are perfectly  right ^_^. I didn't control 
> the acquisition of the data and I know that 5 mm is just a scandal.
> 
>  >
>  >
>  >
>  > About your processing,
>  >
>  > Thanks for posting the Phantom dataset,
>  > it make things much clearer.
>  >
>  >
>  > The cropped region that you posted is unfortunately too
>  > small for performing the processing on it.  The reason
>  > is that the vessels that you are looking at are in the
>  > range of 20 pixels in diameter. The sigma that will capture
>  > those features is such that for your cropped image the
>  > border effects will take over the image before we get
>  > the vessels to show up.
>  >
>  > (unless, of course, I'm looking at the wrong vessels...)
>  >
>  >
> Sorry for the cropping it was for a sake of light data. Moreover I don't 
> know if I can post the whole phantom data, maybe I need some agreement 
> with the company.. However I was more interested  in the tumour 
> eigenvalues ( the ball not the vessels ) whose size is around 7 mm. The 
> purpose was not in the vessels. However you pointed something 
> interesting, i.e. the influence of background in the filter computation, 
> thanks.
> 
>  >
>  >
>  > The application that you want to use is available in
>  >
>  >
>  >        InsightApplications/
>  >                   DuctExtraction
>  >
>  >
>  > It uses the Vesselness filter, (based on the Hessian),
>  > and adds some preprocessing by a gradient anisotropic
>  > filter, and post-processing with a region growing algorithm.
>  >
>  > If you load your images in, you will have an easy test-bed
>  > for your vessel detection program.
>  >
>  >
> 
> I already saw it thanks. Actually I played it a while and cannot make it 
> work successfully. The data I used is not the ugly phantom but a 3D 
> Ultrasound data. Unfortunately you may know that US data is also the 
> devil ^_^. I found that the itkHessian3DToVesselnessMeasureImageFilter 
> was more efficient with this kind of data, so I first created a new 
> DuctExtractor app in order to use the vesselness filter ( BTW if you are 
> interested on it I can post it. ). Seeing that it works not so badly ( 
> but require some improvement ) I wonder whether or not it works 
> correctly with the CT data of the phantom. But I encountered all the pbs 
> I told you before....So I hope that the reason is the bad data sampling...!
> 
> Mmmm....okay I will assume that so badly sampled data cannot be trusted 
> and continue to play with hessian and eigen values.
> 
>  I have written a filter that makes a region growing based on a 
> ImageFunction that returns true or false depending a vesselness w.r.t 
> thresholds. This implies that the imageFunction computes the vesselness 
> on isolated voxels instead of the whole image. Internally I used the 
> RecursiveHessianImageFilter but I am not sure whether it is efficient or 
> not.
> 
> The biggest problem I have encountered is in case of multiscale 
> analysis: by changing the sigma in the RecursiveHessianImageFilter a 
> "reset" occurs ( all the coefficients or things similar of the 
> RecursiveGaussian filter have to be recomputed I guess ) and the 
> performance is extremelly bad. If I stay with the same sigma and call it 
> once on one pixel, after all the calls on other pixels are fast because 
> some "cache" is used. But if I would like to perform a multiscale 
> analysis on one pixel (i.e. SetSigma(1.0), compute, SetSigma(2.0), 
> compute, etc.) the cache cannot be exploited.
> 
> So In order to resume, do you know a way to have an efficient 
> ImageFunction that compute the Hessian on a voxel? Should I write 
> something totally new?
> 
> Thanks and sorry for such a long email...!
> 
> Jerome.
> 
>  > ...
>  >
>  >
>  > The reason why you get better results from the IsoPhantom
>  > (the image that was resampled from the 5mm datasets) is that
>  > the RecursiveGaussian image filter use IIR filters internally
>  > and it bases its computation on 4 consecutive samples. It doesn't
>  > mean, however that the information in the IsoPhantom image is better.
>  >
>  > The solution is to go back to your image acquisition friends,
>  > and show them how their datasets look in orthogonal slices.
>  >
>  >
>  > As Dr. Michael Vannier said once at SPIE:
>  >
>  >
>  >        "Good imaging,
>  >         beats good image processing"
>  >
>  >
>  >
>  > ....
>  >
>  >
>  > Please give it a try to the DuctExtraction application, and let us
>  > know if you find any problem with its usage.
>  >
>  > Note that the application was assuming that the vessels were bright
>  > over a dark background (since most of the time there is a contrast
>  > agent involved). If this is not the case for your images, you may
>  > want to revert the intensities in your image before feeding it into
>  > this program. (The IntensityWindowingImageFilter may be used for this
>  > purpose).
>  >
>  >
>  >
>  >    Best Regards,
>  >
>  >
>  >       Luis
>  >
>  >
>  >
>  > ----------------------
>  > SCHMID, Jerome wrote:
>  >
>  >>
>  >> Luis Ibanez wrote:
>  >>
>  >>  >
>  >>  >
>  >>  > Hi Jerome,
>  >>  >
>  >>  > No, the RecursiveGaussian filter does not requires the voxels
>  >>  > to be isotropic.
>  >>  >
>  >>  > The RecursiveGaussian filter takes the Pixel size into account.
>  >>  >
>  >>  > This means that the computations are performed correctly in
>  >>  > physical units.
>  >>  >
>  >>  > The values reported by the RecursiveGaussian and its associated
>  >>  > filters must be interpreted in terms of derivatives with respect
>  >>  > to millimeters.
>  >>  >
>  >>  > Note that if what you are doing is multi-scale analysis (and it
>  >>  > sounds like that, since you are searching for the sigma that matches
>  >>  > a particular blob...) then you should use the Normalize along scales
>  >>  > option, which will make possible for you to compare the results from
>  >>  > one scale (sigma) with the results of a different scale (sigma).
>  >>
>  >> Hi Luis,
>  >>
>  >> Thanks for the reply. Well I have seen the option of the 
> normalization across scale, but my problem is not related to a 
> mutliscale analysis. Anyway I set the flag to true already so it is not 
> correlated I think. Moreover you confirm to me that the isotropy should 
> not be a pb as the sigma is properly used as a physical unit instead of 
> a voxel one.
>  >>
>  >>  >
>  >>  >
>  >>  > ---
>  >>  >
>  >>  > Just for the record...
>  >>  >
>  >>  > What is the size of the Blobs in your experiments in pixels ?
>  >>  > and what is the size of those pixels in millimeters ?
>  >>  >
>  >>
>  >> I am working on a phantom which has vessels and tumours in a fake 
> liver. My objective is the vessels segmentation and by consequence I 
> would like to detect them and *reject* the tumours. As I told previously 
> by analysing the eigenvalues of the hessian with an adequate scale, one 
> can have a "vesselness" value, as showed in one of the recent filter is 
> the CVS. (BTW which as a little error that I pointed in the bug tracker ).
>  >>
>  >> Unfortunately, my data has an ugly slice thickness: 5 mm....!  I 
> expect in the future more accurate data but in the meantime I do some 
> tests...
>  >> So the diameter of a tumour ( a nice ball ) is 15x15x15 in physical 
> units which correspond roughly to 21x21x3 in voxels.
>  >>
>  >> I guess that you may say me that due to the huge slice thickness 
> some approach may become invalid, but by resampling the volume in iso 
> data I get "better" results....
>  >>
>  >>  >
>  >>  > If you find that the behavior of the filter is incorrect, please
>  >>  > post a minimal piece of code that illustrates your point, and if
>  >>  > necessary post an input/output image on a web accessible place.
>  >>  >
>  >> I will do this:
>  >>
>  >> - Find to this address :  
> http://moscao.free.fr/PrivateDownloads/EigenToRGBSrc.zip,  a zip that 
> contains:
>  >> - an example test that is run like this, "exec InputImage 
> OutputImage [sigma]"
>  >> This example compute the hessian with the sigma value, and the 
> eigenvalues of the hessian matrix of each voxel. An RGB image is 
> eventually created whose each component is an eigenvalue. A useful 
> filter is so added to the archive called 
> "itkHessianEigenValuesToRGBImageFilter.h".
>  >>
>  >> To this address: http://moscao.free.fr/PrivateDownloads/Phantom.zip 
> , a zip archive that contains:
>  >>
>  >> A version of the phantom cropped but *anisotropic* that show a part 
> of the liver with a simple vessel and a tumour. NOTE: I am looking for 
> dark vessel with bright background.
>  >> Its size: 95x81x14 , spacing: 0.7, 0.7, 5.0
>  >> It is an hdr/img Analyze image. Load it in volview with this 
> orientation params: RAS + Posterior + Left. Then in slice 8 (in axial 
> direction) there is a tumour at (46,17,40) (physical units)
>  >>
>  >>
>  >> A version of the phantom cropped *isotropic* ( and BTW scaled 
> previously to avoid huge data) that show a part of the liver with a 
> simple vessel and a tumour.
>  >> Its size: 48x55x44 , spacing: 1.4, 1.4, 1.4
>  >> It is an hdr/img Analyze image. Load it in volview with this 
> orientation params: RAS + Posterior + Left. Then in slice 8 (in axial 
> direction) there is a tumour at (52,31,38) (physical units)
>  >>
>  >>
>  >> Now run the program with the *first* image and set a sigma between 3 
> and 6. Open the image in volview, find the corresponding tumour and look 
> at its eigenvalues: they should be roughly like this:
>  >> 0.005 , 40, 40 ---> NOT a blob...!?
>  >>
>  >> Now run the program with the *second* image and set a sigma between 
> 3 and 6. Open the image in volview, find the corresponding tumour and 
> look at its eigenvalues: they should be roughly like this:
>  >> 3 , 4, 4 ---> Looks more like a blob (i.e. eigenvalues close to each 
> other )
>  >>
>  >>
>  >> Maybe the pb is that one cannot trust the results with such 
> images/spacing but I wonder why by passing to isotropric voxels thre 
> result are more as expected...?
>  >>
>  >> Thanks!
>  >>
>  >> Jerome
>  >>
>  >>
>  >>  >
>  >>  > Best Regards,
>  >>  >
>  >>  >
>  >>  >    Luis
>  >>  >
>  >>  >
>  >>
>  >>
>  >>
>  >>
>  >>
>  >>  >
>  >>  > ----------------------
>  >>  > SCHMID, Jerome wrote:
>  >>  >
>  >>  >> Dear Itk users,
>  >>  >>
>  >>  >> I am playing with vesselness measure based on eigen values of 
> Hessian matrices. I have noticed a quite odd behaviour and I would like 
> to know if it is normal.
>  >>  >>
>  >>  >> Basically the eigenvalues change "a lot" whether or not my input 
> image has or not isotropic voxels. The "a lot" changing is related to 
> the tubes or blobs detection with the Hessian+eigenValues principle, 
> i.e. roughly lambda1 = lambda2 = lambda3 for a blob and lambda1 = 0 and 
> lambda1 << 0 and lambda3 << 0 for a tube.
>  >>  >>
>  >>  >> In my case with a blob and non isotropic image: lambda1 close to 
> 0 and lambda2/3 << 0. -> odd(?)
>  >>  >> In my case with a blob and isotropic image: lambda1 = lambda2 = 
> lambda3 (roughly) -> expected
>  >>  >>
>  >>  >> Knowing that the hessian is computed with the 
> RecursiveGaussianImageFilter that is based on the RecursiveGaussian 
> filter, I wonder whether the pb is here.
>  >>  >>
>  >>  >> A note: my sigma for the detection is close to the radius of the 
> blob, and it is defined in physical units. I expect that the recursive 
> gaussian filter uses sigma in physical units ( which looks to be the 
> case, looking at the distinction variance/sigma, isn'it?).
>  >>  >>
>  >>  >> Any suggestion? Thanks!
>  >>  >>
>  >>  >> Best Regards,
>  >>  >>
>  >>  >> Jerome Schmid
>  >>  >>
>  >>  >> -----------------------------------
>  >>  >> Jerome SCHMID
>  >>  >> Project Manager/ Engineer
>  >>  >> Augmented and Virtual Reality
>  >>  >> MIS Centre
>  >>  >> Prince of Wales Hospital
>  >>  >> Chinese University Of Hong-Kong
>  >>  >> -----------------------------------
>  >>  >>
>  >>  >>
>  >>  >> 
> ------------------------------------------------------------------------
>  >>  >>
>  >>  >> _______________________________________________
>  >>  >> Insight-users mailing list
>  >>  >> Insight-users at itk.org
>  >>  >> http://www.itk.org/mailman/listinfo/insight-users
>  >>  >
>  >>  >
>  >>  >
>  >>  >
>  >>
>  >>
>  >> ------------------------------------------------------------------------
>  >>
>  >> _______________________________________________
>  >> Insight-users mailing list
>  >> Insight-users at itk.org
>  >> http://www.itk.org/mailman/listinfo/insight-users
>  >
>  >
>  >
> 



More information about the Insight-users mailing list