A fast fully connected neural network-based method for random and scattered coincidence correction in polyhedral brain PET

Simulation model construction

Regarding the detector model, polyhedral brain PET offers higher detection efficiency and geometric adaptability compared to conventional ring-shaped brain PET, resulting in the processing of more data within the same timeframe [20,21,22]. To investigate fast random and scatter correction for polyhedral brain PET, a polyhedral brain PET detector model based on the great rhombicuboctahedron was constructed using the GATE simulation platform. The scintillation crystals selected were LYSO, with dimensions of 4 mm × 4 mm × 20 mm, separated by a 0.2 mm gap. The selection of 4 mm × 4 mm × 20 mm LYSO crystals with a 0.2 mm gap is based on established designs for high-resolution PET scanners [23]. This configuration offers an optimal balance between spatial resolution, detection efficiency, and depth of interaction (DOI) effects [20].

For the phantom model, to closely mimic the actual human brain structure, this study utilized a 3D brain phantom derived from 11 human brains developed by Aubert-Broche et al. [24]. This phantom features strict delineations of brain regions and material assignments, allowing for the classification into 12 distinct tissue regions based on pixel values [24]. The corresponding tissue regions and pixel values are listed in Table 1. Each brain phantom model is divided into 362 slices along the axial direction, with each slice measuring 362 pixels × 434 pixels. The voxel size is 0.5 mm × 0.5 mm × 0.5 mm.

Regarding the source model, to construct a high-fidelity digital source, this study employed an anatomy-constrained physiological modeling approach. Taking into account the spatial continuity of tracer uptake and the complexity of microstructures in real biological tissues, we convolved the high-resolution anatomical phantom with a point spread function (PSF). This process effectively simulates the biological partial volume effect at tissue interfaces, generating a radiotracer activity distribution map with natural gradient transitions, thereby avoiding non-physiological artifacts inherent in ideal step-function models [25, 26]. A linear mapping relationship between pixel values and activity was established. The radiation source consisted of γ-photons with an energy of 511 keV, emitted in opposite directions.

In terms of simulation parameters, the Monte Carlo simulation was configured to emulate a modern clinical TOF-PET scanner. The detector time resolution was set to 300 ps to represent state-of-the-art systems. An energy resolution of 15% at 511 keV was assumed. Data were acquired within a 400–600 keV energy window to reject scattered events. The system dead time was modeled as 110 ns in non-paralyzable mode, and a 5 ns coincidence time window was applied, consistent with established simulation protocols [27].

Simulation procedure

Based on the phantom and source parameters set in the previous section, eleven simulations were conducted. During these simulations, only the brain phantom and source models were varied, while all other settings remained constant. The centers of the phantom, source, and the polyhedral brain PET detector were all aligned at the coordinate origin. When γ-photons emitted from the source interacted with the scintillators, the simulation recorded eight specific pieces of information for each photon: (event ID, time T, deposited energy E, interaction coordinates (X, Y, Z), Compton scattering ID, Rayleigh scattering ID).

Neural network architecture

Based on the simulation results, coincidence recombination was performed on the recorded events using a coincidence time window to obtain the initial coincidence events. These coincidence events were subsequently classified into three categories based on Compton scattering occurrence and other information, and corresponding labels were assigned. For the input layer, assuming a coincidence event includes the following information: response event 1 consists of response time T₁, response energy E₁, response coordinates (X₁, Y₁, Z₁), and response event 2 consists of response time T₂, response energy E2, response coordinates (X₂, Y₂, Z₂). Thus, the input layer for a single coincidence event is represented as a 1 × 9 vector (T₁–T₂, E₁, X₁, Y₁, Z₁, E₂, X₂, Y₂, Z₂). The output layer is a 1 × 1 vector, where a value of 10 indicates a true coincidence, 0 indicates a scattered coincidence, and −10 indicates a random coincidence.

To enhance the quality of coincidence data within the training set, a two-stage true coincidence discrimination process was implemented. The first stage involved eliminating coincidence data where the annihilation point was outside a specified cubic volume, determined by calculating the annihilation point for each event. Specifically, the ten digital brain phantoms used for training, testing, and validation were all located within a cubic space defined by X ranging from −90.5 mm to 90.5 mm, Y ranging from −108.5 mm to 108.5 mm, and Z ranging from −90.5 mm to 90.5 mm.

The FCNN training served as the second stage of true coincidence discrimination. For dataset configuration, coincidence data processed from the first eight digital brain simulations were used as the training set, while data from the ninth and tenth simulations were used as the validation set and test set, respectively.

The FCNN was implemented in MATLAB R2025a. To optimize the model architecture and training parameters – including the number of fully connected layers, neurons per layer, initial learning rate, and L2 regularization coefficient – Bayesian optimization was utilized [28]. During this process, 5-fold cross-validation was embedded within every optimization iteration, providing the evaluation feedback necessary for the algorithm to iteratively update and select the optimal hyperparameters [29].

The Bayesian optimization process ran for 59 iterations, at which point the mean squared error (MSE) had converged. The final architecture of the FCNN consists of an input layer, five hidden layers, and an output layer. The input layer accepts a 1 × 9 vector, and the output layer produces a 1 × 1 vector. The number of neurons in the five hidden layers are 511, 217, 126, 18, and 1, respectively. A ReLU activation function is applied between each hidden layer. The optimized initial learning rate is 0.0016658, and the L2 regularization coefficient is 0.000024817.

For the K-fold cross-validation, the final result yielded an MSE of 54.313 ± 0.153.

True coincidence discrimination and random/scatter correction

True coincidence discrimination

As illustrated in Fig. 1, based on the trained network framework, the input layer of the training set is first imported for prediction, yielding a corresponding predicted value M for each coincidence event via FCNN. A lower threshold K is then set to filter coincidence events where M > K, which are presumed to be true coincidences. By comparing these predictions with the actual output layer of the training set and incrementally increasing K, the value K1 is determined when the NECR reaches its maximum. Subsequently, this optimal threshold K1 is applied to the test set to obtain the final predicted NECR.

Fig. 1.

Random and scatter coincidence correction method

Based on the training results, a correction method for RSC is proposed. The eleventh digital brain phantom is selected to evaluate the training performance. The specific methodology is as follows.

In time-of-flight positron emission tomography (TOF-PET), for a specific line of response (LOR), the system matrix element a_ij can be decomposed into the product of a geometric component G_ij and a TOF kernel function K_ij (Δt) [30,31,32]:

Since the types and predicted values M of each coincidence event in the training set are known, the distribution functions of true coincidences, scatter coincidences, and random coincidences with respect to the predicted value M can be denoted as T(M), S(M), and R(M), respectively. The correction function is defined as L(M), expressed by the following equation:

For an unknown coincidence event, its predicted value M can be obtained using the trained FCNN. By inputting this M into the correction curve function L(M), the contribution value L of the true coincidence for this event is derived. Based on this true coincidence contribution value L, RSC for the brain PET image can be performed. For a specific LOR, incorporating Eq. (2), it can be expressed as follows:

To evaluate the training performance, this study employed the single-slice rebinning (SSRB) method to reconstruct images from the coincidence data using the filtered back projection (FBP) algorithm [33, 34]. For 2D voxel assignment along the LOR, a Gaussian kernel function was utilized to distribute values around the annihilation point along the LOR, with a standard deviation of 5 mm. The reconstruction matrix was set to 362 × 434 × 362, with a voxel size of 0.5 mm. The FBP reconstruction matrix image without RSC is denoted as IMG1. The FBP reconstruction matrix image with RSC applied is denoted as IMG2. Furthermore, the FBP reconstruction matrix image generated solely from true coincidence data is denoted as IMG3.

To comprehensively evaluate the performance of the proposed FCNN-based RSC method, both attenuation-corrected (AC) and non-attenuation-corrected (NAC) reconstructions were performed.

For the AC scenario, an analytical attenuation correction was applied. Based on the known anatomical structure and material composition of the digital brain phantom, a theoretical attenuation map (μ-map) was generated. The attenuation factors were then calculated analytically via ray-tracing through this map and applied to the sinogram data prior to reconstruction [35].

It should be noted that the attenuation correction was performed as a standard pre-processing step and was independent of the FCNN model, which focused solely on random and scatter coincidence discrimination.

Evaluation metrics

True coincidence discrimination

The GATE simulation required 354 hours, and the deep learning model completed 1 720 000 iterations. The training results demonstrate high computational efficiency: when predicting 1 761 742 coincidence events, the fully connected neural network required only 9.733550 seconds, averaging 5.52495768 × 10⁻⁶ seconds per event. This indicates the method significantly conserves computational resources and enables rapid prediction for subsequent correction processing.

Figure 2 illustrates the variation of NECR with the lower screening threshold K1. As K1 increases, the NECR gradually rises, reaching a maximum of 353 987.40 when K1 is approximately −0.6. This phenomenon is explained by Fig. 3, which depicts the distribution of predicted values for the three coincidence types and the correction function L(M). The primary reason for setting a lower threshold is that scatter coincidences S(M) and random coincidences R(M) have significant proportions when M< −0.6, whereas L(M) is primarily distributed for M> −0.6. This allows for the easy removal of scatter and random coincidences with predicted values below −0.6. However, in the region where M> −0.6, many random and scatter coincidences remain, as their predicted values largely overlap with those of true coincidences, making them difficult to distinguish. Therefore, a secondary discrimination process was implemented.

Fig. 2.

Fig. 3.

As shown in Table 1, after secondary discrimination, the NECR of the test set increased from 35 485.01 to 43 990.17. The random fraction (RF) decreased from 44.65% to 33.85%, and the scatter fraction (SF) decreased from 40.60% to 29.01%. Meanwhile, 95.60% of the true coincidences were retained, preserving the core imaging data while significantly reducing random coincidences and thereby minimizing imaging noise interference.

Figure 4 displays the annihilation point probability matrix generated from the filtered coincidence data (predicted values M< −0.6) for brain slice 181. The image shows that the filtered annihilation points are predominantly distributed outside the brain slice (outside the red box), with almost no removal of scatter or random coincidences occurring within the brain slice (inside the red box). This demonstrates that the method effectively removes invalid coincidence data outside the brain region, enhancing the quality of the coincidence data for subsequent imaging processing.

Fig. 4.

Random and scatter coincidence correction

Given the overlap of predicted values for true, random, and scatter coincidences when M> −0.6, this paper applied Eq. (3) to perform RSC on this subset of data. As shown in Fig. 5, the CRC, SNR, and CNR were subsequently calculated.

Fig. 5.

The experimental results demonstrate that the proposed RSC method significantly optimizes image quality when synergistically combined with AC. Under AC conditions, the application of RSC elevated the CRC from 0.946 to 0.969, boosted the SNR substantially from 10.24 to 19.13, and improved the CNR from 0.402 to 0.496 compared to the non-RSC baseline. This confirms the method's effectiveness in purifying the signal through accurate modeling and subtraction of scatter and random events, where AC recovers true count magnitude while RSC removes spurious events, resulting in a synergistic enhancement of contrast recovery and noise suppression.

Conversely, under the NAC condition, the differential impact of RSC on separate metrics necessitates a nuanced explanation grounded in the physical definitions of each parameter. Specifically, the CRC saw a slight increase from 0.945 to 0.971; this marked improvement is a direct consequence of RSC's ability to restore true concentration differences by removing spurious signals that mask true lesions and blur image details, making CRC the most sensitive metric to signal purity. In stark contrast, the SNR deteriorated from 11.53 to 9.73 because the removal of scatter and random components leads to a net loss of total counts without AC compensation, thereby increasing statistical noise according to Poisson statistics (noise ∝ 1/√N). Consequently, the CNR exhibited an intermediate behavior, decreasing slightly from 0.569 to 0.521, as it represents the composite effect of competing factors where the gain in contrast is offset by the increase in noise.

Therefore, the distinct performance trends observed under NAC do not indicate algorithmic failure but rather reflect the inherent trade-off between physical accuracy and statistical quality, further underscoring that the proposed RSC method is optimally deployed in conjunction with AC.