[Insight-developers] Can someone explain DiffusionTensor3D::GetRelativeAnisotropy() to me?

[Torsten Rohlfing]

Simon Warfield wrote:

>> The code (purged of numerical safeguards etc) that I don't get is this:
>>
>>  const RealValueType anisotropy = 3.0 * isp - trace * trace;
>>  const RealValueType relativeAnisotropySquared =
>>        static_cast< RealValueType >( anisotropy / ( sqrt( 3.0 ) * 
>> trace ) );
>>  const RealValueType relativeAnisotropy =
>>        static_cast< RealValueType >( sqrt( relativeAnisotropySquared 
>> ) );
>>
>> where "isp" is the inner scalar product and "trace" is the tensor 
>> trace, both as returned by the respective class member functions.
>>
>> My problem is with the definition of "anisotropy" above. I assume 
>> that this is supposed to be the inner product of the anisotropic part 
>> of the tensor, that is, the inner product of  (tensor - (trace / 3) * 
>> identity). This is a double sum over the squares of differences 
>> between elements of the original tensor and the corresponding 
>> elements of (trace/3) times the identity tensor.
>
>
> Check out Equation 27 here:
> http://splweb.bwh.harvard.edu:8000/pages/papers/westin/media2002/Westin_Diffusion_MedAI_2002.pdf 
>
> and the definition of the norm at the top of the column.
>
Simon:

Thank you, but that's exactly my problem - how do you get from Eq. 27 in 
the Westin paper to the implementation in ITK? The Westin equation is 
basically the same (except for the sqrt(2) factor) to Eq. 14 in the 1996 
JMR(B) paper by Basser and Pierpaoli, and their Eq. 10 is the tensor 
norm as described in the Westin paper. But how do you get from the 
numerator in the second line of Westin's Eq. 27 to "3*isp-trace*trace" 
(or isp+1/3*trace*trace, which is basically the same)?

I assume that the sum of the squared elements is broken up, that is, the 
differences that make up each element are separated as (a-b)^2 -> 
(a^2-2ab+b^2), and so sum_i sum_j (a-b)^2 becomes sum_ij a^2 - 2sum_ij 
ab + sum_ij b^2.

Then with a_ij==D and b_ij==1/3 trace(D)I, that would give

sum_ij D_ij^2 - 2/3 trace(D) sum_ij I_ij + (1/3 trace(D))^2 sum_ij I_ij

and with the first term the norm of D and I the identity tensor we would 
have

|D| - 2 trace(D) + 1/3 trace(D)^2

So where is the middle part going in ITK?

Sorry if I am suffering from a mental block here.

Again, thanks for your help!
  Torsten

-- 
Torsten Rohlfing, PhD          SRI International, Neuroscience Program
 Research Scientist             333 Ravenswood Ave, Menlo Park, CA 94025
  Phone: ++1 (650) 859-3379      Fax: ++1 (650) 859-2743
   torsten at synapse.sri.com        http://www.stanford.edu/~rohlfing/

     "Though this be madness, yet there is a method in't"