Single proton LET characterization with the Timepix detector and artificial intelligence for advanced proton therapy treatment planning

Objective. Protons have advantageous dose distributions and are increasingly used in cancer therapy. At the depth of the Bragg peak range, protons produce a mixed radiation field consisting of low- and high-linear energy transfer (LET) components, the latter of which is characterized by an increased ionization density on the microscopic scale associated with increased biological effectiveness. Prediction of the yield and LET of primary and secondary charged particles at a certain depth in the patient is performed by Monte Carlo simulations but is difficult to verify experimentally. Approach. Here, the results of measurements performed with Timepix detector in the mixed radiation field produced by a therapeutic proton beam in water are presented and compared to Monte Carlo simulations. The unique capability of the detector to perform high-resolution single particle tracking and identification enhanced by artificial intelligence allowed to resolve the particle type and measure the deposited energy of each particle comprising the mixed radiation field. Based on the collected data, biologically important physics parameters, the LET of single protons and dose-averaged LET, were computed. Main results. An accuracy over 95% was achieved for proton recognition with a developed neural network model. For recognized protons, the measured LET spectra generally agree with the results of Monte Carlo simulations. The mean difference between dose-averaged LET values obtained from measurements and simulations is 17%. We observed a broad spectrum of LET values ranging from a fraction of keV μm−1 to about 10 keV μm−1 for most of the measurements performed in the mixed radiation fields. Significance. It has been demonstrated that the introduced measurement method provides experimental data for validation of LETD or LET spectra in any treatment planning system. The simplicity and accessibility of the presented methodology make it easy to be translated into a clinical routine in any proton therapy facility.


Introduction
Proton radiation therapy is increasingly used for cancer treatment worldwide (Durante 2019). Currently, 107 proton or combined proton and ion facilities are in operation, and 61are under construction or at the planning stage, according to The Particle Therapy Co-Operative Group website (Particle Therapy Co-Operative Group 2022). The key motivation for applying charged particle radiation therapy are superior physical and radiobiological properties of the proton and ion beams compared to the photon and electron beams used in conventional radiation therapy. Foremost, the characteristic depth-dose distribution of proton or heavier ion beams (the Bragg Peak, BP) enables minimization of the dose deposited in healthy tissues surrounding the tumor while maximizing the physical dose and, thus, radiobiological effects in the tumor region (Durante et al 2017). A number of studies demonstrate reduced toxicity and dosimetric superiority of intensity-modulated proton therapy over intensity-modulated radiation therapy plans in specific cases (Jain et Yu et al 2022). Nowadays, proton beams are often used to treat pediatric tumors when it is critical to limit the harm to healthy growing tissues. Protons and ions are also preferable for tumors located near organs at risk, e.g. ocular cancer, central nervous system cancer, and skull base tumors (Hu et al 2018).
The radiobiological effectiveness of protons and, in particular, heavier ions with respect to photons is greater than one due to an increased ionization density of charged particles on the microscopic scale. The variable ionization density encountered in these beams is typically quantified by linear energy transfer (LET) (ICRU 1970), an average microscopic quantity, or lineal energy (ICRU 1983), a stochastic microscopic quantity. For example, a high-energy proton beam produces in water or tissue a mixed radiation field composed of primary and secondary protons and other secondaries with low-and high-LET components.
In the current practice of treatment planning for proton therapy, one still uses the constant relative biological effectiveness (RBE) of 1.1. The validity of this approach has been questioned and is contradicted by many radiobiological studies-see Paganetti (2014)  The state-of-the-art experimental approach for proton beam commissioning and validation, as well as for treatment plan verification, is the measurement of physical dose using either individual or an array of ionization chambers, scintillating screens, and CCD camera combinations. These approaches originate from conventional photon therapy and integrate the signal calibrated against relative or absolute dose. This is appropriate for physical dose verification in photon therapy but insufficient to verify LET distributions besides physical dose in mixed radiation fields produced by proton or ion beams. The lack of commercial methods for routine clinical LET evaluation in addition to Monte Carlo (MC) simulated LET distributions is a major obstacle to the routine implementation of LET-based quality assurance (QA) and control (QC) procedures that seem necessary based on the generally accepted concept of variable RBE in proton and ion beams.
This contribution implements a method based on a Timepix detector (see figure 1(a)) and artificial intelligence (AI) for evaluating LET spectra and LET D in clinical proton beams. Commercially available as miniaturized radiation cameras (Granja et al 2018b) the hybrid semiconductor pixel detectors Timepix bring (c) A schematic illustration of the Timepix sensor positioned during measurements presenting the incident angle β, which is the angle between the normal to the detector sensor and beam axis. (d) Timepix positioned on step motors of BluePhantom 2 , IBA (Germany) in air. (e) A schematic illustration of the experimental setup for measurements in water phantom including Timepix detector and waterproof cover. (f) Timepix placed inside a waterproof holder in water phantom during measurements of the mixed radiation field. energy-sensitive single particle tracking (Granja et al 2013, Heijne et al 2013, Ballabriga et al 2018 and quantumdosimetry imaging (Granja and Pospisil 2014) into the physical sciences with applications in space (Turecek et al 2011, Kroupa et al 2015, Stoffle et al 2015, Granja et al 2016 and medicine (Bisogni et al 2009, Jakubek et al 2010, Loo et al 2011, Martisikova et al 2011. Timepix was also applied to measure the relative dose with a wide dynamic range in clinical proton radiation fields (Stasica et al 2020b).
In this work, a simple experimental setup consisting of Timepix detector (see figure 1(a)) operated in a water phantom was designed, and measurements of energy deposition patterns of individual protons in air and water were performed. As most RBE models consider only primary and secondary protons disregarding other secondary particles (Wilkens and Oelfke 2004, Chen and Ahmad 2011, Carabe et al 2012, Wedenberg et al 2013, a convolutional neural network (CNN) was developed for efficient particle track recognition resolving protons from photon and electron events. For recognized protons, LET spectra and LET D values were computed and compared with those calculated by MC simulations. The presented work provides the first evidence that this approach fills the aforementioned gap and can provide measured position-and depth-dependent LET spectra and LET D .

Clinical proton beams at CCB Krakow proton therapy center
The Cyclotron Centre Bronowice (CCB) at the Henryk Niewodniczański Institute of Nuclear, Physics Polish Academy of Sciences, in Krakow (IFJ PAN) has been operating in Poland since 2016, offering proton beams for intensity-modulated proton therapy, as well as experimental and industrial purposes. Two treatment rooms are equipped with a rotational gantry with pencil scanning proton beams. By the end of 2022, about 800 patients have undergone radiation therapy in the facility. At CCB, the C-235 (IBA, Belgium) cyclotron accelerates protons up to 235 MeV. The beam is available within the energy range from 70 to 226 MeV, which corresponds to the beam range in water from 4.2 to 31.8 cm. An additional range shifter (RS) made of 4.2 cm thick PMMA can be applied to modulate the beam range. The lateral beam size varies from 3 to 15 mm (σ). The clinical beam intensity ranges from 1 to 300 nA, while it can be reduced to about 1 pA in experimental mode. For treatment planning purposes, Eclipse TPS version 13.6 by Varian Medical Systems (US) equipped with an analytical algorithm is used in CCB. It was commissioned against experimental data and further exploited for the commissioning of research MC systems, including fast FRED MC code (Schiavi et al 2017, Gajewski et al 2021 and GATE/Geant4 (Borys et al 2022). The model of the therapeutic proton beam at CCB was used in this work for MC simulations (see section 2.4).
2.2. Timepix for single particle tracking and quantum-imaging dosimetry The Timepix hybrid semiconductor pixel detector developed at CERN (Heijne et al 2013, Ballabriga et al 2018 consists of a semiconductor sensor bump-bonded to the highly-integrated ASIC readout chip. The detector provides a high-granularity matrix of closely packed pixels with independent signal readout electronics per pixel. The hybrid architecture of the pixel detector allows using radiation-sensitive sensors of different semiconductor materials (e.g. Si, CdTe, GaAs) of varying thicknesses (e.g. 100, 300, 500, 1000, 2000 μm). The Timepix chip (Llopart et al 2007) consists of a 256 × 256 pixels matrix of a total of 65 536 pixels. Each pixel has a size of 55 μm × 55 μm, providing a full sensor sensitive area of 1.98 cm 2 .
Timepix provides highly-sensitive quantum-dosimetry imaging of single particles (Granja and Pospisil 2014). Ionizing particles produce characteristic tracks in the pixelated matrix of the Timepix detector, forming clusters of adjacent pixels (see figures 1(b) and (c)) (Holy et al 2008, Granja et al 2013 which can be visualized online Pospisil 2014, Granja et al 2018b) and registered with spectral (energy loss) and tracking sensitivity in wide field-of-view (Granja et al 2018a). Detailed information is obtained on the deposited energy, track path across the sensor, time, position in the sensor, and direction of incidence of single particles (Granja et al 2018a(Granja et al , 2018b). The energy a particle deposits is given by the sum of the energies deposited in all cluster pixels. Knowing the energy deposition in the cluster and analyzing its morphology, it is possible to recognize the type of the particle In this work, The Minipix Timepix detector with Timepix1 chip generation and 300 μm silicon sensor (Granja et al 2018b) energy calibrated per pixel (Jakubek 2011) was used. The compact radiation detector (see figure 1(a)) has small dimensions (77 mm × 21 mm × 10 mm), weights 25 g and requires a single USB 2.0 connector for integrated power, control and readout. It operates at room temperature without the need for active cooling and, with a waterproof seal, can be immersed in water. The Timepix chip runs in frame mode and pixel-level configuration in energy (so-called time-over-threshold) mode (Llopart et al 2007). The dead time of the detector, i.e. frame readout time, is around 22 ms, and the data rate is up to 45 frames per second. The crossplatform software PIXET (Advacam, Prague) controls the data acquisition and provides online visualization.

Measurements
The measurements were performed at the CCB facility in a water phantom using a motorized BluePhantom 2 system by IBA (Belgium). Timepix was placed inside an in-house designed PMMA waterproof holder and mounted on step motors of the phantom (figures 1(d)-(f)), allowing its precise positioning with respect to the beam.
First, measurements in homogenous fields were performed to collect the training data for CNN proton identification model (see section 2.5) and assess the impact of the incident angle β on the LET measurement (see section 2.6.5). Measurements were performed in air with beam nominal energies of 70 MeV, 100 MeV, 150 MeV, and 200 MeV for the detector set at the isocentre at various angular positions, namely for incident angle β equal to 30°, 45°, 60°and 75°-see figure 1(b). Incident angle β is defined in the figure 1(c). The mean lateral size of the beam (σ) at the isocenter in the air for considered beam nominal energies, i.e. 70 MeV, 100 MeV, 150 MeV, and 200 MeV, is equal to 6.7 mm, 5.3 mm, 4.2 mm, and 3.0 mm, respectively (Gajewski et al 2021).
Second, the holder with the detector was positioned so that β = 45°, and the phantom was filled with water. Such angular position is preferred as it is an optimum between maximizing registration of elongated clusters and minimizing the scattering of protons on an aluminum case of Timepix as it was proved before (Stasica et al 2020b).
Step motors and standard software of the BluePhantom 2 allowed remote control of the detectorʼs position in the phantom. Next, measurements were performed at a depth of 149 mm for a proton beam of nominal energy 171.66 MeV with RS and nominal energy of 150 MeV without the RS. Both nominal energy and RS usage combinations exhibit the same beam range in the water of 156.6 mm, while the beam size at the measurement depth was 8.9 mm and 5.7 mm with and without RS, respectively. The mean lateral size (σ) at the isocenter in the air for the beam nominal energy of 171.66 MeV and with the inserted range shifter is 6.3 mm. Finally, measurements were performed for 150 MeV proton beam along lateral and longitudinal beam profiles. Measurements along lateral beam profiles were performed at a depth of 78.3 mm at distances of 30 mm, 60 mm, and 90 mm, and at a depth of 156.6 mm at distances of 30 mm, 45 mm, and 60 mm (see figure 4, top panel). Measurements along longitudinal beam profiles were performed in the beam axis and at a distance of 37 mm from the beam axis at depths of 30 mm, 120 mm, 156.6 mm, and 161 mm (see figure 5, top panel). The measurement locations were chosen where different LET D values were expected based LET D distributions presented in figures 4 and 5, top panels.
The beam current had to be kept at an optimal level to avoid overlapping of the registered clusters, which could disturb information about the actual clusters' morphology. For the measurements performed out of the beam core, the lowest stable clinical beam current of 1 nA was appropriate. However, to avoid particle overlapping when measuring directly in the beam axes, for both setups in air and the water phantom, the cyclotron was operated in the physical mode to apply unstable cyclotron dark current. The physical operation mode of the cyclotron does not allow for precise registration of the delivered dose rate. Nevertheless, the LET spectra are not expected to be different as the probability that two or more protons are coming through the biological cell simultaneously is close to zero (Oancea et al 2023).

Monte Carlo simulations
GATE in version 9.0, compiled with Geant4 v. 10.6.2 code, was used for MC simulations. The physics list was QGSP_BIC_HP_EMZ (Zacharatou Jarlskog and Paganetti 2008, Winterhalter et al 2020) with the maximum step limiter set to 1 mm. The production cuts of 1 mm were applied for gammas, electrons, and positrons, and 10 μm for protons. An applied beam model was specific for the treatment room and included the energy-dependent parameters of the beam, such as initial energy and energy spread, beam lateral size propagation in air, separately in X and Y directions (emittance parameters), and monitor units to the number of primaries conversion. Moreover, the beam model includes the RS parametrization, i.e. the dimensions, position with respect to the isocentre (fixed nozzle in case of CCB), as well as the composition and density. The beam model in GATE is implemented using polynomials describing the dependence between the parameters and the nominal energy of the beam for the entire range of energies clinically available, i.e. 70-226 MeV. It was prepared by converting the extensively validated beam model implemented in previous work for FRED MC code (Gajewski et al 2021, Borys et al 2022).
The Monte Carlo simulations of the Timepix response were performed in geometry that mimics the experimental setup in water. The geometry consisted of a 14.08 × 14.08 × 0.30 mm 3 silicon box representing the Timepix sensor, a PMMA RS when needed, and a 40 × 40 × 40 cm 3 water phantom box. The ionization potential of water was set to 78 eV. The silicon sensor was positioned in water at β = 45°and various depths and distances from the beam axis corresponding to performed measurements. The GATE Phase Space actor was used to record the type, position, direction, and total deposited energy (including secondaries produced in the sensor) for each particle entering the silicon sensor. Only the particles entering the sensor at its front surface with respect to the beam were considered. The track length of each incident particle was calculated based on the sensor thickness and the incident direction, namely the incident angle β. Based on the total deposited energy and the track length, the LET of each particle was calculated -see section 2.6.1. Next, based on particle type information scored with a Phase Space actor, primary and secondary protons were separated from other secondary particles.
To obtain sufficient statistics of particles depositing energy in the detector, various numbers of primary protons were simulated depending on the detector distance from the beam and ranged from 2 × 10 3 in the beam axis to 2 × 10 7 in lateral regions.
2.5. AI-based particle identification 2.5.1. CNN model Particles were newly identified with AI methods using a deep CNN, because of their great success in handling visual data (Krizhevsky et al 2012, LeCun et al 2015, Kalyani and Janakiramaiah 2021, Sathiyapriya and Shanthi 2022. In particular, we selected a deep CNN (VGG-16 network-a 16-layered model named for the Visual Geometry Group in Oxford) that has shown consistently good performance in broad applications while addressing the problem of network depth in CNN (Simonyan and Zisserman 2015). Network depth refers to the number of layers, and in general, more layers can allow better identification, but at the cost of significantly increasing the parameters that must be solved for, and thus the need for careful training with more data. VGG-16 solves the problem by simultaneously increasing the number of convolutional layers and reducing the size of the convolutions in the layers to a 3 × 3 kernel. The multiple layers of small convolutions keep the data manageable while the multiple small convolutions function like a larger convolution, resulting in an accurate CNN for recognition, classification, and localization (Qassim et al 2018). The model is still large (138 million parameters) requiring it to be initially trained using other images from ImageNet database (Deng et al 2009). Then, the model was retrained on our data (transfer learning) using the Keras software package (Chollet 20152015). The majority of our data was from homogeneous fields, and the particles or photons were known apriori (labeled data). Training and validation were done using 10-fold cross-validation (Browne 2000, Jiang andWang 2017).
The labeled data set was provided by Advacam and was composed of 2899 816 clusters for protons, 853 717 clusters for electrons, and 3143 149 clusters for photons (x-rays, gamma rays). Data were collected in well-defined radiation fields provided by radionuclide sources and particle accelerator beams of light and heavy charged particles: 4-25 MeV electrons from the MT-25 Microtron accelerator at the Nuclear Physics Institute (NPI), Czech Academy of Sciences (CAS), Rez near Prague, 8-35 MeV protons from the U-120M cyclotron of the NPI-CAS Rez, and 70-200 MeV protons from the C-235 cyclotron at CCB IFJ PAN (see section 2.3 and figure 1(b)). X-rays and low-energy gamma rays were measured using 241 Am source at Advacam. Gamma rays were measured using 60 Co and 137 Cs radionuclide irradiators at the Czech Metrology Institute (CMI) Prague. Description of these calibration measurements and the data acquired with a similar Timepix detector with a 300 μm silicon sensor are given elsewhere (Granja et al 2013, 2018a, Stasica et al 2020a, Granja et al 2021, Nabha et al 2022. Data was selectively filtered to remove background and unwanted events such as overlapping clusters, clusters produced by secondary particles, scattered primary particles, or partially collected events, e.g. at the sensor edge and detector artifacts. Since there were more than three times as many protons as electrons in the entire data set, a generative adversarial network (GAN) was used to augment the electron clusters (Goodfellow et al 2014, Salimans et al 2016. To avoid mode collapse in GAN training, instance noise was used (Mescheder et al 2018).
The data set was randomly divided into training and testing using 10-fold cross-validation (Trippa et al 2015), while electron clusters created by GAN were used for training only. The CNN was trained for 100 epochs toward convergence (see figure 2-phase 1). Accuracy and uncertainty were assessed using confusion matrices (Fawcett 2006). After the trained model was established, the mixed radiation field data set was used as the input to identify protons from other particles (see figure 2-phase 2).

Proton recognition uncertainty
The AI uncertainty is the variance of the predictive probability. The predictive uncertainty here is the sum of the aleatoric and epistemic uncertainty (Hora 1996, Kiureghian andDitlevsen 2009). The aleatoric uncertainty gives the inherent randomness of an output (choice of the proton), and the epistemic uncertainty is the variability of the input data (Hullermeier and Waegeman 2021). The variance was calculated from the second to last layer. The chance of the output in percentage was also calculated. It is the chance that a particle could have been something else: if the output is a proton, the chance of it could be an electron or photon, while if it is an electron or photon, the chance of it could be proton only.
2.6. Data analysis 2.6.1. Unrestricted linear energy transfer in silicon The LET of protons traversing the pixelated sensor volume of the Timepix detector was computed using the same formula for both experimental and simulation data. The LET was computed as the ratio of energy deposited by a particle, ò, and its track length, l, in the silicon sensor volume: For measurement data, the track length l was computed assuming that the particle does not stop in the sensor volume based on the cluster length in the sensor plane and sensor thickness which is 300 μm.
2.6.2. LET conversion from silicon to water Typically, the conversion of LET in silicon to water or tissue is carried out using a constant conversion factor (Bolst et al 2017, Bolst et al 2021, Pan et al 2021. Herein, instead of a single value, a fitted function proposed by Benton et al (2010) was applied, which covers a wide range of energies, and suits the purpose of this work.
2.6.3. Dose-averaged linear energy transfer (LET D ) While the full characterization of the mixed radiation field in terms of LET requires information about the entire LET spectrum in a given point, the average value of LET, most often LET D , is commonly provided instead Figure 2. A schematic illustration of AI-based particle identification method presented in this work. In phase 1, a pre-trained model was trained with labeled data measured with Timepix in well-defined radiation fields. In phase 2, a model trained with Timepix data was applied to unlabeled data measured in a mixed radiation field produced by a therapeutic proton beam in a water phantom. The output was labeled data with recognized protons. (Kalholm et al 2021). Herein, LET D was computed on a proton-by-proton basis using the following formula: where i is the particle index.
2.6.4. LET uncertainty LET uncertainty was determined by the error propagation, considering uncertainties of energy deposition, sensor thickness, and particle incident angle β. Energy deposition uncertainty was computed based on the fit using the minimum chi-square method to the experimental data, namely the x-ray fluorescence signal of several materials (e.g. Cu, Fe) and decay peaks of several radionuclides (e.g. 58.4 keV-estimation of Am-241 x-ray decay peak at 59.5 keV, 353.4 keV-estimation of Eu-152 β decay peak at 344.3 keV). Sensor thickness uncertainty estimated by the manufacturer is ±10 μm, and the angle uncertainty estimated based on MC simulations is 2°. The uncertainty of LET conversion from silicon to water is not provided (Benton et al 2010).
LET uncertainty was assumed to be normally distributed and to cover the 3σ range. To account for variable LET uncertainty for individual protons and to transfer it to the LET spectrum, a variable kernel density estimate was approximated by sampling 1000 points for each bin of the LET spectrum. Then, a resulting normal distribution was weighted for each bin by the uncertainty of the number of counts, which is a consequence of applying the particle recognition model.

LET spectra verification in air.
For each measurement in air (see section 2.3 and figures 1(b) and (c)) protons that hit the sensor plane at angles other than β ± 3°were filtered out. Next, LET values in silicon were computed as presented in section 2.6.1 for all protons. A histogram for each measurement was produced and normalized to the maximum value. LET spectra in silicon for measurements in the air for each nominal beam energy and various angles were compared.

LET spectra verification in water.
For verification of the measurement reproducibility, the results of two independent measurements of LET spectra in the mixed radiation field produced by a proton pencil beam of nominal energy of 171.66 MeV with RS and 150 MeV without RS in the beam axis at a depth of 149 mm were compared. As the range of these beams in water is the same, the LET spectra are expected not to differ qualitatively in these points. The particle identification model was applied to separate protons from other particles as presented in section 2.5. Next, LET in water was computed for each proton as presented in sections 2.6.1-2.6.2 and the LET uncertainty was determined as presented in section 2.6.4. Histograms for both measurements were produced, normalized to the maximum value, and compared with spectra obtained from MC simulations.

LET spectra along beam profiles in water.
For the analysis of the LET spectra, the particle identification model and LET silicon-to-water conversion was applied as described above. Also, LET uncertainty was computed the same way. Then, a histogram for each measurement was produced and normalized to 10 8 primary protons corresponding to a physical dose of about 1 Gy in the BP. Obtained spectra were compared to the corresponding MC simulation results. For each measurement, LET D was computed for protons only, as presented in section 2.6.3, and compared to LET D obtained from MC simulations.

Results
3.1. AI-based particle identification Accuracy, sensitivity, and specificity were used to evaluate the performance of the CNN model (Fawcett 2006). Sensitivity and specificity were particularly important as they are considered adequate predictors for imbalance classification (Luque et al 2019). The accuracy achieved for proton recognition is 95.36%. The resulting neural networks have a sensitivity of 92.97% and a specificity of 97.10%. The CNN model was applied for mixed radiation field data sets to identify the particle type, after which LET spectra were computed for recognized protons and compared with MC simulations.

LET spectra verification in air
A comparison of LET spectra measured in the air for proton pencil beams of various energies and different incident angles is presented in figure 3(a). The shift of the LET spectrum was observed at lower beam energies and larger incidence angles. The largest relative difference of the most probable value with respect to the one for 45°was at the level of 8%. Although the difference is considerable, the measurements and simulations in water show that the percentage of protons hitting the sensor at the angle β 75°for all mixed field data measured is lower than 2%.

LET spectra verification in water
A comparison of LET spectra measured and simulated in water for the proton beam nominal energies of 150 MeV without RS and 171.66 MeV with RS at a depth of 149 mm in the beam axis is presented in figure 3(b)) presents . As the beam range is the same in both scenarios, the LET spectra in the considered points should match. As can be seen, measured and simulated LET spectra of recognized protons overlap. For measurements performed at a depth of 78.3 mm LET D values range from 2.5 to 2.8 keV μm −1 . In the case of the measurement location most distanced from the beam core, the largest relative difference in LET D between measurements and MC simulations equal to 12% occurs. For measurements performed at a depth of 156.6 mm LET D values range from 5.5 to 9.1 keV μm −1 . At a distance of 30 mm and 60 mm, there are large discrepancies between LET D obtained from measurements and MC simulations at a level of about 50%-80%. This can result from long tails in the high end of LET spectra for measurements not present in MC simulation results (note log-log scale on both axes).

LET spectra along beam profiles in water
The relative number of protons decreases with distance from the beam axis. Positions of the peak obtained in measurements and simulations are in good agreement for small distances, although a shift occurs in the case of farther points. Also, a wide range of measured LET values for protons from fractions to more than 10 keV μm −1 was observed. For example, at a depth of the beam range 45 mm from the beam core, there is a significant number of protons of LET of 1.5 keV μm −1 , as well as of 10.0 keV μm −1 , while LET D is equal to 5.5 keV μm −1 .
Measured and simulated LET spectra for points along longitudinal profiles of 150 MeV proton beam in water in the beam axis and at a distance of 37 mm from the beam axis are presented in figure 5. The results are presented for measurements and simulations at depths in the water of 30 mm, 120 mm, 156.6 mm and 161 mm for both profiles. Again, vertical lines on LET spectra present LET D values obtained from measurements and MC simulations.
For measurements performed in the beam axis LET D values range from 1.1 to 5.1 keV μm −1 at depths in water from 30 mm to 161 mm respectively. The greatest relative difference between measurements and simulations is at a level of 6%. For measurements performed at a distance of 37 mm from the beam axis LET D values range from 2.7 to 6.5 keV μm −1 . The greatest relative difference between these measurements and MC simulations occurs at a depth of 30 mm and is equal to 20%.
The relative number of protons does not change significantly with depth but is several orders of magnitude smaller at a distance of 37 mm with respect to measurements in the beam axis, as expected. The peak positions in measurements and simulations are overall in good agreement, but discrepancies were observed in the plateau region. The mean relative difference between proton LET D values obtained from all presented measurements and corresponding simulations in the mixed radiation field is 17%.

Discussion
Results of LET spectra measurements and LET D calculations in air and water obtained with therapeutic proton beams were reported and compared with MC simulations. We demonstrated that for high-energy protons, the particle incident angle does not substantially impact LET estimation (see section 3.2 and figure 3, top panel). We also proved that primary proton beam of 150 MeV without RS and 171.66 MeV with RS adjusted to have the same BP depth lead to matching LET spectra (see section 3.3 and figure 3(b)). LET spectra and LET D values distributions measured at different positions and depths in water phantom are presented and compared to MC simulations in section 3.4 and figures 4-5. While LET spectra measured close to the beam axis agree well with MC predictions, the spectra measured at larger distances from the beam axis differ slightly in shape and position with respect to MC simulations.
In these examples, we demonstrated that Timepix technology could provide accurate measurements of single particle LET that generally fit the MC predictions, although some differences in the low and high end of LET spectra were observed (see discussion below).
A wide range of LET values of single protons was observed for most of the investigated measurement locations, especially in the BP region. For intensity-modulated proton radiation therapy, composed of several heterogeneous radiation fields and multiple pencil beams, the LET spectrum will be even wider. In particular, the high-LET component will be responsible for elevated radiobiological effects (Rucinski et al 2021). While overall, high-LET is recognized to be correlated with radiographic image changes, there is no consensus on the impact of LET D on this (Bahn et al 2020, Garbacz et al 2021, Niemierko et al 2021. In such cases, LET D as an averaged quantity may not be an appropriate predictor of the radiation quality (Grün et al 2019). Note that, for example, at a depth of the BP, 45mm from the beam core (see figure 4(E)), the single particle LET values range from 1.5 keV μm −1 to 10.0 keV μ −1 , while calculated LET D equal to 5.5 keV μm −1 may not be representative for such broad spectrum. Thus, the knowledge of the spectrum of proton LET values is essential for LET-based plan optimization.
In this work, the CNN was applied to identify particles of different types, namely protons, photons, and electrons registered with the Timepix detector. This is a novel application of AI for nuclear physics data analysis useful also in proton radiotherapy. Non-AI-based methods are widely used for particle type recognition in Timepix data (Holy et al 2008, Vykydal et al 2009, Bouchami et al 2011a, Hartmann et al 2017, Granja et al 2018b. However, such approaches require manual analysis of cluster features, to assure precise classification of different types of particles. The most time-consuming part of this process must be repeated for different particle types, experimental setups, and sensor thicknesses or materials. Moreover, we risk omitting relations between particle type and cluster features that we are not considering in our analysis. Here we present an AI-based approach as an alternative. The CNN algorithm finds the strongest patterns itself analyzing 2D maps of cluster images. To improve its performance or to include more particle types, more labeled (training) data needs to be provided, and the performance of the model is restricted only by the amount and quality of the training data. Observed differences between measured and simulated LET spectra in water may result from the convolution of several possible causes associated with MC simulations, the measurement method, and the particle identification model. Although not investigated, they may be related to the choice of physical models applied in the Geant4 simulations-e.g. shifts in the figures 4(C), (F) and differences in the high ends of LET spectra in the figures 4(D), (E) and 5(F), (G), (H). Other reasons could be simplifications of the beam model and geometry implemented in the MC simulation code (see section 2.4)-e.g. the geometric details of the gantry nozzle and the aluminum case of the detector were neglected. Also, the pixelated structure and charge-sharing effect , Bouchami et al 2011c, Teyssier et al 2012, Šolc et al 2022 occurring in the Timepix sensor were not included in the MC simulations. However, it is unlikely that these details have a substantial impact on the MC simulationʼs accuracy.
The discrepancies between measurements and MC simulations in the high LET region (see figure 4 panels D and E and figure 5 panels G and H) are the result of different track length estimation approaches. The measured 3D track length is estimated based on the 2D cluster length and the sensor thickness (see section 2.6.1). Contrary, the track length in the MC is estimated based on the particle incident angle and the sensor thickness (see section 2.4). This difference leads to the overestimation of the track length in MC with respect to the measurement, i.e. underestimation of the LET. We are currently investigating possible solutions to overcome this intrinsic limitation of the Timepix detector, by scoring the true track length in MC and the charge-sharing application, similarly to the methodology described in Šolc et al (2022) or by applying detailed simulations of the charge transport inside the detector using Allpix 2 toolkit (Spannagel et al 2018, Spannagel andSchütze 2022).
An increased LET uncertainty in the low and high ends of measured LET spectra and further from the beam axis is associated with lower statistics for these data and applying a logarithmic scale. For most of the data, we observe a low or negligible uncertainty originating from the AI model, energy, and angle measurement described in section 2.6.4. The impact of the uncertainty resulting from AI is limited due to large labeled data sets applied to train the CNN. The high sensitivity and specificity parameters also confirm the excellent performance of the AI model. A limitation of the presented approach is the LET computation method used when analyzing the Timepix data. This method relies on the assumption that protons cross the Timepix sensor, ignoring stopping protons near the end of the beam range (Nabha et al 2022), which leads to an overestimation of the track length and the associated underestimation of LET, particularly at the end of beam range and lateral regions. This issue will be further investigated with MC simulations and could be addressed by, e.g. applying a correction factor. Even though the water phantom allows very accurate positioning of the Timepix in water, positioning uncertainties might occur at a level of 1 mm. The future goal is to estimate LET D robustness to LET spectra measurements with Timepix, which can be, for instance, assessed by analyzing the difference between LET spectrum obtained from Timepix and MC simulations.
Further experimental uncertainty is associated with LET calculation for low-energy particles impinging the detector at large angles. The measurements in the air show a slight discrepancy of the LET spectra, which increases for lower beam energies and larger incident angles (see figure 3.2). Although a contribution from particles with this characteristic is minor (2% of protons hit the detector at 75°), further measurements in the air are planned in CCB to estimate the impact of particle incident angle β on LET measurements, which is expected to be greater for lower proton energies.
The limitations of AI in this work include the class imbalance problem which occurs when the distribution of classes in the training data set is unequally distributed and negatively affects the CNN (Buda et al 2018). In this work, the number of electron clusters in the training data set was significantly less than for the other two particle types. This problem was addressed using GAN to generate more electron clusters for training the CNN model. There is also a lower number of protons compared to the MC simulations at the low end of the LET spectra (see figures 4(A), (B), (D), (F) and 5), which could have been misinterpreted as electrons and incorrectly excluded by CNN model. Another challenge of the training data set was that the beam nominal energy range differed for each type of particle in the training data set and the mixed radiation field.
For future development of the CNN-based particle recognition model, more data will be supplied to balance the training data set. Since the current study is supervised, investigating a fusion of supervised and semiunsupervised methods can improve and strengthen the current model's performance while still taking advantage of the learned features from the labeled data set.
Future objectives also include further experimental work to collect sufficient ion (Z > 1) data sets. This will allow verifying our assumption that the current model, which is based on image recognition, identifies clusters from ions as protons as clusters from ions are visually more similar to protons than electrons or photons (Bouchami et al 2011a, Granja et al 2018a. Developing a CCN model to recognize heavy secondaries, such as carbon and oxygen, will allow assessing their contribution to the mixed radiation field (Hartmann et al 2017). In the case of light ions, such as deuterons or alphas, the separation from protons might be challenging due to the possible similarity of tracks produced by low-energy protons and high-energy light ions. Although this would provide more detailed knowledge of the composition of the mixed field, we believe the separation of light ions from protons would have a minor impact on the estimation of biological effects as long as microscopic ionization patterns and LET of protons and light ions are similar.
Of note, it is to be discussed how the LET measurement method presented in this manuscript relates to other non-commercial approaches based on microdosimetry. Microdosimetry Rosenzweig 1955, Rossi 1959) considers stochastic analogs to dose and LET, namely specific energy, lineal energy, and their spectra in microscopic structures such as cell nuclei (ICRU 1983). It plays a key role in linking measurable physical quantities describing the radiation and their biological effects (Rossi 1979). Tissue equivalent proportional counters are the most common tools for microdosimetric measurements and are considered the reference (Rossi and Rosenzweig 1955). To address the challenge of experimental characterization of the mixed radiation field, also passive detector technology, e.g. optically stimulated luminescent detectors Mckeever 2011, Parisi et al 2022b), as well as active detectors based on diamond (Verona et al 2018(Verona et al , 2020a  . However, none of these methods have been commercialized and implemented in the clinical practice of proton radiotherapy yet. The presented approach, despite its limitation related to the silicon sensor thickness, aims, in turn, at developing a method for fast systematic validation of MC simulation-based treatment planning systems in terms of LET of individual protons as well as average LET distributions in individual voxels defined by the volume of the Timepix sensor.
Experimental characterization of LET in silicon for treatment planning is an essential advancement towards moving beyond absorbed dose measurements in the clinical routine. We demonstrated here that LET spectra and LET D values computed and measured in silicon generally agree, thus, it can be assumed that simulations of LET in water will also be correct.
The presented approach has the potential to be applied for the commissioning of MC-based TPS in terms of LET spectra and LET D , which is a necessary step toward the clinical implementation of proton treatment planning that includes LET as a radiation quality factor. Currently, there is a lack of clarity on whether LET spectra and LET D measurements should be included in patient-specific QA procedures. The results of this work suggest that a wide range of proton LET values may occur in the mixed radiation fields produced during proton therapy. While the impact of LET differences in the low-dose penumbra for individual pencil beams might be less relevant, in the case of treatment plans consisting of multiple pencil beams the effect can accumulate. Thus, future work also includes verification of LET spectra and LET D for clinical treatment plans with both MC simulations and measurements with Timepix, to investigate the high-LET component in the mixed radiation fields.

Conclusions
In this work, we proposed a method for experimental verification of LET spectra and average LET distributions for mixed radiation fields produced by the therapeutic proton beam. For this purpose, a single-particle sensitive semiconductor pixel detector Timepix operated in a water phantom and CNN for particle track recognition were applied. For recognized protons, LET was computed, and the resulting LET spectra and LET D values were compared with MC simulation results.
Currently, treatment planning in proton therapy disregards LET of protons, and only physical dose measurements are performed for QA and QC purposes. Implementation of more advanced LET-based treatment planning is obstructed by the lack of commercially available LET spectra or LET D measurement methods. Timepix detector can be a potential solution offering both spectral information and single particle tracking, contrary to conventional detectors that integrate deposited energy.

Acknowledgments
Work in CCB IFJ PAN was supported by the National Center for Research and Development (NCBiR) within grant No. LIDER/43/0222/L-12/20/NCBR/2021. Work in Advacam was performed in the frame of Contract No. 4000130480/20/NL/GLC/hh from the European Space Agency. The calculations presented in this work were partially performed using the Ziemowit computer cluster in the Laboratory of Bioinformatics and Computational Biology at the Biotechnology Centre, Silesian University of Technology, created in the EU Innovative Economy Programme POIG.02.01.00-00-166/08 and expanded in the POIG.02.03.01-00-040/13 project. The work performed by CO was partially supported by the 18HLT04 UHDpulse funded by the EMPIR programme.