Prostate cancer is the second most common cancer in men in the United States¹. Early diagnosis of prostate cancer allows for a number of treatment options. Current screening processes include Prostate Specific Antigen (PSA) level in blood. In case of abnormal results, Transrectal Ultrasound (TRUS)-guided biopsy is often recommended. Although the PSA level has high sensitivity, it suffers from low specificity and may cause unnecessary anxiety and cost. Interest in MRI for diagnosis has increased over recent years; it is observed that the tissue perfusion changes due to angiogenesis, the pathological condition often found in cancer tissue. We are interested in accurate estimation of the perfusion parameters: (i): $$$K_{trans} \left[min^{-1}\right]$$$ (volume transfer constant between blood plasma and extravascular extracellular space or transfer constant) and (ii): $$$K_{ep} \left[min^{-1}\right]$$$ (rate constant)² derived from DCE-MRI of the prostate.The purpose of this work is to develop a robust and efficient numerical optimization technique to find the (nonlinear) least squares estimates of $$$K_{trans}$$$ and $$$K_{ep}$$$ from time-resolved 3D DCE-MRI. The performance of this technique is compared to the conventional fitting strategy, Levenberg Marquardt Iteration (LM)³, as well as to the commercially available DynaCAD (DC) package.

For each voxel in the region of interest (ROI), the acquired tissue enhancement curve, $$$C(t)$$$, is fitted to the perfusion model, $$$P(t)$$$, by minimizing a cost function. The Tofts² model shown in Eq. 1 and Eq. 2 estimates the tissue perfusion. In this model, $$$I\left(K_{trans},K_{ep},t\right)$$$ is the tissue impulse response function, $$$P(t)$$$ is the estimated tissue perfusion, $$$AIF$$$ is the Arterial Input Function selected from the L or R iliac artery close to the prostate, and $$$\otimes$$$is the convolution operator. This model is restated in linear algebraic form, where $$$A$$$ is a convolution matrix and $$$H\left(K_{ep}\right)$$$ is the discrete-time exponential part of $$$I$$$. The cost function $$$J$$$ in Eq. 3 is the difference between the acquired and estimated perfusion, $$$C(t)$$$ and $$$P(t)$$$ respectively, which is minimized to find the optimum $$$K_{trans}$$$ and $$$K_{ep}$$$. $$$J$$$ depends on both $$$K_{trans}$$$ and $$$K_{ep}$$$, and thus we have a 2D optimization problem for each voxel. As shown in Eq. 4 a closed-form expression for $$$K_{trans}$$$ (dependent on $$$K_{ep}$$$) can be derived. By replacing $$$K_{trans}$$$ from Eq. 4 in Eq. 3, we reduce the dimensionality of $$$J$$$ from 2D to 1D which is shown in Eq. 5. This technique is called Variable Projection (VARPRO, or VP)⁴. The optimum $$$K_{ep}$$$ is efficiently determined via Golden Section search⁵. The average runtimes per iteration for one voxel using Matlab 2015b are about 0.53 msec and 1.22 msec for VARPRO and LM method, respectively.

$$I=I\left(K_{trans},K_{ep},t\right)=K_{trans}e^{K_{ep}t}=H\left(K_{ep}\right)K_{trans}\hspace{7cm}\text{(Eq.1)}$$

$$P(t)=AIF\otimes I=AI=AH\left(K_{ep}\right)K_{trans}\hspace{9.1cm}\text{(Eq.2)}$$

$$J\left(K_{trans},K_{ep}\right)=\|AIF\otimes I\left(K_{trans},K_{ep},t\right)-C(t)\|_{2}^{2}=\|AH\left(K_{ep}\right)K_{trans}-C(t)\|_{2}^{2}\hspace{1.2cm}\text{(Eq.3)}$$

$$\nabla_{K_{trans}}(J)=0\Rightarrow K_{trans}=\left(AH\left(K_{ep}\right)\right)^\dagger C(t)=B^\dagger C(t),where, B^\dagger=\left(B^*B\right)^{-1}B^*\hspace{1.7cm}\text{(Eq.4)}$$

$$J\left(K_{ep}\right)=\|BB^\dagger C(t)-C(t)\|^{2}_{2}\equiv-C^{*}BB^\dagger C\hspace{8.8cm}\text{(Eq.5)}$$

We acquire 55 time frames of the prostate using an accelerated 3D time-resolved DCE-MRI⁶ sequence (256×384×38, 0.86×1.15×3.0 mm³, frame time 6.6 sec). Perfusion parameters from pixels within the prostate were estimated with the: (i) the above VARPRO, (ii) LM, and (iii) DynaCAD. The residual defined as root mean square difference between $$$C(t)$$$ and $$$P(t)$$$, is calculated using a subset of the time frames, (7 to 55). In our preliminary study, we have performed the VARPRO technique on 15 cases which were matched with the results of LM and DC.

Fig. 1A shows a plot of $$$J$$$ as a function of $$$K_{ep}$$$ (Eq. 5), for a representative voxel, demonstrating how the optimum $$$K_{ep}$$$ can be readily found. Fig. 1B shows the fitted perfusion estimate, $$$P(t)$$$ to the acquired data, $$$C(t)$$$. Fig. 2A shows the residual, for VP and LM plotted versus iteration. Fig. 2B shows the residual versus absolute computation time. Fig. 3 shows the perfusion parameter maps for a representative prostate study (9x4 mm² abnormality in left mid peripheral zone suggests presence of significant prostate carcinoma; Brachytherapy performed) using VARPRO, LM, which are similar to maps produced with DynaCAD. Fig. 4A shows a plot of the normalized residual versus the number of time frames used as input data. As shown, after about 30 time frames the residual reaches a plateau. The pharmacokinetic maps according to a few of the subsets are shown in Fig. 4B with $$$K_{trans}$$$ on top and $$$K_{ep}$$$ at the bottom row.

To attempt to differentiate benign versus malignant tissue, optimum pharmacokinetic mapping of $$$K_{trans}$$$ and $$$K_{ep}$$$ from 3D DCE-MRI images of the prostate can be performed quickly using the VARPRO technique described here; the results are consistent with those generated using the traditional LM fitting strategy and DynaCAD. Using about 30 time frames (6.6 sec) seems to be sufficient for accurate estimation of perfusion parameters.