ABSTRACT

Conclusion:

By combining the sparse Legendre representation of the renogram curves and the Moore-Penrose matrix inverse techniques, we obtained improved noise reduction in the deconvoluted data, leading to better elimination of non-physiological signals -as negative values- and the avoidance of the smear effect of conventional smoothing on the vascular peak, which both influenced the meanTT measurement.

Results:

For simulated data, the root mean square error (RMSE) between the theoretical non-noisy renal retention curve (RRC) and the results of the deconvolution methods applied to the noisy RRC were up to two times lower with LG (p<0.001). The RMSE of the reconvoluted renogram and the theoretical one was also lower for LG (p<0.001) and showed better preservation of the original signal. The SF was neither improved nor degraded by the proposed method. For patient data, no statistically significant difference was found between the SF for the LG method compared with the database values, and the meanTT better agreed with the physician’s diagnosis than the matrix or clinical software (Hermes) outputs. A visual improvement of the RRC was also observed.

Methods:

The stability and accuracy of the proposed method were tested using a renal database containing Monte Carlo-simulated studies and real adult patient data. Two clinical parameters, namely, split function (SF) and mean TT (meanTT), obtained with LG were compared with values calculated with the established method that combines matrix deconvolution and a three-point linear smoothing (F121) as recommended by the 2008 International Scientific Committee of Radionuclides in Nephrourology consensus on renal TT measurements.

Objectives:

This study aimed to introduce an improved deconvolution technique for Tc-99m-mercaptoacetyltriglycine renograms based on the combination of a sparse Legendre polynomial representation and the Moore-Penrose inversion matrix (LG). This method reduces the effect of noise on the measurement of renal retention function transit time (TT).

Introduction

Deconvolution methods were initially applied in the early 1940s (1) to blood flow measurements and later to gastroenterology kinetic tracers (2). The first quantitative renogram analysis by deconvolution was proposed in the early 1970s (3). Various renography tracers have been developed since then. Tc-99m-mercaptoacetyltriglycine (MAG3), which is the focus of present research, has become the radiopharmaceutical of choice in various clinical contexts (4) and has been routinely used for years. Furthermore, increased interest in renography is expected following the recent development of novel PET radiotracers (5) that open the door to three-dimensional renography, while scintigraphy remains limited to two-dimensional renograms.

Deconvolution is a mathematical process that can be viewed as an imitation of the renal retention curve (RRC), which could be obtained if the radiotracer activity was instantly and directly injected into the renal artery. The RRC provides information about the quantity of radiotracer retained in the kidney over time. It is possible to extract information that relates to the spread of tubular transit times (TTs) and has physiological significance. For example, urinary tract obstruction can lengthen the TT of the radiotracer through the kidney. The renal system can be described as a series of pathways with different lengths and therefore with different TTs. The parameters with the highest clinical relevance are the mean TT (meanTT) and the left-to-right ratio or split function (SF). SF has been proven (6) to be proportional to the plateau height of the individual RRC. After removing the early vascular phase by back extrapolation of the plateau, the meanTT is computed as the integral of the RRC divided by the plateau height. For MAG3, typical normal values ranged from 2 to 4 min for the whole kidney meanTT (4,7) and are between 45% and 55% for SF.

Mathematically, deconvolution can be seen as the solution of the integral equation:

where R (t) is the renogram obtained in the clinical routine, B (t) the input function (blood input), and H (t) the unknown RRC. Different methods have been intensively investigated to solve this equation, such as Laplace transform, matrix systems or constrained least square in the time domain and fourier deconvolution in the frequency domain (7). A study revealed that these deconvolution techniques performed differently in clinical cases (7). For example, the constrained least square was more reliable in assessing the meanTT in the presence of very noisy data. Conversely, the matrix method performs better when noise is very low (7).

Unfortunately, the deconvolution process is one of the mathematical inverse problems that are classified as ill-posed. This means that small errors or variations in the input data lead to large errors in the output, i.e., the deconvoluted curve. On some occasions, the process can lead to a solution RRC, which has no clinical meaning, while the RRC that is reconvoluted with the input function still gives the correct initial renogram. None of the previously mentioned deconvolution methods are exempted to this problem. Moreover, they are all very sensitive to small variations in the inputs, renogram, and blood input curve and are therefore very sensitive to data noise. For the matrix method, which is the most frequently implemented method, the initial signal value is crucial. This value determines not only whether the matrix can be inverted, but any slight change or error in this value will propagate throughout the deconvolution process (7).

Given the presence of a high level of noise in the majority of scintigraphy data, it is not possible to obtain an RRC by deconvolution without a reduction in these statistical fluctuations. Without appropriate noise reduction, large oscillations are often present in the solution, sometimes with negative values that have no physiological meaning. Simple filtering of the input data, a common method of noise reduction, is not very effective. Excessive filtering may modify the characteristics of the underlying physiological information present in the curves. On the contrary, applying a minimum filtering value will still result in an unstable and inaccurate output. Previous studies have shown that the accuracy of the RRC depends not only on the degree but also on the type of filtering (8). Moreover, filtering must be applied before (and sometimes also on) the output of the deconvolution process, mainly to remove negative values. The RRC contains a mix of renal and non-renal phases. It follows the idea that the initial shape of the RRC should be cautiously considered since the first time points (vascular stage) are higher than those that follow. One of the drawbacks of excessive filtering is a spread of the vascular component in addition to a greater difficulty in determining the end of the RRC plateau. A general conclusion from the literature, as summarized in the ISCORN consensus (4,7), is that linear filtering with a 1-2-1 kernel and with a variable number of passes, as proposed by Fleming (9,10), is the recommended option.

This study was conducted to determine whether sparse Legendre polynomials (LP) can be used to create a representation of time-activity curves (TACs) as a Poisson noise removal tool, as was recently shown to be effective for standard renography processing (11), which in combination with an adapted deconvolution technique could be suitable for improving nuclear renogram analysis by deconvolution in routine practice.

Materials and Methods

Representation of LP and the Pseudo-inverse Method

This study relied on the representation of a noisy TAC by sparse Legendre expansion (11,12). The description of the noisy TAC ƒ (t) can be written as a limited sequence of the LP of K_max +1 terms giving a denoized curve ƒ_L(t):

where  P_k (t) is the LP of order k and the Legendre coefficients L (k) are computed from:

The $\frac{2\mathrm{k}+1}{\mathrm{N}}$  term is a normalization factor linked to the LP.

As the tracer is injected into the blood, the renogram function R (t) (which is the TAC of the kidney obtained from scintigraphy images) is the mathematical convolution of the plasma or impulse function B (t) and the kidney impulse response function H (t). The relationship between these functions is given by the following equation:

R(t) = H(t) o B(t)   (4)

where o denotes the convolution operator. This process is schematically represented in Figure 1. B (t) is usually obtained from a region of interest (ROI) over the heart. H (t) describes the TAC that would be obtained if the tracer was injected in bolus (δ mathematical function) directly into the kidney artery and is also called the impulse response function. H (t)  does not depend on the shape of the blood input function and truly characterizes the fundamental behavior of the kidney.

Figure 1

In the first step, we applied the finite Legendre transform to the raw blood input B (t) and renogram R (t) curves to perform noise removal.                       

where ${\mathrm{L}}_{\mathrm{j}}^{\mathrm{B}}$ and ${\mathrm{L}}_{\mathrm{ı}}^{\mathrm{R}}$ are component vectors of dimension j and ι and KB_max, KR_max are the individual number of Legendre coefficients for each curve. ${\mathrm{P}}_{\mathrm{j}}^{\mathrm{B}}$ and ${\mathrm{P}}_{\mathrm{ı}}^{\mathrm{R}}$ are LP components of a P^B ∈ ℝ^jxm and a P^R ∈ ℝ^jxm matrix, respectively, where m had the dimension of the number of time points in the renogram.

In a second step, the unknown function H (t) was also developed as a Legendre expansion:

where the ${\mathrm{P}}_{\mathrm{i}}^{\mathrm{H}}$ are the LPs of a P^H ∈ ℝ^jxm matrix and ${\mathrm{L}}_{\mathrm{i}}^{\mathrm{H}}$ are the unknown Legendre coefficients of the renal retention function of dimension q. In our previous publication (11), the maximum number of Legendre coefficients for ${\mathrm{L}}_{\mathrm{i}}^{\mathrm{B}}$ and ${\mathrm{L}}_{\mathrm{ı}}^{\mathrm{R}}$ was 2* $\sqrt{\mathrm{m.}}$ Here, we found the limit for q to be 4* $\sqrt{\mathrm{m.}}$

Substituting equation (7) in (4) gave the following equation:

Convolution of B (t) with the ${\mathrm{P}}_{\mathrm{i}}^{\mathrm{H}}$ gave the matrix equation:

where (${\mathrm{P}}_{\mathrm{B}}^{\mathrm{T}}$ (t)) is the transposed convoluted matrix and L contains the unknown coefficients of the RRC.

To solve this system, we needed the inverse of the non-square matrix ${\mathrm{P}}_{\mathrm{B}}^{\mathrm{T}}$ (t). This was performed using the Moore-Penrose inverse. This generalized inverse is calculated using its singular-value decomposition and is noted (${\mathrm{P}}_{\mathrm{B}}^{\mathrm{T+}}$). The solution of equation (7) was obtained by the calculation of the Legendre coefficient vector L’.

L’ could be different from  L due to the Moore-Penrose inversion process, which acts as a least square minimizer (13).

From the L’, we obtained the filtered response function H’ (t) as a sparse Legendre expansion with q terms:

Kidney Database

For this study, we used a freely accessible online kidney database (www.dynamicrenalstudy.org) that contains Monte Carlo (MC) (14)-simulated studies and real adult patient data. The pharmacokinetics of the MC database is driven by a multi-compartmental model for Tc-99m-MAG3 and based on the first-order differential equation:

where  C_i(t) is the tracer concentration in compartment i at time t, V_i the volume of the ith compartment, and  r_j⟶_i the transfer rate constant from compartment j to i. The simulation incorporated up to 69 phantom structures to provide a realistic patient representation in terms of anatomical and pharmacokinetic characteristics. The renal cortex and medulla were modeled by delay functions to generate time distributions close to real patient TTs. The MC dataset comprised six studies based on the same phantom with two clearances and a total of 30 simulations. Each study represented a specific SF, and simulations were available for two levels of simulated injected activity (50 and 100 MBq), anterior and posterior views, and a reference study (RS) posterior view. RS was a simulation without any physical image-degrading effect such as noise, tissue background, attenuation, and scattering in both phantom and camera, giving the actual time variant tracer distribution (ground truth). Each study was based on the characteristics of its RS, and a difference in the kidney-to-skin distance was taken into account. Each simulated dataset consisted of a dynamic renogram acquisition of 120 frames of 10 s and 128×128 pixels. A summary of the MC study characteristics is presented in Table 1. For each study, we extracted the RRC from the RS to obtain a perfect theoretical noise-free shaped retention curve. The noise-free RCC and input function (B) from an ROI over the heart of the RS were convoluted to obtain the theoretical renogram (R).

Table 1

Renogram Processing

Original ROIs were drawn on studies with the Hermes renogram analysis (Hermes Medical Solutions AB, Stockholm Sweden) and were copied onto our software (11). The kidney background ROIs were drawn from the lateral direction, going from the lower to the upper pole to avoid the pelvi-ureteric activity. An ROI was drawn over the left ventricle to obtain the input function and another between the heart and kidneys for input function background subtraction.

For the LG and matrix deconvolution (15), we always applied the same method for the starting point of the input curve and the determination of the plateau of the RRC. The input curve peak time was taken as the zero time (7). The plateau of the retention function was calculated from the mean value (P_Mean) of the RRC curve between 1.2 and 2.0 min. This PMean replaced all the values from zero to the time point where the RRC fell below 85% of the PMean (16).

For the matrix method and theoretical data, we combined pre- and post-filtration (MPPF) using the F121 with a progressively rising number of passes up to 12 for pre-filtering to cover the range of passes determined by Fleming’s formula (9,10). Post-filtering was applied up to six passes, for each pre-filtering value, starting at the third point of the deconvoluted curve to avoid the influence of the vascular part. For the simulated data, we only considered F121 pre-filtering using the number of passes as determined by Fleming’s formula and an F121 post-filtering with one pass, again excluding the first three points. This method is denoted MFF. The blood input function was filtered by F121 with a variable number of passes.

We applied the LP on the input function and on the noisy renogram before LG. We used a method based on an autocorrelation of the Legendre coefficients (11) to determine the optimum number of coefficients to be used. We also varied the number of Legendre coefficients to check if there was no other best solution.

For the first part of this study, the theoretical renogram R and input function B were added with different Poisson noise realizations to obtain noisy curves R’ and B’. Using R’ and B’, the two deconvolution processes (LG and Matrix) were applied to obtain RRC’. In the first step, to obtain a global estimate of the accuracy of the deconvolution methods, the root mean square error (RMSE) between the results of the two deconvolution processes (matrix and LG) on the noisy curves RRC’ and the theoretical RRC was calculated. A convolution of the computed RRC’ with the filtered input function was used to calculate the RMSE with the theoretical noise-free renogram R to assess the preservation of the original signal. In a second step, the SF and whole kidney meanTT were computed from the RRC’ obtained by both deconvolution methods.

The second part of the study used the MC-simulated studies with the presence of all physical-degrading effects and with respect to the different noise levels. The SF and TTs obtained from the RRC with both deconvolution methods were compared. Thus, it was not possible to determine the best F121, and only Fleming’s formula was used with one F121 post-filtration (MFF).

The last part of this study was the confrontation of the outputs of the two deconvolution methods (LG and MFF), to which we added the outputs of the Hermes kidney analysis clinical software, applied to 31 patient studies from the database. These studies consist of MAG3 renograms recorded in a single acquisition (30 min) of 180 frames of 10 s in 128×128 pixels. The selection of 31 studies from the database was based on the clinical diagnostics and aimed to have a ratio close to 50% between normal and pathological kidneys. The clinical data, SF, and diagnosis were available from the database. The SF provided in the database was compared with the values obtained with the deconvolution and analysis methods used in this work. Three specialists for nuclear medicine (BW, IM, and FH) performed a blind test analysis of the patient studies with the Hermes renogram analysis program to establish if the whole kidney meanTT was pathological or not, and the results were compared with the outputs of the LG method and MFF to determine if a correlation exists between the clinical diagnostic and calculated meanTT. A bias was noted in the Hermes results of the meanTT that was discovered at the beginning of the study with the help of the MC-simulated studies. It was confirmed by the Hermes support team that the system starts the analysis at 20 s. Moreover, the Hermes software does pre- and post-linear filtering (F121) and an apodization of the results to avoid negative values. Other details of the matrix deconvolution process in Hermes were unavailable. The determination of the plateau was computed in Hermes using the first and second derivatives. While the 20-s difference for the starting point is not significant for diagnosis, it has an influence on the matrix deconvolution output. For patient data, we considered this difference when classifying the kidney TT as pathological or normal. This feature of the Hermes software was part of the motivation for developing our own matrix deconvolution software.

CHU UCLouvain Namur Site De Sainte-Elisabeth Hospital Ethics Comitte approval (number: 08/21) was obtained, and the requirement to obtain informed consent was waived for this study because of its retrospective design.

Statistical Analysis

Data were divided into three subsets (RS, noisy MC, real patient data), and two parameters were analyzed (SF and whole kidney meanTT). The RMSE was used to measure the deviation between observed and theoretical values. The mean and standard deviation (SD) were used to describe the spread of measurements. The slope, intercept, and coefficient of determination were obtained from linear regression for a pairwise comparison of the outcomes of the processing methods. A Bland-Altman plot was used to compare pairwise methods, including the 95% confidence interval for the limits of agreement. Statistical t-tests were performed at a 5% level of significance (p<0.05) using XLStat (version 2019.1.3, France).

Results

Theoretical Renogram and Input Function Added with Poisson Noise

Considering the matrix-deconvoluted theoretical noisy RRC’, the composite number of pre- and post-filtering passes for F121 that led to the lowest RMSE between the deconvoluted RRC’ and the non-noisy theoretical RRC was considered optimal. This was consistent with previous findings (8). This optimal composite number of pre- and post-filtering F121 passes was selected for all further comparison with the LG method. The RMSE between theoretical and computed values for the retention function and renogram curves are presented in Table 2. The LG method showed an RMSE up to two times lower than the MPPF for the RRC’ (p<0.001). In addition, the reconvolution error was the lowest for LG (p=0.001). The MFF was slightly worse in both cases but nevertheless close to the lowest MPPF RMSE. The better performance of the LG can be visualized in Figure 2a, where the MPPF method still had ripples and negative values, which have no physiological interpretation. The LG method preserved most of the original RRC. In particular, the vascular peak and subsequent points were smeared out when MPPF was applied, a phenomenon that gradually was increasingly present as the number of F121 passes increased in pre-filtering (Figure 2b).

Table 2

Figure 2

For the simulated data, SF values are given in the database, and the meanTT was computed from non-noisy RRC by Riemann integration of the curve. An excellent agreement was noted between the database SF values, and the results computed from the three methods (Table 3) with no significant difference (p>0.30).

Table 3

The meanTT calculated from the non-noisy RRC curves and the results from the RRC’ obtained with LG and the matrix method -with MPPF or MFF- were compared using Bland-Altman plots (Figure 3a,b,c). LG resulted in a bias closest to zero and the smallest SD, while the MPPF with the lowest RMSE and the MFF performed very similarly.

Figure 3

MC-simulated Studies

For the MC-simulated studies, we only considered the MFF for the matrix deconvolution. For each study, the simulated data were created with the same kinetic parameters (SF and meanTT) based on their respective RS. For the SF, we did not find any significant difference in the mean between the LG and MFF methods when compared with the expected values (p=0.63). The average value of the meanTT extracted from the RS used as reference was 4.0 min for the left kidney and 3.9 min for the right kidney. When comparing the mean of the meanTT for all posterior studies and all anterior studies (Table 4), a difference was observed between the two methods. This is illustrated in Figure 4a, b. LG was always higher and closer to the expected meanTT value than MFF. Considering the same view and same kidney, the value dispersion was also reduced with LG. In Figure 4, the posterior and anterior studies were separated to take into account the difference in the patient configuration, in the background subtraction for the kidney and the input function, and to highlight the signal-to-noise ratio worsening in anterior images.

Table 4

Figure 4

Patient Studies

For the SF, no significant difference was found between LG and the database values (p=0.48). A comparison of the whole kidney meanTT was performed between the Hermes clinical software, LG, and MFF (Figure 5a, b, c). A general good agreement was noted in the normal range of meanTT (3.5±1 min for hydrated patients following ISCORN), while discrepancies between methods were more frequently observed for kidneys with meanTT outside this normal range. A systematic bias was observed between LG and Hermes (bias: 1.1, SD: 1.4). The LG meanTT values were higher than that in Hermes (p<0.001). The same trend was observed between MFF and Hermes (bias: 0.6, SD: 1.1, p<0.0001). The difference between LG and MFF was less marked (bias: 0.1, SD: 0.7, p=0.18).

Figure 5

The blind test realized by the three physicians gave a discordant diagnosis on 4 of 62 kidneys for the determination of an abnormal TT. The discordance appeared for patients with results in the “gray zone” where the final diagnosis remains physician dependent (17). For the remaining 58 kidneys, we obtained discrepancies for only 10.3% with LG, 12.1% with MFF, but 24.1% with the Hermes system between diagnosis based on prolongated TT values and physician diagnosis.

Discussion

A common result for all the data -theoretical, MC-simulated, and patient- is the absence of significant differences between the LG and database values for the determination of the SF, a clinically important parameter.

On the theoretical data, the LG showed better preservation of the original information in the renogram and a better recovery of the theoretical RRC (Figure 2). For MC simulations, the meanTT was closer to the expected values with a smaller dispersion (Table 4).

Figure 2

Table 4

The patient database does not contain any information about TT contrary to simulated data where the ground truth is known. This was part of the motivation to include the blinded physician diagnoses in the study, although the diagnosis does not always give evidence of a slowed TT or obstruction, and it remains difficult to define since no gold standard exists (18,19). The results obtained with patients had the same trend as the simulation data, suggesting that real-world data processing with LG behaves in the same manner. Moreover, the diagnoses were more in agreement with the values obtained with LG. The display of the RRC curve was also visually improved by the LG (Figure 6) method when reconvolution of the RRC with the input function is used as a quality control process, this is an interesting property.

Figure 6

The filtration of the input curve and renogram by the LP instead of F121 before the matrix deconvolution was also investigated (data not shown). This globally resulted in an improvement of the deconvolution over F121 filtering. Nonetheless, the LG method used in this study presents even more improvements in the results due to the Moore–Penrose inversion matrix, which acts as a least square fit. Even if the matrix method is dependent on the first point, the filtration of the curves remains crucial (7). As pointed out by previous studies (9), F121 is not the best technique, but still remains the most used. Moreover, there is no assumption in the LG method on the form of the input function in contrast to some methods where the input function is modeled by mathematical functions (20). This LG method should not be seen as limited to kidney deconvolution but should be also applied to other types of dynamic study (work in progress) in nuclear medicine.

The LG method can be automated and has no extra computation needs compared with the matrix method.

In summary, the LP-based deconvolution appears to be a clinically feasible alternative to the classic matrix deconvolution for renogram analysis.

Study Limitations

The two main limitations of the study are the low number of variable parameters (SF and clearances) in the MC simulations and the lack of independently determined TT in the patient database. However, the results showed a similar trend in both simulated and real data to the benefit of our method.

Conclusion

The study demonstrates that the LG method is more stable on simulated data and further preserves the original information than the use of a traditional linear filter, like the usual and ISCORN-recommended 1-2-1 filter, combined with the matrix deconvolution. LG gave the best and a near-perfect correlation with the expected values for the SF determination of simulated and real patient data. With LG, the kidney meanTT was less influenced by the noise in simulated data. For patient data, we observed a better correlation between the medical diagnosis and the values obtained by the LG and a better visual rendering of the RRC curve.