Investigation of the influence of sampling schemes on quantitative dynamic fluorescence imaging.

Dynamic optical data from a series of sampling intervals can be used for quantitative analysis to obtain meaningful kinetic parameters of probe in vivo. The sampling schemes may affect the quantification results of dynamic fluorescence imaging. Here, we investigate the influence of different sampling schemes on the quantification of binding potential (BP) with theoretically simulated and experimentally measured data. Three groups of sampling schemes are investigated including the sampling starting point, sampling sparsity, and sampling uniformity. In the investigation of the influence of the sampling starting point, we further summarize two cases by considering the missing timing sequence between the probe injection and sampling starting time. Results show that the mean value of BP exhibits an obvious growth trend with an increase in the delay of the sampling starting point, and has a strong correlation with the sampling sparsity. The growth trend is much more obvious if throwing the missing timing sequence. The standard deviation of BP is inversely related to the sampling sparsity, and independent of the sampling uniformity and the delay of sampling starting time. Moreover, the mean value of BP obtained by uniform sampling is significantly higher than that by using the non-uniform sampling. Our results collectively suggest that a suitable sampling scheme can help compartmental modeling of dynamic fluorescence imaging provide more accurate results and simpler operations.


Introduction
Fluorescence imaging (FI), as one of important branches of molecular imaging technology, has attracted much more attentions and provided an important tool for early diagnosis and staging of cancer [1][2][3]. Traditional static FI that acquires fluorescent images at the time point of the probe reaching a state of equilibrium can qualitatively describe the location and distribution of the tumor. Especially, combing with tomographic technique, it can provide accurate quantitative information of tumor, such as the three-dimensional localization and volume size. However, static FI cannot provide the dynamic status of the tumor. Combing with kinetics analysis, dynamic FI can estimate various parameters associated with metabolic rate and receptor density. Further combining with singular value decomposition (SVD) analysis, dynamic FI can separate specific signals from non-specific ones, which further improves the accuracy of quantifying the binding potential of a targeted probe [4].
However, sampling schemes may affect the quantification results of dynamic FI, including the sampling starting point, sampling sparsity, and sampling uniformity. To be exact, the timing of the probe injection should be used as the sampling starting point, which can be accomplished with the help of caudal venous catheterization [5][6][7]. This would significantly increase the complexity of experimental operation, especially for a larger number of animals used in the longitudinal study. It is necessary to affirm whether the acquisition at the timing of the probe injection is needed. Especially when we are only interested in some robust kinetic parameters, such as the binding potential (BP), a receptor density-related parameters. If we can get relatively robust dynamic parameters in case of missing some time points immediately after the injection of the probe, the experimental operations can be much simpler. In current dynamic FI, the sampling sparsity and sampling uniformity were usually empirically determined [8][9][10][11][12][13]. Too small a sampling interval results in heavier camera burden and data processing burden, but too large a sampling interval results in an unsatisfactory recovery of the time activity curve (TAC), which may affect the accuracy of the determination of kinetic parameters. Considering the decay of radionuclides, nonuniform sampling was usually employed in dynamic PET imaging to ensure the accuracy of the TAC and the quality of dynamic PET images [14]. In such scheme, intensive sampling points were adopted in the primary stage of probe metabolism, and sparse sampling points were used in the latter stage. However, such characteristic becomes different in dynamic FI, because the decay factor can be neglected due to the stability of fluorescent dye [15]. From this point of view, the uniformity and the scale of the sampling interval are probably worth investigating.
Here, we take the quantification of the BP values as an example to investigate the influence of different sampling schemes on quantification of dynamic FI. Logan graphical analysis with reference tissue model (GARTM) combined with SVD are used to calculate the BP values. Based on three groups of sampling schemes including the sampling starting point, sampling sparsity, and sampling uniformity, we design the corresponding investigation experiments. All of the investigations are conducted with the theoretically simulated and experimentally measured data.

Materials
Athymic female nude BALB/c mice (4-6 weeks, 18-22 g) were obtained from the Laboratory Animal Center of Fourth Military Medical University (FMMU). Cy5.5-GX1 synthesized by FanBo Co. Ltd. (Beijing, China) and its characterization had been done in our previous work [13].

Xenografted mouse model
All the animal studies were performed in accordance with the FMMU animal protocol. The tumor was introduced by injecting 2 * 10 6 BGC-823 gastric cancer cells into the subcutaneous area on the right limb of mice. Then it was allowed to grow for about two weeks before imaging. At that time, the tumor size was about 8-10 mm in diameter [16]. Ten mice were randomly divided into two groups (blocked group and unblocked group) according to whether the free GX1 was pre-injected before the injection of imaging probe (Cy5.5-GX1).

Dynamic fluorescence images
Dynamic fluorescence images were acquired using the IVIS Kinetic imaging system (PerkinElmer, USA). After collecting five background images, the Cy5.5-GX1 was intravenously injected with a dose of 2 nmol [17]. Then, we collected dynamic fluorescent images for 63 minutes immediately after the injection of the probe. The imaging interval and camera integration time were set to be 18 s and 1 s. During data acquisition process, the mice were anesthetized with isoflurane and a thermostat-controlled thermal heater kept their body temperature. In addition, the fluorescent images of each mice were acquired once again after 24 hours of the probe injection.

Kinetic models and parameter estimation
A two-tissue (three-compartment) model and a one-tissue (two-compartment) model were used to describe tracer kinetics in the tumor region and reference region, respectively [18]. Detailed description of kinetic models and related parameters can be found in [13]. Here, the GARTM model was applied to calculate the BP values, which can be expressed as Eq. (1) [19]: where ROI(t) and ref(t) represent the fluorescent activity in the tumor and reference regions respectively, k 2 ref denotes the washout rate from interstitial space into plasma in the reference region, DVR represents volume ratio of distribution and int ' is a constant. According to the definition of the GARTM model, the value of The BP values, as a macro-parameter reflecting the binding affinity and the available receptor density, could be calculated as BP = DVR-1 [20].

Simulations
Simulations were performed to simulate the uptake process of Cy5.5-GX1. By assuming that the plasma input function is similar to that of IRDye800CW-EGF, we chose the plasma input function of Cy5.5-GX1 from blood sampling experiments published previously [21]. Considering the molecular weight of Cy5.5-GX1 was about 2k Da [22,23], the kinetic rates constants K 1 and k 2 were set to be 0.013 min −1 and 0.13 min −1 . The k 3 and k 4 were simply set to be 0.07 min −1 , 0.1 min −1 in the unblocked group, and 0.02 min −1 , 0.1 min −1 in the blocked group. To better simulate the actual experimental data, 3% Gaussian noise was added to the simulated tracer uptake curves. The detailed simulation process can be found in the work of Hamzei et al [23]. The number of simulations for each group is 10,000 times. All the simulated time-activity curves (TACs) of fluorescence tracer were interpolated to 9 s intervals from 0 min to 63 min.

Sampling schemes
Three categories of sampling schemes were designed according to sampling starting point, sampling uniformity and sampling sparsity. In the scheme of sampling starting point, there are two ways for data acquisition, the deferred data acquisition and the truncated data acquisition. By assuming that the probe injection time was 0 min and the actual data acquisition time was T0 min, the deferred data acquisition ignored the data in the period of 0 min to T0 min and directly used the data in the period of T0 min to (60 + T0) min for quantification. However, the truncated data acquisition used the data in the period of 0 min to 60 min for quantification, and the uncollected data in the period of 0 min to T0 min was recovered using the linear interpolation method. For simplicity, the symbol D represents the deferred data acquisition and T denotes the truncated data acquisition. In this work, we set up six sampling starting points, including 0 s, 36 s, 72 s, 108 s, 144 s, and 180 s. In the scheme of sampling uniformity, there are also two ways for data acquisition, uniform acquisition and non-uniform acquisition. In the way of non-uniform acquisition, the data acquisition process was segmented into three sections. Similarly, the symbol U represents the uniform acquisition and N denotes the non-uniform acquisition. In the scheme of sampling sparsity, we designed four sparsity levels that used different number of sampling points, including 200, 100, 50, and 40. The four sparsity levels were simplified as S1 to S4. The detailed definitions of sampling schemes are summarized in Table 1. For better understanding, we take an example to illustrate the meanings of abbreviations. For example, the sampling scheme of UD-S1 (54 s) refers to the uniform acquisition in informality and the first level in sparsity (Table 1, 18 s*  200), and the deferred data acquisition was adopted with the sampling starting point of 54 s.

Image and data analysis
Identification of the regions of interest (ROIs) with the same surface area in the tumor and muscle region was accomplished with Living Image 4.5 software (PerkinElmer, USA). Dynamic fluorescence data pre-processing and analysis were performed using MATLAB 2015b (The MathWorks, USA). The region of muscle contralateral to the tumor was selected as the reference region because of its negligible specific binding [24]. The mean values of fluorescent signals (*10 6 photon/cm 2 /s) in each ROI were calculated after the average of five background images was subtracted, which were used to determine the TACs within the tumor and muscle ROIs. The singular value decomposition (SVD) method was used to process the acquired fluorescence data to un-mix the specific signals before applying GARTM [25]. In addition, all curves were interpolated to 9 s temporal resolution in order to match the sampling scheme.

Statistic
Pharmacokinetic parameters were expressed as their mean values and standard deviation (means ± SD, N = 5). SPSS software (IBM, USA) was employed to conduct all of the statistical analyses. Statistical significance was evaluated using one-way ANOVA, LSD, and Student's t-test. P values smaller than 0.05 were considered to be of statistically significant difference.

Influence investigation of sampling starting point
We first used simulation data to investigate the influence of sampling starting point on quantification of BP values. Considering the fact that the sampling uniformity and the sampling sparsity would not affect the conclusion, we here only summarized the statistical results of UD-S1 and UT-S1 groups. Figure 2(A) and 2(B) show the estimated BP values for different sampling starting points in the unblocked group and the blocked group respectively. We find that the BP value increased with an increase of T0 for both data acquisition ways. Furthermore, such increasing trend was much faster for the deferred data acquisition than that of the truncated data acquisition. We further statistically analyzed the simulation data, and find that significantly statistical differences can be observed from the second sampling starting point (LSD, P < 0.05) in the UD-S1 group. However, such significantly statistical differences were not observed until the third sampling starting point (LSD, P < 0.05) in the UT-S1 group. The estimated BP values for different sampling starting points in the real experiment were presented in Fig. 2(C) (unblocked) and 2(D) (blocked). The same conclusion can be addressed as that in the simulation. However, there were no significantly statistical differences among all sets of data in both the unblocked and blocked groups (ANOVA, P > 0.05). Moreover, we calculated the parametric maps that displayed the BP values at each point in ROI. Figure 2(E) and 2(F) present such parametric maps in the unblocked (2(E)) and blocked (2(F)) groups for the sampling scheme of UD-S1. No obvious differences were observed with the change of sampling starting point. To estimate the difference quantitatively, we calculated the image structural similarity index (SSIM) between two parametric maps. There was still no difference and the SSIM was 100%.

Influence investigation of sampling uniformity
Similarly, simulations were first conducted to investigate the influence of sampling uniformity on quantification of BP values. Figure 3(A) and 3(B) show the estimated BP values of different sampling starting points for the sampling schemes of UD-S1 and ND-S1. We find that the BP value in the UD-S1 scheme was slightly higher than that in the ND-S1 scheme. However, there was no statistical difference between each pair of data (t-test, P > 0.05). In other words, at the same conditions of data acquisition way, sampling sparsity, and sampling starting point, the BP value estimated by using the uniform sampling scheme was slightly higher than that obtained by using the non-uniform sampling scheme. Quantification results from the real experiment further confirmed this conclusion. Figure 3

(C) and 3(D)
show the estimated BP values from the real experiment for the unblocked and blocked groups respectively. Almost the same conclusion can be addressed from the real experimental data except that the standard deviations are a little larger. We also calculated the parametric maps for the sampling schemes of UD-S1 (54 s) and ND-S1 (54 s), as shown in Fig. 3(E) (unblocked) and 3(F) (blocked). No obvious differences were observed with the change of sampling uniformity, with the SSIM being of 100%.

Influence investigation of sampling sparsity
As shown in Fig. 4(A) (unblocked) and 4(B) (blocked), we first used simulation data to investigate the influence of sampling sparsity on quantification of BP values. Here, we selected quantification results from the sampling schemes of UD (54 s) and ND (54 s) as examples for results presentation. Two conclusions can be addressed. First, the standard variations of the BP values had a tendency to increase with the increase of sampling sparsity in both the unblocked (4(A)) and blocked (4(B)) groups. Second, the mean value of the BP values had a tendency to increase with the increase of sampling sparsity in UD (54 s) scheme. However, a trend of rise and fall with the increase of sampling sparsity was observed in ND (54 s) scheme, with the peak point at S3. Further statistical analysis proved that significantly statistical difference can be observed both in the UD (54 s) and ND (54 s) schemes (LSD, P < 0.05). Then, real experimental data was further used to verify these conclusions, as shown in Fig. 4(C) (unblocked) and 4(D) (blocked). Statistical analysis results showed that no significantly statistical difference can be observed (ANOVA, P > 0.05). Finally, the parametric maps (Fig. 4(E) for the unblocked group, and Fig. 4(F) for the blocked group) showed that no obvious differences were observed with the change of sampling sparsity, with the SSIM being of 100%. It was also obvious that when the other conditions were identical, the BP value obtained by the uniform sampling interval was higher than that of the nonuniform sampling interval and the BP value obtained by the deferred data acquisition way was higher than that obtained by the truncated data acquisition way. To further investigate the trend of quantitative BP values with the change of sampling sparsity level, we designed a series of sampling schemes with different number of sampling points. For the sake of simplicity, uniform sampling scheme was used here. This means that the selected sampling points were uniformly distributed on the sampling interval. The quantitated BP value as a function of the number of sampling points were shown in Fig. 5, where Fig. 5(A) and 5(B) show the results from the simulation data, and Fig. 5(C) and 5(D) are those from real experiment. We find that the variation of the mean BP value became flat when the number of sampling points was greater than 60, which was equivalent to a sampling interval of 1 min. On the other hand, the standard deviation of the BP value had a tendency to decrease with the increase of the number of sampling points. Results of the statistical analysis showed that statistical difference can be observed from the simulation data (LSD, all P < 0.05) for both the unblocked (Fig. 5(A)) and blocked (Fig. 5(B)) groups, but cannot be found from the real experimental data (ANOVA, all P > 0.05).

Discussion
The binding potential (BP) of the targeted probe is a very important pharmacokinetics parameter in tumor diagnosis and treatment, because it can reflect the expression level of the available receptors in the tumor [26]. With the ability of analyzing real-time changes of targeted probes in vivo, dynamic fluorescent imaging (FI) allows non-invasive quantification of the BP values. Incorporating the SVD analysis to un-mix the specific signal component of the tumor from un-specific component, dynamic FI allows in vivo quantitative results similar to those of dynamic PET imaging [4]. In this study, based on the GARTM model, we investigated the effects of sampling starting point, sampling uniformity and sampling sparsity on quantification of BP values. Currently, there are some limitations in the SVD analysis. Some factors may affect the quantification results, such as the size of ROI, imaging resolution, imaging depth, and signal-to-noise ratio. They are open issues that need to be further studied. Considering the fact that the TACs recovered by the SVD analysis were more consistent with the in vivo pharmacokinetics properties of the targeted probe and more similar to the statistical results of the BP values of raw data, only the statistical results of the BP values obtained by the recovered TACs were displayed.
When investigating two ways of data acquisition, no statistical difference can be observed for the real experimental data. The main reason may be that the variance of the data in the real test is too large. However, the quantitated BP values using the truncated data acquisition way is still higher than those using the deferred data acquisition are. This proved that the tail vein catheterization can be avoided to simplify the experimental operations without loss of accuracy of quantitating BP values. We only need to record the time of the probe injection and data acquisition of the first frame. This may lead to another problem. The location of the mice may change slightly before and after the probe injection, which may cause errors in quantification. However, such errors were far less than the signal of the tracer (the difference is multiple magnitude), so that it could be neglected.
When the other conditions were the same, the BP values obtained by uniform sampling scheme were higher than those obtained by the non-uniform sampling scheme were. For the GARTM model, the BP value is estimated by linear regression of sampling time points discretized in the TACs of the tumor and muscle ROIs. In theory, the greater the difference between the tumor and muscle TACs, the larger the estimated BP values. According to the in vivo pharmacokinetics properties of the probe, the difference between the tumor and muscle TACs over time will become increasingly large. Thus, compared with the non-uniform sampling scheme, more points were used in the uniform sampling scheme in the late stage when the total sampling points are identical. This is equivalent to further highlighting the proportion of tumor TACs, leading to a higher BP value. Currently, there are no gold standard methods for quantifying in vivo receptor concentrations. Each receptor quantitative model has their own advantages and disadvantages. For example, previous studies had found that the BP values estimated by the GARTM model was a little small and obtained by the Lammertsma simplified reference tissue model (SRTM) was slightly larger [23,[27][28][29]. Therefore, when the BP value is quantitatively analyzed, it may be possible to improve the accuracy of the BP value estimation by choosing a suitable combination of sampling scheme and kinetic quantitative model. Subsequently, the characteristics of several commonly used receptor quantitative model are summarized in Table 2 for reference purposes only [23]. It not only includes the two models abovementioned, but also summarizes the other two models of the reduced full reference tissue model (RFRTM) and the basis function method (BFM). Similar to the above reasons, the decrease of the number of sampling points makes the BP values become larger. The muscle TACs-dominant time quantum is much shorter than that of tumor TACs in the entire 63 minutes of dynamic imaging. The proportional reduction of sampling points is more detrimental to the muscle TACs, leading to the increase of the BP values. In the non-uniform sampling scheme, the number of sampling points in the earlier period is more than that used in the uniform sampling scheme. As a result, the growth rate of the BP values caused by the proportional reduction of sampling points in the non-uniform sampling scheme is significantly smaller than that in the uniform sampling scheme. When the number of sampling points is the same, the non-uniform sampling scheme in the Sparsity 4 cannot be shrunk with the same proportions as those used in the previous sparsity levels. Thus, we adjusted the segmented nodes of the sampling scheme (equivalent to a lot more in the early stages), which led to the decrease of the BP value. This result further proved the above explanation. On the other hand, the standard deviation of the BP values became larger as the number of sampling points decreases. This might be caused by the fact that the TACs were more affected by the noise in the integration. As shown in Fig. 5, the variation of the mean BP value became very flat when the number of sampling points exceeded 60. This indicated that we should use as many sampling points as possible to improve the estimation accuracy of the BP values. From this point of view, this is an advantage of dynamic fluorescent imaging compared to the current PET dynamic imaging. Higher sampling rate can be used for the dynamic fluorescent imaging to improve the estimation accuracy of the BP values while ensuring the quality of dynamic images.
The parametric maps, as shown in Fig. 2-4, were insensitive to the investigated sampling schemes. No obvious differences can be observed with the SSIM being 100%. These results are reasonable, because the influence of the investigated sampling schemes on each pixel of the image is consistent, which will not cause any structural changes. This indicated that we might not consider the investigated sampling schemes if we are only interested in parametric maps. In addition, it was also worth noting that the standard deviation of the BP values in the simulation was far less than that in the real experiment. This may be caused by the fact that the individual differences of animals were not considered in the simulations.
In this study, we investigated the influence of different sampling schemes on quantification of BP with theoretically simulated and experimentally measured data. The investigated sampling schemes included the sampling starting point, sampling uniformity, and sampling sparsity. Results showed that it is necessary to choose a suitable combination of the kinetic quantitative model and sampling scheme, which may improve the accuracy of the quantitated BP values and simplify the experimental operation in dynamic fluorescent imaging. For using the GARTM model to quantify the BP values, the preferred sampling schemes are the combination of the truncated data acquisition, uniform sampling, and the sampling sparsity level with the number of sampling points being more than 60. Under such sampling schemes, the procedure of tail vein catheterization can be avoided to simplify the experimental operation in dynamic fluorescent imaging. In summary, our results will provide a valuable reference for the quantification of the BP values with the dynamic fluorescence imaging.