Using the PKmodeling extension once the dynamic contrast data has been fit and the individual voxel data visually assessed a portion of the voxels are not fit. This doesn’t appear to be related to noise or data quality. To show this I attached images from 2 voxels with similar concentration profiles with only 1 having been successfully fit by the extension. Is this a expected feature of this kind of analysis? Or should I expect all voxels to be successfully fit?
@bmb777 thank you for investigating this and sharing your feedback!
It is hard to tell what went wrong in this case. As a general answer, no - you cannot expect all voxels to be fit. But I agree with you that from the look of it, this specific voxel you identified should probably work.
To investigate what happened, as the first step, you should toggle generation of the diagnostic image (see Advanced options). The diagnostic image will have an error code explaining why the fit failed, as described in the documentation: https://www.slicer.org/wiki/Documentation/Nightly/Modules/PkModeling (better formatted here for readability!). Unfortunately, there are many reasons why it could have failed …
Output map with the optimizer diagnostics. The code is encoded in 2 hex numbers. Lower 4 bits encode the optimizer errors are as follows:
- 0: OIOIOI – failure in leastsquares function
- 1: OIOIOI – lmdif dodgy input
- 2: converged to ftol
- 3: converged to xtol
- 4: converged nicely
- 5: converged via gtol
- 6: too many iterations
- 7: ftol is too small. no further reduction in the sum of squares is possible.
- 8: xtol is too small. no further improvement in the approximate solution x is possible.
- 9: gtol is too small. Fx is orthogonal to the columns of the jacobian to machine precision.
- 10: OIOIOI: unknown info code from lmder.
- 11: optimizer failed, but diagnostics string was not recognized.
Upper 4 bits encode other non-optimizer errors or notifications:
- 16 (0x10): Ktrans was clamped to [0…5].
- 32 (0x20): Ve was clamped to [0…1].
- 48 (0x30): BAT detection failed.
- 64 (0x40): BAT at the voxel was less than AIF BAT.
You can also see some illustrations in these PRs (of course, it will not look as colorful and pretty for the default color map assigned to the fit diagnostic map):
Let us know what you find!
Thank you for providing the diagnostic codes. However, I may be misunderstanding the interpretation of the output diagnostics, as the value at the pixels that are not fit but appear as if the should has a diagnostic code of 50, which isn’t provided as a error in the numbers above. Is there something I am misunderstanding?
If you have 50, this means 48 (BAT failure) + 2 (converged to ftol).
You can try using
UseConstantBAT mode in the advanced options, and set it to 7 (the frame number corresponding to the initial sharp uptake of contrast, based on the curves you share earlier).
Although, looking at the curve, I do not see why BAT detection would fail there. Would be good to investigate this. If you can share a sample dataset, please submit an issue and attach the dataset here: https://github.com/QIICR/PkModeling/issues/new.
For that particular fit I had selected the constant batch option and set it to 7.
I am currently doing a fairly intensive DCE data analysis for a conference next week. I’ll submit a sample dataset with the issue mentioned once that is completed.
Thanks for the information.
Did I understand correctly - when you set BAT to 7 fit works for the affected voxels?
Whenever you have time to respond. Good luck with the submission!
No, the fits I attached in post one above were obtained using a constant BAT. One of the reasons I was confused why the software returned a diagnostic code of failure to calculate a BAT.
As a secondary question, is the AUC parameter calculated based of the PK fits or the raw converted concentration versus time data? For example in the above scenario will the C v. T profile that was not fit for a voxel return a AUC value of 0 due to the failed fit?
Sorry for the delayed reply, I know you are working on a deadline, but I was away Thu-Fri.
Yes, indeed AUC is calculated for the fitted signal and also relies on BAT for aligning the AIF and signal curves. Just since I had to find this to answer your question, AUC is calculated in this line: https://github.com/QIICR/PkModeling/blob/master/CLI/itkConcentrationToQuantitativeImageFilter.hxx#L394. and
shiftedVectorVoxel is the fitted vector shifted to compensate for the BAT delay: https://github.com/QIICR/PkModeling/blob/master/CLI/itkConcentrationToQuantitativeImageFilter.hxx#L342.
Hope this helps.