Neighbor correlation logic

[cookpa]

Hi @frankyeh,

I was looking at this code with @mattcieslak

https://github.com/frankyeh/DSI-Studio/blob/8b5b13c2e5e05b83d5ebe8a8234cc715e66242a2/libs/dsi/image_model.cpp#L438C1-L482C2

with a view to implementing the DWI neighbor correlation in Python for inclusion in MRIQC. I am having a little trouble understanding the pairing logic

    std::vector<std::pair<size_t,size_t> > corr_pairs;
    for(size_t i = 0;i < src_bvalues.size();++i)
    {
        if(src_bvalues[i] == 0.0f)
            continue;
        float min_dis = std::numeric_limits<float>::max();
        size_t min_j = 0;
        for(size_t j = i+1;j < src_bvalues.size();++j)
        {
            tipl::vector<3> v1(src_bvectors[i]),v2(src_bvectors[j]);
            v1 *= std::sqrt(src_bvalues[i]);
            v2 *= std::sqrt(src_bvalues[j]);
            float dis = std::min<float>(float((v1-v2).length()),
                                        float((v1+v2).length()));
            if(dis < min_dis)
            {
                min_dis = dis;
                min_j = j;
            }
        }
        corr_pairs.push_back(std::make_pair(i,min_j));

This seems to make two pairs for each pair (i,j), considering only candidate pairs j > i. So for i=0, vector src_bvectors[i] will be compared to vectors j in the range (1,src_bvalues.size()). So we may find a pair (0,a) where a > 0. But when a=i, we will compare vectors j in the range (a+1, src_bvalues.size()). So it will be paired with some other vector.

I feel like I'm missing something here, because it seems that as i increases, we are picking pairs from an increasingly small pool of candidates. It also seems to that I would get different pairs if I shuffled the acquisition order of the bvals and bvecs.

Apologies if I've got this totally wrong!

[frankyeh]

You do make a good point! This implementation may give a false underestimate of the neighboring DWI correlation. My mistake. I am correcting it now.

On Wed, Mar 27, 2024 at 10:20 AM Philip Cook @.***> wrote:

Hi @frankyeh https://github.com/frankyeh,

I was looking at this code with @mattcieslak https://github.com/mattcieslak

https://github.com/frankyeh/DSI-Studio/blob/8b5b13c2e5e05b83d5ebe8a8234cc715e66242a2/libs/dsi/image_model.cpp#L438C1-L482C2

with a view to implementing the DWI neighbor correlation in Python for inclusion in MRIQC. I am having a little trouble understanding the pairing logic

std::vector<std::pair<size_t,size_t> > corr_pairs;
for(size_t i = 0;i < src_bvalues.size();++i)
{
    if(src_bvalues[i] == 0.0f)
        continue;
    float min_dis = std::numeric_limits<float>::max();
    size_t min_j = 0;
    for(size_t j = i+1;j < src_bvalues.size();++j)
    {
        tipl::vector<3> v1(src_bvectors[i]),v2(src_bvectors[j]);
        v1 *= std::sqrt(src_bvalues[i]);
        v2 *= std::sqrt(src_bvalues[j]);
        float dis = std::min<float>(float((v1-v2).length()),
                                    float((v1+v2).length()));
        if(dis < min_dis)
        {
            min_dis = dis;
            min_j = j;
        }
    }
    corr_pairs.push_back(std::make_pair(i,min_j));

This seems to make two pairs for each pair (i,j), considering only candidate pairs j > i. So for i=0, vector src_bvectors[i] will be compared to vectors j in the range (1,src_bvalues.size()). So we may find a pair (0,a) where a > 0. But when a=i, we will compare vectors j in the range (a+1, src_bvalues.size()). So it will be paired with some other vector.

I feel like I'm missing something here, because it seems that as i increases, we are picking pairs from an increasingly small pool of candidates. It also seems to that I would get different pairs if I shuffled the acquisition order of the bvals and bvecs.

Apologies if I've got this totally wrong!

— Reply to this email directly, view it on GitHub https://github.com/frankyeh/DSI-Studio/issues/84, or unsubscribe https://github.com/notifications/unsubscribe-auth/AACDI4LI7ASPQQOFR4IY4ULY2LIR7AVCNFSM6AAAAABFK7GY7OVHI2DSMVQWIX3LMV43ASLTON2WKOZSGIYTAOJTG43TQOA . You are receiving this because you were mentioned.Message ID: @.***>

[frankyeh]

Here's the commit to correct this issue: https://github.com/frankyeh/DSI-Studio/commit/3596f8ff81859dabc1d0e9165fd02014f6d338ef

I also remove the mask requirement and added other code to handle some exceptional conditions.

Thank you again for reporting this issue.

[mattcieslak]

Hi Frank, thanks for the quick fix! I was wondering if it might be possible to include a masked and unmasked NDC? Without the mask it's likely that the size of the field of view will affect NDC

[frankyeh]

Good point. I am adding it now.

On Wed, Mar 27, 2024 at 10:55 AM Matt Cieslak @.***> wrote:

Hi Frank, thanks for the quick fix! I was wondering if it might be possible to include a masked and unmasked NDC? Without the mask it's likely that the size of the field of view will affect NDC

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were mentioned.Message ID: @.***>

[cookpa]

Thanks so much!

[frankyeh]

https://github.com/frankyeh/DSI-Studio/commit/9aeb41dfaeeaf0dc762ce4ffdcac788c5629bdf0 now both masked and unmasked version are reported in the QC

[frankyeh]

Here I am also testing new QC metrics for DWI contrast: https://github.com/frankyeh/DSI-Studio/commit/6c53317d0b21b525c06f440b48818efcab87eebc and https://github.com/frankyeh/DSI-Studio/commit/d2a06b5e58ef8c2b736e019ac999e7edf77fa152

It is defined as the ratio between neighboring DWI correlation and orthogonal DWI correlation (NDC/ODC). The "orthogonal DWI" is the DWI nearest to the orthogonal ring of an evaluating DWI in the q-space. It works for single-shell, multi-shell, and grid data.

Good DWI contrast will have a ratio value larger than 1.2

For DWI acquired with a very low b-value, the NDC will be very good, but the DWI contrast will be close to 1.

[mattcieslak]

this is good to know, I am in QC mode right now and will be adding them

[frankyeh]

It is open to better suggestions. I am also testing using different data. Frank

On Wed, Mar 27, 2024 at 2:54 PM Matt Cieslak @.***> wrote:

this is good to know, I am in QC mode right now and will be adding them

— Reply to this email directly, view it on GitHub https://github.com/frankyeh/DSI-Studio/issues/84#issuecomment-2023716465, or unsubscribe https://github.com/notifications/unsubscribe-auth/AACDI4IBVM6IKYUEKS3MEQLY2MIWLAVCNFSM6AAAAABFK7GY7OVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAMRTG4YTMNBWGU . You are receiving this because you were mentioned.Message ID: @.***>