An International Study of Factors Affecting Variability of Dosimetry Calculations, Part 4: Impact of Fitting Functions in Estimated Absorbed Doses

DISCUSSION

Previous publications on factors affecting variability in the dosimetry calculations of the SNMMI ¹⁷⁷Lu Dosimetry Challenge (5,13,18) showed combined variabilities from all elements of the dosimetry workflow. Here, we focused on analyzing variability in AD determination from time–activity data using submissions from tasks 4 and 5.

In the first step of this analysis, we identified the most used fitting methods for modeling different time–activity curve regions. Some models produced outliers and deserve some discussion. Specifically, the 4-parameter biexponential model (Eq. 4; Fig. 1, method 11) in some cases resulted in outliers in AD (Fig. 3). The fitting problem is ill-conditioned, because the number of parameters equals the number of time points. Therefore, it is not recommended to apply a 4-parameter model to fit 4 data points. Supplemental Figure 3A, where we replicated the plot of one of the participants, shows the model fit the data well, but it is not physiologically reasonable because of the high activity at time 0.

From Figure 4, we see that most organs had an almost instantaneous uptake; therefore, monoexponential fits were commonly used for organs (Fig. 1). However, some participants used monoexponential fits for lesions; for example, 3 participants used monoexponentials for lesion 2 of patient B, despite visual evidence (Fig. 4; Supplemental Fig. 3B) showing poor fit. Visual inspection is essential to ensure appropriate fitting functions (19).

The trapezoidal approach was frequently used to interpolate time–activity data, combined with several extrapolation methods for I3. When physical half-life was assumed (Fig. 1, method 1), AD results were higher than when the effective half-life was calculated (Figs. 3, 5, and 6). Although we recognize the benefits of simple models, they should be reserved for instances in which exponential models fail to fit the data adequately, and this was not the case for any region in the challenge dataset (Fig. 4). In addition, assuming physical decay after the last time point is a conservative approach, appropriate for radiation-protection purposes, but it may not be ideal for optimizing RPTs, because both overestimating and underestimating the AD can affect treatment optimization.

The second step of the analysis showed higher variability in ADs for lesions than for organs at risk, based on QCD differences between tasks and IQR of . This may be attributed to the more complex pharmacokinetics observed in lesions than in organs, suggesting that TIA and AD estimation for lesions could inherently have more variability because of uncertainties in time–activity data fitting and integration. Lesions typically showed an uptake phase followed by a washout phase (i.e., a biexponential behavior); a monoexponential function, which assumes instantaneous uptake, does not model this behavior. In addition, there was higher variability in the kidneys of patient B than in those of patient A. As before, the higher variability can be attributed to the slower uptake phase in patient B’s kidneys (Fig. 4). Similarly, among all lesions, ones with slower uptake, such as lesion 2 of patients A and B, had higher variability. For regions with slower uptake where the assumption of instantaneous uptake is invalid, a model that emphasizes robust curve fitting is essential, particularly during the slow washout phase, because it contributes most significantly to TIA (20). If data from multiple pharmacokinetic phases, both uptake and washout, are included, the model should account for these phases to prevent inaccuracies in later phases of the curve.

The observed higher variability in total kidney ADs than in individual kidney ADs (Fig. 2) may be explained by the likelihood that some participants added the ADs for each kidney separately, without accounting for differences in kidney mass. In addition, it can be caused by the estimation of the AD to 1 kidney without accounting for the mass of the combined kidneys in the dosimetry software. The accurate approach to combining the ADs of the separate kidneys is to perform the scaling of AD by the mass of each kidney to obtain a weighted average; an alternative approach is to ensure the mass of each kidney is used when calculated separately.

The paired analysis from Figure 3 shows that most outliers were observed with in-house fitting implementations. In addition, as noted in Supplemental Table 1, more outliers were reported for patient A than for patient B, potentially because of a failure to perform decay correction, which was required for patient A but not for patient B. The original spreadsheets for task 4 lacked activity at each time point, limiting our ability to verify whether decay correction was applied.

Step 3, conducted by the dosimetry expert, revealed a narrower AD range than that of the participant results. Outliers from participants’ analyses could not be reproduced using standard practices, even with the same fitting methods, suggesting these outliers stemmed from errors (e.g., software issues) or deviations from best practices, such as improper decay correction, rather than inherent variability in the fitting process.

This baseline analysis showed that simplified approaches, such as trapezoidal interpolation with physical or linear decay after the last time point (Fig. 1, methods 1 and 4), as well as models that do not apply additional constraints on the fit (Fig. 1, method 6), resulted in higher ADs in organs and lesions than with other models (Figs. 5 and 6). ADs calculated with RTStructs were higher than those calculated with binary masks, despite masks showing higher activity levels (Fig. 4), because RTStructs represented smaller volumes.

The variability in TIA calculation was highest in I3 (Fig. 7), emphasizing the importance of appropriate and consistent extrapolation after the last time point to decrease the variability in TIA and AD estimates. To mitigate errors from the contributions within this interval, the European Association of Nuclear Medicine guidelines recommend that the fractional contribution of the extrapolations from I3 to the TIA should be less than 20% (21). In a comparison of I3 results for patient A (last SPECT at 124 h after injection) and patient B (last SPECT at 193.3 h after injection), outliers in patient A’s data led to a TIA value double the median. This suggests that including later scans, as with patient B, could help mitigate such outliers. Variability in I1 was small because of its short duration and limited contribution to the total area under the curve, although assumptions during this interval still influenced overall TIA variability.

During the analysis, we identified several issues in participants’ submissions that affected activities, which serve as the input for pharmacokinetic modeling. These included improper decay correction and the inclusion of lesions in the assessment of healthy liver.

For patient A, images were decay-corrected to the injection time, whereas for patient B, they were corrected to the start of the SPECT scan. Despite instructions, some participants did not take this into account correctly, resulting in TIA overestimation ranging from 43% to 204%. This suggests that a quality-control check is needed to ensure that assumptions about decay correction are consistent with the data.

When estimating AD to healthy tissue for ¹⁷⁷Lu toxicity predictions, it is essential to clearly define VOIs. Specifically, the uptake of activity by lesions located within the tissue (e.g., lesions in the liver) should not be included in calculating the total organ TIA (13); including this activity could significantly increase the whole-organ TIA and result in an overestimation of AD to healthy tissue. For the challenge data, the TIA for patient A’s liver, calculated with Equation 3, was 93.4% larger when lesions were included. Patient B, who had only 2 small liver lesions, showed a 6% increase. This highlights the need for agreement on the structures on which AD calculations are performed.

We acknowledge certain limitations in this study. First, the small number of patients (only 2) limits the generalizability of the findings, which are specific to ¹⁷⁷Lu-DOTATATE and may not apply to other therapies. In addition, this study focused on variability in dosimetry results, rather than assessing their accuracy. Task 4 aimed to minimize the impact of segmentation-related variability. However, quantities calculated in VOIs differed when using RTStructs and binary masks, depending on how the software participants used calculated VOI quantities from RTStructs, which define VOIs by polygons. It does point out the need to harmonize not just VOI interchange formats but also how VOI quantities are calculated from them. Another limitation was the spreadsheet design for task 4, which did not request activities or fitting parameters, data that would have helped identify errors and issues. Although the baseline analysis was conducted by a single dosimetry expert, the results and methods were extensively discussed within the team to minimize potential errors. We recognize the limitation of relying on 1 analyst but believe this was mitigated through Dosimetry Challenge team discussions. Voxel-level fitting was not explored because of the larger number of complexities that must be addressed during its implementation. Voxel-level noise and the requirement for precise voxel-level registration can introduce significant variability, but it was a challenge to pursue it at this stage when participants were asked to report organ-level values (although they could have calculated them with voxelized approaches). However, the TIA map provided in task 5 was created using voxel-level fitting. We observed that this approach resulted in artifacts in some voxels because of alignment issues. Also, to decide on the best model, the Akaike Information Criterion was used instead of the corrected Akaike Information Criterion, despite the small amount of data. Lastly, in all tasks, many submissions lacked sufficient detail on the fitting methods, parameters, and results, making them difficult to interpret or replicate. We emphasize the importance of systematic and detailed reporting of dosimetry methods to ensure reproducibility and clarity in the field.

Nuclear medicine societies have been actively engaged in ongoing efforts to establish recommended best practices for dosimetry workflows (22). There is also substantial existing literature discussing the performance of time–activity data fitting methods and their impact on TIA and AD variability (6–8). This study highlights an ongoing need for improvement in adherence to current best practices in fitting and integration in clinical settings. Considering the outcomes of this work and the findings from these referenced studies, we have identified 4 critical areas in which the field could improve by establishing harmonized practices that are feasible even in busy clinical environments.

The first area, assess fit quality, would ensure proper fitting functions through visualization and knowledge of the shape of the time–activity curve for organs and lesions for the particular agent (15,19). Although quantitative metrics of goodness of fit should be used and reported, relying solely on metrics such as the coefficient of determination to accept or reject a fit would seem unwise. Some outliers observed in our analysis would have been avoided by following this recommendation (Supplemental Figure 3).

For the second area, optimize pharmacokinetic modeling, the use of fitting and exponential functions, rather than trapezoidal integration, should be preferred when they reasonably fit the data. This is primarily because of the need to perform extrapolation from the last time point to infinity. Use of a model that prioritizes robust curve fitting, particularly during the slow washout phase, contributes most significantly to dosimetry (20). If data from multiple pharmacokinetic phases, including both uptake and washout, are included, the model should account for these phases to prevent inaccuracies in later phases of the curve. There might be cases in which acquiring an extra imaging time point would improve the calculations (e.g., a lesion with the uptake phase would benefit from the collection of points that allow modeling of that phase).

With regard to the third area, carefully validate time–activity curve estimation software, most outliers observed in this study were associated with in-house software, which may lack the robustness or validation of commercial packages. It is crucial to benchmark all software, especially in-house tools, using well-established datasets before the software is applied to new data. We recommend using literature-based datasets, such as ground truth synthetic data (15) or challenge data from patient A in Supplemental Table 2, as gold standards to test the fitting and integration performance. Alternatively, new software can be validated by comparing its results with those from well-established software.

The fourth area, implement sanity checks, would incorporate sanity checks into the dosimetry workflow to ensure data are decay-corrected to scan time, input and output units are accurate, and results are presented systematically and clearly. For more on structured reporting, see Lassmann et al. (23) and Sgouros et al. (24).