Open AccessArticle

The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres

Edoardo d’Andrea

^1,2

Andrea Politano

^1,2,

Bartolomeo Cassano

³,

Nico Lanconelli

⁴,

Marta Cremonesi

⁵,

Vincenzo Patera

^6,7 and

Massimiliano Pacilio

^2,*

Postgraduate School of Medical Physics, Azienda Ospedaliero-Universitaria Policlinico Umberto I, 00161 Rome, Italy

Department of Medical Physics, Azienda Ospedaliero-Universitaria Policlinico Umberto I, 00161 Rome, Italy

Department of Research, Diagnosis and Innovative Technologies, Medical Physics Unit, IRCCS Regina Elena National Cancer Institute, 00144 Rome, Italy

⁴

Department of Physics and Astronomy, Università di Bologna, 40126 Bologna, Italy

⁵

Medical Physics Unit, European Institute of Oncology, IRCCS, 20141 Milan, Italy

⁶

Department of Basic and Applied Science for Engineering, University of Roma “Sapienza”, 00161 Rome, Italy

⁷

Sezione di Roma I, Istituto Nazionale di Fisica Nucleare, 00161 Rome, Italy

Author to whom correspondence should be addressed.

Appl. Sci. 2025, 15(2), 958; https://doi.org/10.3390/app15020958

Submission received: 18 November 2024 / Revised: 14 January 2025 / Accepted: 16 January 2025 / Published: 19 January 2025

(This article belongs to the Special Issue Internal and/or External Radiotherapy and Dosimetric Evaluation: From Cellular Level to the Patient)

Download

Browse Figures

Versions Notes

Abstract

This study compares dosimetric approaches for lung dosimetry in ¹⁶⁶ radioembolization (Ho-TARE) with direct Monte Carlo (MC) simulations on a voxelized anthropomorphic phantom derived from a real patient’s CT scan, preserving the patient’s lung density distribution. Lung dosimetry was assessed for five lung shunt (LS) scenarios with conventional methods: the mono-compartmental organ-level approach (MIRD), voxel S-value convolution for soft tissue (kST, ICRU soft tissue with 1.04 g/cm³) and lung tissue (kLT, ICRU lung tissue with 0.296 g/cm³), local density rescaling (kST_L and kLT_L, respectively, for soft tissue and lung tissue), or global rescaling for a lung mean density of 0.221 g/cm³ (kLT₂₂₁). Significant underestimations in the mean absorbed dose (AD) were observed, with relative differences with respect to the reference (MC) of −64% for MIRD, −93% for kST, −56% for kST_L, −76% for kLT, −68% for kLT₂₂₁, and −60% for kLT_L. Given the high heterogeneity of lung tissue, standard dosimetric approaches cannot accurately estimate the AD. Additionally, MC results for ¹⁶⁶Ho showed notable spatial absorbed dose inhomogeneity, highlighting the need for tailored lung dosimetry in Ho-TARE accounting for the patient-specific lung density distribution. MC-based dosimetry thus proves to be essential for safe and effective radioembolization treatment planning in the presence of LS.

Keywords:

lung shunt; internal dosimetry; Monte Carlo; radioembolization; SIRT; TARE; lung dosimetry; molecular radiotherapy; liver cancer

1. Introduction

The use of Transarterial Radioembolization (TARE) for cancerous liver lesions is increasingly taking on the role of a curative treatment, moving beyond its traditional use as a palliative option, which has been the norm until recent times. For instance, TARE is gaining interest in the clinical landscape of liver diseases as a neoadjuvant approach to hepatic lobectomy [1] or as a bridge to transplantation [2]. In these new clinical scenarios, personalized dosimetry becomes essential to assess the tolerability of the treatment, especially in the presence of a lung shunt (LS), as the absorbed dose (AD) is the only reliable risk indicator for the preliminary assessment of the potential occurrence of tissue reactions.

The complications occurring after radioembolization treatment are classified into two main categories [3,4]: post-radioembolization syndrome (PRS) and complications related to non-target deposition of the microspheres and cross-irradiation to the structures surrounding the target regions. Clinical evidence of PRS includes symptoms like nausea, fatigue, vomiting, anorexia, fever and abdominal pain, caused by internal radiation and microembolization [3].

The off-target distribution of microspheres to non-tumor tissues may result in complications like radioembolization-induced liver disease (REILD), pulmonary complications (i.e., radiation pneumonitis), and gastrointestinal (GI) tract complications.

REILD [4,5,6] arises following normal hepatic parenchyma radiation exposure as a consequence of both the deposition of microspheres in normal parenchyma and the irradiation of surrounding structures. REILD leads to aberrations and hepatic dysfunction with the appearance of jaundice and ascites.

Lung complications arise from the non-target deposition of microspheres in patients with high LS fractions. Radiation pneumonitis symptoms are dry cough, fever, and ventilatory dysfunction manifesting within 6 months after treatment, with an incidence rate of <

1 %

[3,5,6].

The distribution of microspheres in the GI tract can occur in the case of variant anatomy and circulation, changes in flow dynamics during microsphere infusion, and microsphere reflux. Even low doses of radiation can compromise gastric acid production with consequent mucosal thinning, edema, and chronic inflammation, leading to gastroenteritis or ulceration [3,5,6].

Other side effects for the off-target deposition of microspheres are pancreatitis and cholecystitis, which could also arise from unexpected cutaneous and gallbladder embolization.

In this scenario, internal dosimetry evaluation plays a crucial role in preventing complications and effectively managing the clinical outcomes of radioembolization treatment.

The traditional eligibility criteria based on generic thresholds of pulmonary activity [7,8] are highly conservative, not patient-specific, and are no longer adequate, as they do not provide an appropriate level of personalization, which may lead to the exclusion of patients more easily, thereby limiting the therapeutic potential of TARE. Recently, a dosimetric study on ⁹⁰Y-labeled microspheres was conducted to evaluate the accuracy of conventional and simplified methods for personalized lung dosimetry, in comparison with the Monte Carlo (MC) direct radiation transport simulation technique [9]. The results demonstrated the inadequacy of simplified methods due to the high heterogeneity of lung tissue (LT) and the range of values skewed toward low tissue densities, highlighting the need for MC simulation of radiation transport.

In TARE with ⁹⁰Y-labeled microspheres, pre-treatment dosimetry is performed using a gamma-emitting radiotracer, typically ^99mTc-labeled macroaggregated albumin (^99mTc-MAA), and evaluated via single-photon emission computed tomography (SPECT) imaging [7,8]. This choice historically stems from the lack of suitable imaging channels for ⁹⁰Y, necessitating the use of a different radiotracer. Even if ⁹⁰Y possesses a small channel of internal pair creation that allows positron emission tomography (PET) imaging [10], the use of ⁹⁰Y microspheres as a pre-treatment radiotracer, given the low branching ratio (about

3 \times 10^{- 5}

pairs/decay [10,11]), would require high activities to perform imaging acquisition with sufficient resolution to be used for dosimetric evaluation. Nevertheless, the use of ^99mTc-MAA to simulate the biodistribution of therapeutic microspheres is intrinsically biased, being two different physical agents and not accurately reproducing the actual biodistribution of the treatment microspheres [12,13,14].

In recent years, a new therapeutic device based on poly-L-lactic acid (PLLA) microspheres labeled with ¹⁶⁶Ho (Ho-MS) has been approved for clinical use (Quiremspheres^®, Terumo Therapeutic Interventional Oncology, Deventer, Netherlands) [15]. ¹⁶⁶Ho is a

β^{-}

emitter with a half-life of approximately 26.8 h, exhibiting multiple

β

-decay pathways. The two most probable decay channels have average energies of 0.651 MeV and 0.693 MeV, with probabilities of 49.9% and 48.8%, respectively, decaying either to the first excited state of 166Er (followed by a characteristic gamma emission of approximately 81 keV) or directly to its ground state [16]. The primary advantage of Ho-MS lies in its multiple imaging-friendly properties: the intrinsic gamma emission allows for SPECT imaging, enabling pre-treatment biodistribution evaluation with a small amount of Ho-MS as a scout dose, which increases the predictive accuracy of pre-treatment dosimetry and may improve treatment precision and efficacy [17]. Additionally, the paramagnetic nature of Ho-MS enables magnetic resonance imaging (MRI), significantly enhancing the spatial resolution of the biodistribution, and potentially improving the dosimetric assessments [18]. Another promising approach is dual-isotope imaging, where both Ho-MS and ^99mTc-stannous phytate are administered simultaneously, allowing for the visualization of both healthy liver tissue and lesions. However, this approach presents technological and timing challenges [19].

According to current European guidelines [8], lung dosimetry in TARE using Ho-MS is addressed through a mono-compartmental model based on the well-established MIRD schema, which evaluates the LS fraction using planar imaging of ^99mTc-MAA and calculates the absorbed dose assuming a standard lung mass of 1 kg. The upper limit is set at a 30 Gy mean absorbed dose (

\bar{A D}

) in the lungs, associated with a 20% LS, based on this assumed lung mass. In the case of TARE with Ho-MS, investigating the accuracy of traditional, simplified approaches for lung dosimetry is also advisable, especially due to differences in the decay spectrum compared to ⁹⁰Y, as ¹⁶⁶Ho exhibits multiple low-energy

β^{-}

decay channels, along with more numerous and intense electron-capture (EC) decay emissions.

This study aims to compare the most common dosimetric approaches for calculating the

\bar{A D}

in the lungs in ¹⁶⁶Ho TARE, using MC simulation as a reference. The MC simulation models direct radiation transport on a CT-based anthropomorphic voxelized phantom, with realistic lung density distributions across various LS scenarios. This work serves as Part B of a previously published paper focused on ⁹⁰Y TARE [9]. In the final section, comparisons are reported between the MC results for both radionuclides, along with an analysis of the equivalent dose at 2 Gy per fraction (

E Q D_{2}

) [20] distributions for each.

2. Materials and Methods

2.1. Anthropomorphic Voxelized Phantom and Common Dosimetric Methods

The same phantom described and used in a previously published paper [9] was employed in the present study, following the same rationale. Here is a brief summary of its characteristics: the CT scan (voxel size of 1.37 × 1.37 × 3.27 mm³) of a healthy, normotype patient (weighing 72 kg and 166 cm tall) was segmented with ITK-Snap software (version 3.8.0) [21] into six tissue classes (air, lungs, adipose tissue, soft tissue, liver, and bone) and resampled to cubic voxels with a 2.21 mm side (converting a volume of about 6.1 mm³ to 10.8 mm³) using 3D-Slicer software (version 5.2.2)’s [22] built-in module “Resample Scalar Volume” with the “Nearest Neighbor” algorithm [23].

Three spherical lesions with a density higher than the liver parenchyma (1.2 g/cm³ and 1.04 g/cm³, respectively [24,25]) were introduced with diameters of 10, 20, and 30 mm. For lung tissue, the actual spatial density distribution of the patient was assumed, based on the HU-density calibration of the CT scanner. Table 1 summarizes the main density features of the lungs.

Five activity distributions were considered, with specific LS conditions (0%, 10%, 20%, 30%, 40%) and a fixed activity concentration ratio between lesions and healthy liver of 5:1, assuming a uniform distribution of microspheres within each compartment. Numerical values representing the activity concentration in each voxel were determined using the “relative calibration” method [12].

The use of an anthropomorphic phantom created from a CT scan allows us to define a realistic and controlled framework to test several approaches of internal dosimetry.

Several common dosimetric methods were tested in this study (Table 2), namely the mono-compartmental organ-level model of the MIRD Committee (hereinafter referred to as MIRD), convolution with voxel S-values (VSVs) using the soft tissue (ST) (kST) and the lung tissue (LT) kernel (kLT), both computed with dedicated MC simulations (further details in the following section).

The MIRD model [7,8,26,27] is a widely used approach for dosimetry at the organ level with beta-emitter radiopharmaceuticals. It allows the computation of the mean absorbed dose within a specific 3D ROI, based on the activity in that ROI and the ROI’s mass, under the assumption of a uniform radionuclide distribution in a homogeneous medium with complete energy deposition within the compartment. While the MIRD model is routinely employed for liver dosimetry, it has also been proposed for lung dosimetry with additional approximations, such as assuming a lung mass of 1 kg and neglecting the contribution from the hepatic dome to the homolateral lung. In this study, the MIRD approach was applied using the specific mass (727 g) of the reference phantom.

The convolution with the VSV [28,29,30] employs a 3D kernel of values representing the absorbed dose per decay in a target voxel with its centroid at position

(i, j, k)

, due to the source voxel positioned at

(0, 0, 0)

. The VSV kernel is computed using dedicated Monte Carlo (MC) simulations of electron and photon transport in a lattice of uniform biological tissue (e.g., soft tissue or lung tissue) with well-defined characteristics. This kernel is then convolved with the time-integrated activity distribution map to derive the absorbed dose map. In this study, VSV kernels for both soft tissue and lung tissue were considered for ¹⁶⁶

Ho

. The absorbed dose was subsequently corrected for local density variations relative to the uniform tissue density assumed in the VSV calculations, with rescaling based on the CT-derived voxel density.

Local energy deposition (LED) was not considered suitable for dosimetry with Ho-MS due to the presence of

γ

emissions, which result in unavoidable radiation transport outside the source voxels. First-order (local) density corrections to the convolution results from both tissues (kST_L and kLT_L, respectively) were applied according to Dieudonné et al. [31]. A global density correction for the convolution results obtained from the lung kernel was also performed, rescaling for the phantom-specific mean lung density (kLT₂₂₁).

Table 2. Summary of the “common” dosimetric methods tested in this study. Each method is reported together with its corresponding abbreviation (Tag) used throughout the paper, the general mathematical formulation, and the specific strategy of personalization of the calculations (phantom-specific feature).

Method	Tag	General Formulation	Phantom-Specific Feature
Mono-compartmental	MIRD	$15.87 \times \frac{A_{0} \cdot L S}{M_{l u n g s}}$ ^(†)	Organ Mass
Convolution kernel ST	kST	$\tilde{A} \otimes V S V_{S T} (‡)$	Activity Distribution
Convolution kernel ST + Local Rescale	${kST}_{L}$	$\frac{\tilde{A} \otimes V S V_{S T}}{ρ (x, y, z)} \times ρ_{S T} (‡)$	Activity and Density Distributions
Convolution kernel LT	kLT	$\tilde{A} \otimes V S V_{L T} (*)$	Activity Distribution
Convolution kernel LT + Global Rescale	${kLT}_{221}$	$\frac{\tilde{A} \otimes V S V_{L T}}{{\bar{ρ}}_{l u n g s}} \times ρ_{L T} (*)$	Activity Distribution & Mean Organ Density
Convolution kernel LT + Local Rescale	${kLT}_{L}$	$\frac{\tilde{A} \otimes V S V_{L T}}{ρ (x, y, z)} \times ρ_{L T} (*)$	Activity and Density Distributions

^(†) The numerical factor [27] accounts for the physical decay of the specific radionuclide and the conversion factors needed to compute the absorbed dose to the whole organ in Gy. It requires the total administered activity in GBq (named

A_{0}

), the lung shunt fraction percentage (LS), and the mass of the lungs in kg (

M_{l u n g s}

). ^(‡)

\tilde{A}

= time-integrated activity map; ⊗ = convolution product;

V S V_{S T}

= soft tissue S-voxel kernel. ^(*)

\tilde{A}

= time-integrated activity map; ⊗ = convolution product;

V S V_{L T}

= lung tissue S-voxel kernel.

2.2. Monte Carlo Simulations

All MC simulations of radiation transport were performed using the GATE/Geant4 code, version 8.1 [32], on a desktop computer equipped with an Intel Core i7 Extreme Edition 8-Core Processor and 16 GB of Random Access Memory (RAM). The simulations were set up as described in the previously published work [9]. A brief summary of the simulation parameters is reported in Table 3.

The simulation input files included activity and density map distributions in Metaimage format, a file with tissue composition information, and “macro” files with all the simulation parameters. The density description of the lungs was obtained from the calibration curve of the CT image, whereas the ICRU lung tissue composition [24,25] was assigned to all the lung voxels based on user-defined HU ranges to identify the tissue class.

Dosimetric calculations using MC simulations were performed for all the activity maps defined for the reference phantom (see Section 2.1).

Additionally, a dedicated simulation was conducted to compute the VSV kernels for ST and LT (available as Supporting Material of this paper in “.CSV” format), assuming ICRU soft tissue (density of 1.04 g/cm³) or ICRU lung tissue (density of 0.296 g/cm³) as medium [24,25,33], each uniformly distributed over a cube with a 600 mm side and cubic voxels of 2.21 mm. The ¹⁶⁶Ho VSV kernels were validated by comparing the

\bar{A D}

to the liver and lungs obtained via convolution with the values from a direct MC simulation. In this validation, the reference phantom had activity only in the lungs (LS = 10%), or in the liver (without lesions), with uniform liver and lung densities of 1.04 g/cm³ and 0.296 g/cm³, respectively.

2.3. Evaluation of Voxel-Based Dosimetric Approaches

Each voxel-based dosimetric approach was evaluated in terms of

\bar{A D}

in the lung region, as shown in Equation (1).

\bar{A D} = \frac{1}{N_{Lungs}} \sum_{x, y, z} A D (x, y, z)

(1)

where the sum is taken over all voxels within the lung region of interest (ROI),

N_{Lungs}

is the total number of lung voxels, and

A D (x, y, z)

is the absorbed dose in the voxel with its centroid at the position

(x, y, z)

The relative difference (RD) with respect to the reference MC simulation was evaluated as defined in Equation (2). Correlation plots compare the reference MC

\bar{A D}

on the x-axis with respect to the tested method’s

\bar{A D}

on the y-axis for each LS condition.

R D (%) = \frac{Tested - MC}{MC} \times 100

(2)

As regards the dosimetric comparison between ⁹⁰Y and ¹⁶⁶Ho, MC simulations of both radionuclides have been evaluated as AD rate (ADr) maps for LS = 10% only. The ADr maps were compared using cumulative dose–volume histograms (DVH) and homogeneity index (Equation (3)) as defined by Semerenko et al. [34].

H I = \frac{D_{5 %}}{D_{95 %}}

(3)

In targeted radiotherapy (TRT), radiation is delivered continuously with a decreasing dose rate, causing both lethal and sub-lethal damage, the latter competing with the tissue repair mechanisms. Since the treatment extends over hours or days, the radiation-induced damage and the repair of sub-lethal damage occur simultaneously. The biological effective dose (BED) in the Dale formulation [35,36] quantifies the biological effects of the radiation delivery patterns taking into account the repair of sub-lethal damage as a competitive process with the radiation damage. The Dale BED formulation depends on two key parameters: the tissue’s

α / β

ratio, which describes tissue sensitivity to radiation, and the sub-lethal damage repair rate (

μ

), related to the rate at which single-strand DNA breaks are repaired. Repair half-times (

T_{μ}

) of typical organ tissues range from 0.5 to 3 h [36,37], and the repair rate is given by

μ = log (2) / T_{μ}

. The BED for continuous radiation delivery with an exponentially decreasing dose rate (Equation (4)) depends on the initial dose rate (

R_{0}

), the radionuclide’s decay constant (

λ = log (2) / T_{1 / 2}

, where

T_{1 / 2}

is the physical half-life), the repair rate (

μ

), and the tissue’s

α / β

ratio [36].

B E D_{T R T} = \frac{R_{0}}{λ} \cdot (1 + \frac{R_{0}}{(μ + λ) \cdot (α / β)})

(4)

Assuming an administered activity of ⁹⁰Y and ¹⁶⁶Ho so that the

\bar{A D}

in the whole liver region was 60 Gy [38], resulting from the mono-compartmental model for the liver region (

A_{m a x} (GBq) = \frac{60 Gy \cdot M_{l i v e r}}{D F}

, where

M_{l i v e r}

is the mass of the whole liver in kg and DF is a dosimetric factor that numerically is 15.87

Gy \cdot kg / GBq

for ¹⁶⁶Ho [27] and 49.87

Gy \cdot kg / GBq

for ⁹⁰Y [7]), 1.8 GBq and 5.7 GBq were obtained, respectively. The BED maps were computed by converting the MC voxel AD values using Equation (4).

T_{μ}

of 0.5 and 1 h for LT were considered [37], along with an

α / β

ratio of 3 Gy [39] and the radionuclide-specific physical half-life (

T_{1 / 2}

From the BED maps, the corresponding

E Q D_{2}

maps and the DVHs were calculated using the linear quadratic model (Equation (5)) [40], considering the values of

T_{μ}

for LT.

E Q D_{2} = \frac{B E D}{1 + \frac{2}{α / β}}

(5)

Additionally, the volumetric lung absorbed dose constraints from external beam radiotherapy (EBRT) [41], listed in Table 4, have been evaluated on the

E Q D_{2}

DVHs for the LS = 10% case.

Further analysis of the probability of incidence of radiation pneumonitis (RP) for the LS = 10% case was carried out using the normal tissue probability control (NTCP) curves extracted from the QUANTEC publication [42]. In particular, a fitted parametrized NTCP curve with a logistic model (Equation (6)) relative to the

\bar{A D}

was used to compute the incidence of RP for the specific

E Q D_{2}

of the simulated cases, using the fit parameters reported by the same publication (

b_{0}

= −3.87 and

b_{1}

= 0.126

{Gy}^{- 1}

, corresponding to an AD correlated to a complication probability of 50% of 30.75 Gy).

Probability of RP = \frac{e^{b_{0} + b_{1} \cdot \bar{A D}}}{1 + e^{b_{0} + b_{1} \cdot \bar{A D}}}

(6)

3. Results

3.1. VSV Kernel for Soft and Lung Tissue for ¹⁶⁶Ho

The LT and ST VSV kernels used for the convolution method are reported in Figure 1 as the AD to the target voxel per unit decay plotted against the source–target distance.

Both kernels were calculated considering the complete decay spectrum of ¹⁶⁶Ho, and considering an extended lattice to properly account for the photon contribution at high source–target distances. As previously reported [9], it is essential to crop the liver region activity when performing convolution calculations with the LT kernel, to prevent abnormally high liver-to-lung cross-irradiation caused by excessive radiation transport from the liver.

The two kernels were validated using dedicated MC simulations on the corresponding uniform tissue with uniform activity distributions (see also Section 2.2). The internal validation results for ¹⁶⁶Ho showed an RD of 1.6% for LT and −0.5% for ST (Table 5), well within the corresponding uncertainty.

The cross-irradiation contribution to the lungs from the liver region, considering liver-only activity for the LT kernel, was 0.12 Gy/GBq, compared to 0.04 Gy/GBq as obtained from the MC simulation under the same conditions.

3.2. ¹⁶⁶Ho Monte Carlo Simulations vs. “Classical” Dosimetric Approaches

The dosimetric results for ¹⁶⁶Ho obtained from the MC simulation versus the tested approaches (Table 2) are shown in the correlation plot in Figure 2.

The relative uncertainty of the MC absorbed dose values at a 95% confidence level was analogous to the previously analyzed paper on ⁹⁰Y [9], with a relative uncertainty of the mean absorbed dose of less than 5% for all the simulated LS cases.

All tested methods severely underestimate the MC results, with an RD of approximately −64% for MIRD, −93% for kST, −56% for

k S T_{L}

, −76% for kLT, −68% for

k L T_{221}

, and −60% for

k L T_{L}

. The RD results are reported in Table 6, while the

\bar{A D}

of each approach is reported in Table 7.

For

{kST}_{L}

and

{kLT}_{L}

, which showed a lower RD with respect to the reference MC simulation, an example slice in the coronal view of the AD distributions and the DVHs were compared with the MC data (Figure 3) for LS = 10% only.

The homogeneity of the distribution varies for each method due to the different radiation transport computations performed by each approach, with HI values of 4.6 for

{kST}_{L}

, 5.6 for

{kLT}_{L}

, and 6.1 for MC.

3.3. ¹⁶⁶Ho Versus ⁹⁰Y Monte Carlo Simulations

To obtain a first comparison of the performance of the two radionuclides, only the MC simulations for each device have been considered for LS = 10%. To make the AD distributions of ⁹⁰Y and ¹⁶⁶Ho independent by the differences in the physical half-life of the two radionuclides, the ADr maps and the corresponding dose rate–volume histograms (DrVHs) (Figure 4) were computed.

The ADr distribution of ¹⁶⁶Ho appears significantly more inhomogeneous compared to that of ⁹⁰Y, a fact supported by the HI, which is 3.5 for ⁹⁰Y and 6.3 for ¹⁶⁶Ho. This difference is due to the physical characteristics of the decay spectra of the two radionuclides, primarily the low-energy electron component of the ¹⁶⁶Ho spectrum, which is characterized by less energetic

β^{-}

decay spectra and discrete electrons [16]. The stronger inhomogeneity of the ADr distribution of ¹⁶⁶Ho implies a stronger correlation of the AD distribution with the anatomical density heterogeneity of the lung tissue.

3.4. Impact on Clinical Decision-Making: From the $\bar{A D}$ to the EBRT Dosimetric Constraints for DVH

All the results hereinafter presented are not intended to be for general application, given the strong dependence of the computations on the specific activity and density distributions of the evaluated case.

Table 8 summarizes the results for the MC

\bar{A D}

, the activity as high as safely administrable (AHASA), and the maximum lung activity (MLA) to ensure that the

\bar{A D}

to the lungs does not exceed the 30 Gy threshold [7,8] for the specific simulated case.

The MLA was found to be nearly constant across all LS cases, at approximately 500 MBq, due to a negligible contribution of cross-irradiation from the hepatic region to the lungs.

The AD maps were converted to BED maps using Equation (4), then to

E Q D_{2}

maps, from which the DVHs were extracted. In this section, we attempt to simulate a clinical evaluation of the LS = 10% reference phantom simulation as if treated with the prescribed activity (see Section 2.3).

The

\bar{A D}

from the MC simulations was found to be 21.4 Gy for ⁹⁰Y at 1.8 GBq, and 35.3 Gy for ¹⁶⁶Ho at 5.7 GBq. This suggests that the treatment with ⁹⁰Y is safely administrable, while for ¹⁶⁶Ho, a reduction in administered activity would be necessary.

Moreover, the mean

E Q D_{2}

values to lungs were between 13.5 Gy and 14.3 Gy for ⁹⁰Y and between 25.7 Gy and 30.1 Gy for ¹⁶⁶Ho. The ¹⁶⁶Ho results significantly exceed the reference constraint (see Table 4), confirming the need to reduce the administered activity.

The DVHs of the

E Q D_{2}

maps are shown in Figure 5 along with the volumetric constraints from Table 4.

For the simulated administered activities, ⁹⁰Y fully meets two constraints (

V_{30}

and

V_{40}

), while

V_{20}

compliance depends on the tissue

T_{μ}

. Conversely, ¹⁶⁶Ho partially meets the

V_{40}

constraint, depending on the specific

T_{μ}

. Neither radionuclide meets the

V_{5}

and

V_{10}

requirements.

The probability of RP onset is graphically represented in Figure 6, along with the points for the LS = 10% case of ⁹⁰Y and ¹⁶⁶Ho for the maximum and minimum

T_{μ}

of LT. The probability of developing RP for the ⁹⁰Y case is between 10% and 11%, while for ¹⁶⁶Ho, it is between 35% and 48%.

No incidence of RP associated with the volumetric constraints was reported due to the variability in the NTCP curves of the different datasets shown in the QUANTEC publication [42].

4. Discussion

The MC simulation seems to represent the most valuable dosimetric tool in TRT of low-density strongly heterogeneous tissues such as lung tissue, especially for radionuclides with a low-energy decay spectrum.

The results of the MC simulations with GATE/Geant4 code have been validated in previously published works for ⁹⁰Y [9,43]. Moreover, the results for ¹⁶⁶Ho were validated for both LT and ST in homogeneous lungs (with LS = 10% and a density of 0.296 g/cm³) and liver (with LS = 0% and density of 1.04 g/cm³), respectively, proving the consistency of the MC simulations with respect to the VSV convolution (see Section 3.1). A more robust validation, but beyond the scope of the present work, to assess the reliability of the MC simulation of strongly heterogeneous low-density lungs, would involve a comparison of different MC simulation codes.

To ensure that resampling from the original CT resolution to cubic voxels did not introduce significant density changes leading to notable variations in AD, an additional MC simulation was performed exclusively on the heterogeneous lung tissue. This simulation maintained the original CT resolution

(1.37 \times 1.37 \times 3.27 {mm}^{3})

while assuming a uniform activity distribution and an LS of 10%. The results were compared with the same heterogeneous lungs resampled using the procedure outlined in Section 2.1. The comparison showed an RD below 1% with respect to the simulation using the original CT resolution.

The analysis of the S-values plot over the source–voxel distance (Figure 1) allows us to test the main radiation–matter interaction features of ¹⁶⁶Ho in ST and LT, primarily the particle range in Continuous Slowing Down Approximation (CSDA), corresponding to the abrupt slope change in the shown plots. The CSDA range for the most energetic end-point electrons of the

β^{-}

spectrum is approximately 9 mm in ST and 30 mm in LT, well correlated with the corresponding values for the slope change of about 9.7 mm and 29.8 mm, respectively.

As for ⁹⁰Y [9], the use of the lung kernel brings the necessity of an image crop of the activity map to extrapolate the lungs’ activity only to overcome an unwanted cross-irradiation contribution to the lungs’ AD from the activity in the liver region.

Comparison of the “classical” dosimetric approaches (MIRD, kST,

{kST}_{L}

, kLT,

{kST}_{221}

and

{kST}_{L}

) with the reference MC simulations using ¹⁶⁶Ho show a remarkable underestimation of the lungs’

\bar{A D}

. The largest discrepancies arise for kST, resulting in −93% RD, and an analogous effect is obtained for ⁹⁰Y [9], due to the higher density of the soft tissue used in the VSV calculations compared to the lung density. On the other hand, the best performance was obtained by the convolution with local density corrections for both ST and LT, resulting in approximately

- 56 %

for

{kST}_{L}

and

- 60 %

for

{kLT}_{L}

. Notably, the use of a specific kernel for lung tissue with an ICRU tissue description does not yield better results in terms of

\bar{A D}

, which results in an

\bar{A D}

about 12% higher with respect to the

{kST}_{L}

. This last result is remarkably higher than the analogous previously analyzed results on ⁹⁰Y [9], strengthening the hypothesis that, for highly heterogeneous low-density tissues, the description and characterization of the tissue’s density distribution is essential for accurate lung dosimetry in MRT, especially when the radionuclide involved has abundant low-energy decay emissions.

It is noteworthy that the LED method is not applicable to ¹⁶⁶Ho due to the gamma emission of ¹⁶⁶Ho that makes the assumption of local deposition unreasonable.

Examining the AD spatial distributions of the reference MC simulation in comparison to the two voxel-based dosimetric methods with stronger physical assumptions (

{kST}_{L}

and

{kLT}_{L}

) for LS = 10%, a pronounced heterogeneity can be easily observed, underscoring the fundamental importance of quantifying the correlation between the density distribution and the AD distribution. This is particularly crucial for radionuclides such as ¹⁶⁶Ho, which exhibit sufficiently abundant low-energy

β^{-}

decay electrons and internal conversion electrons, making the AD distribution more sensitive to morphological heterogeneities.

It is necessary to consider that all results obtained in this paper should not be generalized as a global trend for data obtained by each method; thus, all results are specific to the simulated phantom, and the assumption of the same RD in clinical practice should be discouraged. However, the general performance of each method is representative of its limitations, with the major differences arising from the physical assumptions of each approach.

Looking at the ADr spatial distribution differences in the MC simulations of ⁹⁰Y and ¹⁶⁶Ho for the LS = 10% case, it is clear that ¹⁶⁶Ho exhibits a more heterogeneous distribution than ⁹⁰Y, with a doubled HI. These remarkable spatial ADr differences, together with variations in the absolute ADr values, arise because the decay spectrum of ¹⁶⁶Ho includes a greater variety of low-energy emissions, leading to a highly inhomogeneous dose rate distribution that is strongly dependent on the anatomy of the individual patient. This makes the choice of methods that account for radiation transport particularly important, especially for settings like LT.

To our knowledge, there are no prior studies on ¹⁶⁶Ho lung dosimetry; thus, this study serves as a groundbreaking MC-based investigation into lung dosimetry for ¹⁶⁶Ho-labeled microspheres liver radioembolization, highlighting the need and specific challenges of lung dosimetry for this therapeutic device.

Studying the potential clinical impact of the results presented in this work, it is evident from Table 8, as previously demonstrated for ⁹⁰Y [9], that the primary determining factor for AD in the lungs is the activity present in the pulmonary region. For the examined case, to prevent the lungs from reaching the 30 Gy limit advised by international guidelines [7,8], the maximum activity that lung tissue can tolerate is approximately 500 MBq, corresponding to a maximum administrable activity ranging from 4.88 GBq to 1.26 GBq along with increasing LS.

The clinical impact becomes even more significant when attempting to simulate a treatment for both ¹⁶⁶Ho and ⁹⁰Y for the case of LS = 10%. Planning the treatment according to the guidelines of the Ho-MS manuals [15], it is necessary to administer an activity sufficient to achieve 60 Gy to the whole liver using the mono-compartmental MIRD model. Using the same criterion for both radionuclides, a treatment activity of 1.8 GBq for ⁹⁰Y and 5.7 GBq for ¹⁶⁶Ho was assumed. According to these prescriptions, the

\bar{A D}

to the lungs for ¹⁶⁶Ho exceeds the 30 Gy limit. Converting the AD map with these administrations into

E Q D_{2}

by applying the corresponding BED values, the corresponding

\bar{A D}

E Q D_{2}

again exceeds the constraints for EBRT listed in Table 4. Evaluating the volumetric distribution of

E Q D_{2}

with the DVH for each radionuclide, the five volumetric constraints in Table 4 are mostly met for ⁹⁰Y, while for ¹⁶⁶Ho, almost all of them are violated, except for

V_{40}

depending on the considered

T_{μ}

. However, in both cases, the constraints associated with low isodoses are violated for both devices. This is primarily due to the fact that the constraints considered pertain to partial irradiations of lung tissue in EBRT treatments, not suitably matching the simulated irradiation conditions, i.e., uniformly distributed activity distribution. This characteristic also limits the strength and validity of the evaluation of all the volumetric constraints considered.

An additional point is brought by evaluating the probability of radiation-induced pneumonitis based on known NTCP curves for partial lung irradiation from EBRT [42], indicating a probability between 10% and 11% for ⁹⁰Y and between 35% and 48% for ¹⁶⁶Ho according to the interpolation of the BED concept. However, again, these values should be considered with caution, as they refer to datasets based on partial lung irradiations, where the probability of RP is strongly dependent not only on the absorbed dose but also on the irradiated area.

The high dose values recorded in this study should be regarded as entirely case-specific and remarkably influenced by the morphological characteristics of the lung tissue considered, which in this specific case has a significantly low density and a small volume, representing a particular case where dosimetric verification is particularly crucial.

To emphasize the importance of the results already discussed, a patient treated with ¹⁶⁶Ho with a LS of approximately 11% was evaluated. This patient was characterized by a lung mass of 908 g and a mean density of 0.28 g/cm³, compared to the reference phantom lung mass of 727 g and a mean density of 0.22 g/cm³, and a liver lesion located superficially on the hepatic dome at the liver–lung interface. The SPECT imaging revealed a difference in uptake between the ipsilateral lung and the contralateral one, with the former showing slightly higher activity, primarily attributable to the combined effect of the imaging acquisition protocol and the superficial position of the liver lesion. The prescribed activity for the patient was 4.33 GBq, resulting in an average lung dose of 19.3 Gy from the MC simulation, which can be compared to the analogous value for the phantom with an LS of 10%, which is 26.6 Gy. This significant difference in

\bar{A D}

is mainly due to the mass difference with respect to the phantom’s lungs.

Evaluating the corresponding probability of RP occurrence according to the previously outlined model and within the discussed limits, it is found to be between 10% and 11% for the actual patient and between 18% and 23% for the reference phantom, depending on the considered

T_{μ}

Lung dosimetry remains a crucial issue in TARE in the presence of LS. The need to pay major attention to the dosimetric analysis of patients with considerable levels of LS, even when below the treatment exclusion levels emphasized by international guidelines, becomes highly significant given the new clinical role that this therapeutic option is assuming, transforming from a simple palliative tool to a curative approach as a bridge to liver transplantation [2] or for downstaging lesions for subsequent surgical intervention [44]. In this context, personalized treatment planning becomes a fundamental tool at the service of clinical practice, shifting the focus from a mere assessment of the amount of activity present in the pulmonary compartment to a specific dosimetric evaluation, both in terms of the mean dose absorbed by lung tissue and in terms of the spatial distribution of the absorbed dose. This evaluation, combined with the use of appropriate constraints on the corresponding DVHs, could help manage the specificities of patients in a targeted manner, potentially expanding the horizon of treatable patients with RE.

A useful approach to the decision-making process is certainly the definition of a probability of adverse effects’ occurrence at a certain level of exposure of healthy tissue through NTCP curves. Therefore, it is beneficial for radioembolization therapy to have models specifically tailored to this irradiation modality, along with dosimetric constraints on DVHs derived from future studies on patient cohorts.

Currently, the most reliable approach for lung dosimetry appears to be the direct simulation of radiation transport using MC code, allowing for a complete consideration of the specificity of patient characteristics from both morphological and functional perspectives.

We do not claim that the results presented in this research are a generalizable guideline; each treated case requires independent evaluation. Moreover, a real case is also affected by the specifics of emission imaging, a topic neither addressed nor considered in the presented results.

Future studies could analyze the impact of various clinical scenarios to better assess the benefits of direct MC simulations and their role in the clinical landscape of lung dosimetry in liver radioembolization.

5. Conclusions

Lung dosimetry tailored for ¹⁶⁶Ho radioembolization is a fundamental step in the treatment workflow for patients undergoing these procedures, especially in cases where involvement of the pulmonary compartment is significant enough to create a risk for the development of adverse effects such as radiation pneumonitis.

This study demonstrates that, due to the significant heterogeneity and low density of lung tissue, the radiation transport from radionuclide decay is substantial, making conventional dosimetric calculation methods inadequate for accurate dosimetric evaluation under these conditions.

Furthermore, this study has shown that the complexity of the ¹⁶⁶Ho spectrum demands a proper dosimetric evaluation and morphological framework to correctly manage possible correlations between the heterogeneity of the absorbed dose distribution in lung tissue and potential adverse effects.

The MC simulations comparing ¹⁶⁶Ho and ⁹⁰Y showed a significant difference in the spatial distribution of ADr between the two radionuclides, with greater ADr heterogeneity for ¹⁶⁶Ho due to its generally lower energy and complex decay spectrum.

In this context, internal dosimetry using direct radiation transport simulation with Monte Carlo code has become an applicable approach in the clinical scenario, even with common hardware tools, with computation times fully compatible with everyday practice. Moreover, the advent of Monte Carlo codes on GPUs is already making the use of this type of dosimetric calculation tool agile and straightforward, with the potential to further implement personalized planning, reaching a complexity comparable to that of EBRT treatments, thereby positively impacting the future treatment and management of patients undergoing these therapies.

Further understanding of the clinical role of direct MC simulations in lung dosimetry of liver radioembolization treatments could be improved by analyzing the impact of different clinical situations, evaluating patients with relatively high and various LS values.

In assessing the probability of adverse effects, it is essential to adopt a robust and reliable method, such as direct MC simulation, to develop models for evaluating NTCP curves, potentially paving the way for more advanced and consistent planning strategies of the TRT treatments.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/app15020958/s1, ¹⁶⁶Ho S-value kernel for the soft tissue ICRU with a density of 1.04 g/cm³ and lung tissue ICRU with a density of 0.296 g/cm³ over 2.21 mm voxel size. Data are available in a csv-formatted file as position indices

(i, j, k)

followed by the absorbed dose values expressed in mGy/MBq·s.

Author Contributions

Conceptualization, E.d. and M.P.; methodology, E.d., N.L. and M.P.; software, E.d.; validation, E.d., A.P. and M.P.; formal analysis, E.d., A.P., B.C. and M.P.; investigation, E.d., A.P., N.L. and M.P.; resources, E.d., M.C., V.P. and M.P.; data curation, E.d., A.P., B.C. and M.P.; writing—original draft preparation, E.d.; writing—review and editing, A.P., B.C., N.L., M.C., V.P. and M.P.; supervision, B.C., V.P. and M.P.; project administration, M.P. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The original contributions presented in this study are included in the article and in the Supplementary Material; further inquiries can be directed to the corresponding author.

Conflicts of Interest

The authors declare no conflicts of interest.

References

Ahmed, A.; Stauffer, J.A.; LeGout, J.D.; Burns, J.; Croome, K.; Paz-Fumagalli, R.; Frey, G.; Toskich, B. The use of neoadjuvant lobar radioembolization prior to major hepatic resection for malignancy results in a low rate of post hepatectomy liver failure. J. Gastrointest. Oncol. 2021, 12, 751. [Google Scholar] [CrossRef] [PubMed]
Tohme, S.; Sukato, D.; Chen, H.W.; Amesur, N.; Zajko, A.B.; Humar, A.; Geller, D.A.; Marsh, J.W.; Tsung, A. Yttrium-90 radioembolization as a bridge to liver transplantation: A single-institution experience. J. Vasc. Interv. Radiol. 2013, 24, 1632–1638. [Google Scholar] [CrossRef] [PubMed]
Riaz, A.; Lewandowski, R.J.; Kulik, L.M.; Mulcahy, M.F.; Sato, K.T.; Ryu, R.K.; Omary, R.A.; Salem, R. Complications following radioembolization with yttrium-90 microspheres: A comprehensive literature review. J. Vasc. Interv. Radiol. 2009, 20, 1121–1130. [Google Scholar] [CrossRef] [PubMed]
Riaz, A.; Awais, R.; Salem, R. Side effects of yttrium-90 radioembolization. Front. Oncol. 2014, 4, 198. [Google Scholar] [CrossRef] [PubMed]
Laidlaw, G.L.; Johnson, G.E. Recognizing and managing adverse events in Y-90 radioembolization. In Proceedings of the Seminars in Interventional Radiology; Thieme Medical Publishers, Inc.: Stuttgart, Germany, 2021; Volume 38, pp. 453–459. [Google Scholar]
Sangro, B.; Martínez-Urbistondo, D.; Bester, L.; Bilbao, J.I.; Coldwell, D.M.; Flamen, P.; Kennedy, A.; Ricke, J.; Sharma, R.A. Prevention and treatment of complications of selective internal radiation therapy: Expert guidance and systematic review. Hepatology 2017, 66, 969–982. [Google Scholar] [CrossRef] [PubMed]
Dezarn, W.A.; Cessna, J.T.; DeWerd, L.A.; Feng, W.; Gates, V.L.; Halama, J.; Kennedy, A.S.; Nag, S.; Sarfaraz, M.; Sehgal, V.; et al. Recommendations of the American Association of Physicists in Medicine on dosimetry, imaging, and quality assurance procedures for ⁹⁰Y microsphere brachytherapy in the treatment of hepatic malignancies. Med. Phys. 2011, 38, 4824–4845. [Google Scholar] [CrossRef] [PubMed]
Weber, M.; Lam, M.; Chiesa, C.; Konijnenberg, M.; Cremonesi, M.; Flamen, P.; Gnesin, S.; Bodei, L.; Kracmerova, T.; Luster, M.; et al. EANM procedure guideline for the treatment of liver cancer and liver metastases with intra-arterial radioactive compounds. Eur. J. Nucl. Med. Mol. Imaging 2022, 49, 1682–1699. [Google Scholar] [CrossRef] [PubMed]
d’Andrea, E.; Lanconelli, N.; Cremonesi, M.; Patera, V.; Pacilio, M. The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization with ⁹⁰Y Microspheres. Appl. Sci. 2024, 14, 7684. [Google Scholar] [CrossRef]
Kao, Y.H.; Steinberg, J.D.; Tay, Y.S.; Lim, G.K.; Yan, J.; Townsend, D.W.; Takano, A.; Burgmans, M.C.; Irani, F.G.; Teo, T.K.; et al. Post-radioembolization yttrium-90 PET/CT-part 1: Diagnostic reporting. EJNMMI Res. 2013, 3, 1–13. [Google Scholar] [CrossRef] [PubMed]
Dryák, P.; Šolc, J. Measurement of the branching ratio related to the internal pair production of Y-90. Appl. Radiat. Isot. 2020, 156, 108942. [Google Scholar] [CrossRef] [PubMed]
Chiesa, C.; Sjogreen-Gleisner, K.; Walrand, S.; Strigari, L.; Flux, G.; Gear, J.; Stokke, C.; Gabina, P.M.; Bernhardt, P.; Konijnenberg, M. EANM dosimetry committee series on standard operational procedures: A unified methodology for 99m Tc-MAA pre-and 90 Y peri-therapy dosimetry in liver radioembolization with 90 Y microspheres. EJNMMI Phys. 2021, 8, 1–44. [Google Scholar] [CrossRef] [PubMed]
Elschot, M.; Nijsen, J.F.; Lam, M.G.; Smits, M.L.; Prince, J.F.; Viergever, M.A.; van den Bosch, M.A.; Zonnenberg, B.A.; de Jong, H.W. 99m Tc-MAA overestimates the absorbed dose to the lungs in radioembolization: A quantitative evaluation in patients treated with 166 Ho-microspheres. Eur. J. Nucl. Med. Mol. Imaging 2014, 41, 1965–1975. [Google Scholar] [CrossRef] [PubMed]
Kokabi, N.; Webster, L.A.; Elsayed, M.; Switchenko, J.M.; Chen, B.; Brandon, D.; Galt, J.; Sethi, I.; Cristescu, M.; Kappadath, S.C.; et al. Accuracy and safety of scout dose resin yttrium-90 microspheres for radioembolization therapy treatment planning: A prospective single-arm clinical trial. J. Vasc. Interv. Radiol. 2022, 33, 1578–1587. [Google Scholar] [CrossRef] [PubMed]
Quirem: Quirem-SpheresTMMicrospheres—QSuite v2.1 Instructions for Use. Available online: https://www.quirem.com/wp-content/uploads/2022/05/LC-80095-01-Q-Suite-2.1-IFU-Multilingual-1.pdf (accessed on 10 October 2024).
NuDat 3: Nuclear Structure and Decay Data Provided by National Nuclear Data Center (NNDC) at Brookhaven National Laboratory. Available online: https://www.nndc.bnl.gov/nudat3/decaysearchdirect.jsp?nuc=166Ho&unc=NDS (accessed on 10 October 2024).
Smits, M.L.; Dassen, M.G.; Prince, J.F.; Braat, A.J.; Beijst, C.; Bruijnen, R.C.; de Jong, H.W.; Lam, M.G. The superior predictive value of 166 Ho-scout compared with 99m Tc-macroaggregated albumin prior to 166 Ho-microspheres radioembolization in patients with liver metastases. Eur. J. Nucl. Med. Mol. Imaging 2020, 47, 798–806. [Google Scholar] [CrossRef]
Smits, M.L.; Elschot, M.; van den Bosch, M.A.; van de Maat, G.H.; van het Schip, A.D.; Zonnenberg, B.A.; Seevinck, P.R.; Verkooijen, H.M.; Bakker, C.J.; de Jong, H.W.; et al. In vivo dosimetry based on SPECT and MR imaging of ¹⁶⁶Ho-microspheres for treatment of liver malignancies. J. Nucl. Med. 2013, 54, 2093–2100. [Google Scholar] [CrossRef]
van Rooij, R.; Braat, A.J.; de Jong, H.W.; Lam, M.G. Simultaneous 166 Ho/99m Tc dual-isotope SPECT with Monte Carlo-based downscatter correction for automatic liver dosimetry in radioembolization. EJNMMI Phys. 2020, 7, 1–12. [Google Scholar] [CrossRef]
Sgouros, G.; Bolch, W.E.; Chiti, A.; Dewaraja, Y.K.; Emfietzoglou, D.; Hobbs, R.F.; Konijnenberg, M.; Sjögreen-Gleisner, K.; Strigari, L.; Yen, T.C.; et al. ICRU REPORT 96, dosimetry-guided radiopharmaceutical therapy. J. ICRU 2021, 21, 1–212. [Google Scholar] [CrossRef]
Yushkevich, P.A.; Piven, J.; Hazlett, H.C.; Smith, R.G.; Ho, S.; Gee, J.C.; Gerig, G. User-guided 3D active contour segmentation of anatomical structures: Significantly improved efficiency and reliability. Neuroimage 2006, 31, 1116–1128. [Google Scholar] [CrossRef] [PubMed]
Fedorov, A.; Beichel, R.; Kalpathy-Cramer, J.; Finet, J.; Fillion-Robin, J.C.; Pujol, S.; Bauer, C.; Jennings, D.; Fennessy, F.; Sonka, M.; et al. 3D Slicer as an image computing platform for the Quantitative Imaging Network. Magn. Reson. Imaging 2012, 30, 1323–1341. [Google Scholar] [CrossRef]
Meijering, E.H.; Niessen, W.J.; Pluim, J.P.; Viergever, M.A. Quantitative comparison of sinc-approximating kernels for medical image interpolation. In Proceedings of the Medical Image Computing and Computer-Assisted Intervention–MICCAI’99: Second International Conference, Cambridge, UK, 19–22 September 1999; pp. 210–217. [Google Scholar]
ICRU. Tissue substitutes in radiation dosimetry and measurement. In ICRU Report No 44; International Commission on Radiation Units and Measurements: Bethesda, MD, USA, 1989. [Google Scholar]
ICRU. Photon, electron, proton and neutron interaction data for body tissues. In ICRU Report No 46; International Commission on Radiation Units and Measurements: Bethesda, MD, USA, 1992. [Google Scholar]
Bolch, W.E.; Eckerman, K.F.; Sgouros, G.; Thomas, S.R. MIRD pamphlet no. 21: A generalized schema for radiopharmaceutical dosimetry—standardization of nomenclature. J. Nucl. Med. 2009, 50, 477–484. [Google Scholar] [CrossRef] [PubMed]
Smits, M.L.; Nijsen, J.F.; van den Bosch, M.A.; Lam, M.G.; Vente, M.A.; Huijbregts, J.E.; van het Schip, A.D.; Elschot, M.; Bult, W.; de Jong, H.W.; et al. Holmium-166 radioembolization for the treatment of patients with liver metastases: Design of the phase I HEPAR trial. J. Exp. Clin. Cancer Res. 2010, 29, 1–11. [Google Scholar] [CrossRef] [PubMed]
Bolch, W.E.; Bouchet, L.G.; Robertson, J.S.; Wessels, B.W.; Siegel, J.A.; Howell, R.W.; Erdi, A.K.; Aydogan, B.; Costes, S.; Watson, E.E.; et al. MIRD pamphlet no. 17: The dosimetry of nonuniform activity distributions—Radionuclide S values at the voxel level. J. Nucl. Med. 1999, 40, 11S–36S. [Google Scholar]
Lanconelli, N.; Pacilio, M.; Meo, S.L.; Botta, F.; Di Dia, A.; Aroche, L.T.; Pérez, M.C.; Cremonesi, M. A free database of radionuclide voxel S values for the dosimetry of nonuniform activity distributions. Phys. Med. Biol. 2012, 57, 517. [Google Scholar] [CrossRef] [PubMed]
Pacilio, M.; Amato, E.; Lanconelli, N.; Basile, C.; Torres, L.A.; Botta, F.; Ferrari, M.; Diaz, N.C.; Perez, M.C.; Fernández, M.; et al. Differences in 3D dose distributions due to calculation method of voxel S-values and the influence of image blurring in SPECT. Phys. Med. Biol. 2015, 60, 1945. [Google Scholar] [CrossRef] [PubMed]
Dieudonné, A.; Hobbs, R.F.; Lebtahi, R.; Maurel, F.; Baechler, S.; Wahl, R.L.; Boubaker, A.; Le Guludec, D.; Sgouros, G.; Gardin, I. Study of the impact of tissue density heterogeneities on 3-dimensional abdominal dosimetry: Comparison between dose kernel convolution and direct Monte Carlo methods. J. Nucl. Med. 2013, 54, 236–243. [Google Scholar] [CrossRef]
Jan, S.; Santin, G.; Strul, D.; Staelens, S.; Assié, K.; Autret, D.; Avner, S.; Barbier, R.; Bardiès, M.; Bloomfield, P.; et al. GATE: A simulation toolkit for PET and SPECT. Phys. Med. Biol. 2004, 49, 4543. [Google Scholar] [CrossRef]
Van Dyk, J.; Keane, T.; Rider, W. Lung density as measured by computerized tomography: Implications for radiotherapy. Int. J. Radiat. Oncol. Biol. Phys. 1982, 8, 1363–1372. [Google Scholar] [CrossRef] [PubMed]
Semenenko, V.A.; Reitz, B.; Day, E.; Qi, X.S.; Miften, M.; Li, X.A. Evaluation of a commercial biologically based IMRT treatment planning system. Med. Phys. 2008, 35, 5851–5860. [Google Scholar] [CrossRef]
Dale, R. Dose-rate effects in targeted radiotherapy. Phys. Med. Biol. 1996, 41, 1871. [Google Scholar] [CrossRef]
Dale, R.; Carabe-Fernandez, A. The radiobiology of conventional radiotherapy and its application to radionuclide therapy. Cancer Biother. Radiopharm. 2005, 20, 47–51. [Google Scholar] [CrossRef] [PubMed]
Pacilio, M.; Betti, M.; Cicone, F.; Del Mastro, C.; Montani, L.; Chiacchiararelli, L.; Monaco, A.; Santini, E.; Scopinaro, F. A theoretical dose-escalation study based on biological effective dose in radioimmunotherapy with 90 Y-ibritumomab tiuxetan (Zevalin). Eur. J. Nucl. Med. Mol. Imaging 2010, 37, 862–873. [Google Scholar] [CrossRef] [PubMed]
Stella, M.; Braat, A.J.; van Rooij, R.; de Jong, H.W.; Lam, M.G. Holmium-166 radioembolization: Current status and future prospective. CardioVascular Interv. Radiol. 2022, 45, 1634–1645. [Google Scholar] [CrossRef] [PubMed]
Fowler, J. The first James Kirk memorial lecture. What next in fractionated radiotherapy? Br. J. Cancer. Suppl. 1984, 6, 285. [Google Scholar] [PubMed]
Dale, R.; Plataniotis, G.; Jones, B. A generalised method for calculating repopulation-corrected tumour EQD2 values in a wide range of clinical situations, including interrupted treatments. Phys. Medica 2024, 118, 103294. [Google Scholar] [CrossRef] [PubMed]
Bisello, S.; Cilla, S.; Benini, A.; Cardano, R.; Nguyen, N.P.; Deodato, F.; Macchia, G.; Buwenge, M.; Cammelli, S.; Wondemagegnehu, T.; et al. Dose–volume constraints fOr oRganS At risk In Radiotherapy (CORSAIR): An “All-in-One” multicenter–multidisciplinary practical summary. Curr. Oncol. 2022, 29, 7021–7050. [Google Scholar] [CrossRef] [PubMed]
Marks, L.B.; Bentzen, S.M.; Deasy, J.O.; Bradley, J.D.; Vogelius, I.S.; El Naqa, I.; Hubbs, J.L.; Lebesque, J.V.; Timmerman, R.D.; Martel, M.K.; et al. Radiation dose–volume effects in the lung. Int. J. Radiat. Oncol. Biol. Phys. 2010, 76, S70–S76. [Google Scholar] [CrossRef] [PubMed]
Milano, A.; Gil, A.V.; Fabrizi, E.; Cremonesi, M.; Veronese, I.; Gallo, S.; Lanconelli, N.; Faccini, R.; Pacilio, M. In Silico Validation of MCID Platform for Monte Carlo-based voxel dosimetry applied to ⁹⁰Y-radioembolization of liver malignancies. Appl. Sci. 2021, 11, 1939. [Google Scholar] [CrossRef]
Yu, Q.; Khanjyan, M.; Fidelman, N.; Pillai, A. Contemporary applications of Y90 for the treatment of hepatocellular carcinoma. Hepatol. Commun. 2023, 7, e0288. [Google Scholar] [CrossRef]

Figure 1. VSV kernels (available as Supporting Material of this paper) represented as a plot of the AD per unit decay to the target voxel (y-axis) versus the source–target voxel distance (x-axis) for the soft tissue (a) and lung tissue (b) for ¹⁶⁶Ho on a square voxel of 2.21 mm side.

Figure 2. Correlation plot of the mean absorbed dose in the lungs (

\bar{A D}

) per GBq of administered activity, obtained from MC simulations with the reference phantom (x-axis), compared with those obtained using the methods listed in Table 2 (y-axis). Each point in the data series represents an increasing LS value (10%, 20%, 30%, and 40%), with a line representing the linear interpolation of each dataset, provided as a qualitative visual guide only.

Figure 2. Correlation plot of the mean absorbed dose in the lungs (

\bar{A D}

Figure 3. Example slices in the coronal view of the AD spatial distributions (left) for

{kLT}_{L}

(a),

{kST}_{L}

(b), and MC (c), along with their respective color scales, are shown, whereas the plot (right) reports the corresponding DVH (d). All dosimetric approaches show significant heterogeneity in the AD spatial distribution. The values of the AD maps and the DVH are given in Gy per GBq of administered activity.

Figure 3. Example slices in the coronal view of the AD spatial distributions (left) for

{kLT}_{L}

(a),

{kST}_{L}

Figure 4. Example slices in the coronal view of the ADr maps (left) of ¹⁶⁶Ho (a) and ⁹⁰Y (b), along with the corresponding DrVHs (c), are shown for LS = 10%. The data demonstrate a different degree of inhomogeneity between the two radionuclides, due to the distinct physical characteristics of their decay spectra.

Figure 5. Cumulative DVHs of

E Q D_{2}

for ⁹⁰Y (blue band) and ¹⁶⁶Ho (yellow band) for the 0.5–1 h range of

T_{μ}

along with the volumetric constraints (red dots) listed in Table 4.

Figure 5. Cumulative DVHs of

E Q D_{2}

for ⁹⁰Y (blue band) and ¹⁶⁶Ho (yellow band) for the 0.5–1 h range of

T_{μ}

along with the volumetric constraints (red dots) listed in Table 4.

Figure 6.

\bar{A D}

NTCP model for RP incidence in partial lung irradiation treatments from EBRT as reported in QUANTEC [42]. The black solid line is the logistic model according to Equation (6) with parameters

b_{0} = - 3.87

and

b_{1} = 0.126 {Gy}^{- 1}

, the blue dots refer to ⁹⁰Y cases, and the yellow ones to ¹⁶⁶Ho cases, each for the LS = 10% case and for the labeled

T_{μ}

Figure 6.

\bar{A D}

NTCP model for RP incidence in partial lung irradiation treatments from EBRT as reported in QUANTEC [42]. The black solid line is the logistic model according to Equation (6) with parameters

b_{0} = - 3.87

and

b_{1} = 0.126 {Gy}^{- 1}

, the blue dots refer to ⁹⁰Y cases, and the yellow ones to ¹⁶⁶Ho cases, each for the LS = 10% case and for the labeled

T_{μ}

Table 1. Summary of the main characteristics of the density distribution of the phantom’s lungs, including mean, median, and most probable value (MPV). All reported values are in g/cm³. The distribution is positively skewed, with most of the voxels belonging to the low-density region. The lungs’ volume is 3292 mL with a corresponding mass of 727 g.

Minimum	Maximum	Mean	Median	MPV
0	1.06	0.221	0.179	0.130

Table 3. Summary of the simulation characteristics used in this paper. Further details can be found in a previously published paper [9].

Parameter	Value
Volume Parametrization	Nested sampling algorithm
Physics interactions model	emstandard-opt3
Range Cut	0.1 mm
Number of Primaries	$10^{9}$ for $L S \neq 0 %$ $10^{10}$ for $L S = 0 %$

Table 4. Lungs’ volumetric absorbed dose constraints for conventional fractionation from the CORSAIR summary paper [36].

Absorbed Dose Constraint
$V_{40 Gy} < 10 %$
$V_{30 Gy} < 15 %$
$V_{20 Gy} < 20 %$
$V_{10 Gy} < 40 %$
$V_{5 Gy} < 50 %$
$D_{mean}$ < 20 Gy

Table 5. Summary of the internal validation of the VSV kernels. The

\bar{A D}

(Gy/GBq) is reported for the liver, assuming a uniform density of 1.04 g/cm³ for the ST kernel, and for the lungs, assuming a uniform density of 0.296 g/cm³ for the LT kernel. Results from the MC simulations are compared with those from the convolution method using the calculated kernels, along with the corresponding RD between the two methods.

Table 5. Summary of the internal validation of the VSV kernels. The

\bar{A D}

Tissue	MC	Convolution	RD
Soft Tissue	10.41	10.36	−0.5%
Lung Tissue	1.39	1.41	1.6%

Table 6. Lungs’

\bar{A D}

relative differences for all the tested methods with respect to the reference MC simulations for each LS case.

Table 6. Lungs’

\bar{A D}

relative differences for all the tested methods with respect to the reference MC simulations for each LS case.

Lung Shunt	MIRD	kST	${kST}_{L}$	kLT	${kLT}_{221}$	${kLT}_{L}$
10%	−65%	−92%	−56%	−77%	−69%	−61%
20%	−64%	−93%	−56%	−76%	−68%	−60%
30%	−64%	−93%	−56%	−76%	−68%	−59%
40%	−64%	−93%	−56%	−76%	−68%	−59%

Table 7. Lungs’

\bar{A D}

of the reference MC simulation along with the tested methods’

\bar{A D}

for all the LS cases in units of Gy per GBq of administered activity, as visually represented in Figure 2.

Table 7. Lungs’

\bar{A D}

of the reference MC simulation along with the tested methods’

\bar{A D}

for all the LS cases in units of Gy per GBq of administered activity, as visually represented in Figure 2.

Lung Shunt	MC	MIRD	kST	${kST}_{L}$	kLT	${kLT}_{221}$	${kLT}_{L}$
10%	6.2	2.2	0.5	2.7	1.4	1.9	2.4
20%	12.0	4.4	0.9	5.3	2.8	3.8	4.9
30%	17.9	6.5	1.3	7.9	4.3	5.7	7.3
40%	23.9	8.7	1.8	10.5	5.7	7.6	9.7

Table 8. MC

\bar{A D}

results for all the LS cases and the corresponding activity as high as safely administrable (AHASA) and the maximum lung activity (MLA). All the results are reported in Gy per GBq of administered activity.

Table 8. MC

\bar{A D}

Lung Shunt	MC $\bar{AD}$ (Gy/GBq)	AHASA (GBq)	MLA (GBq)
10%	6.2	4.88	0.49
20%	12.0	2.50	0.50
30%	17.9	1.67	0.50
40%	23.9	1.26	0.50

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

© 2025 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

d’Andrea, E.; Politano, A.; Cassano, B.; Lanconelli, N.; Cremonesi, M.; Patera, V.; Pacilio, M. The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres. Appl. Sci. 2025, 15, 958. https://doi.org/10.3390/app15020958

AMA Style

d’Andrea E, Politano A, Cassano B, Lanconelli N, Cremonesi M, Patera V, Pacilio M. The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres. Applied Sciences. 2025; 15(2):958. https://doi.org/10.3390/app15020958

Chicago/Turabian Style

d’Andrea, Edoardo, Andrea Politano, Bartolomeo Cassano, Nico Lanconelli, Marta Cremonesi, Vincenzo Patera, and Massimiliano Pacilio. 2025. "The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres" Applied Sciences 15, no. 2: 958. https://doi.org/10.3390/app15020958

APA Style

d’Andrea, E., Politano, A., Cassano, B., Lanconelli, N., Cremonesi, M., Patera, V., & Pacilio, M. (2025). The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres. Applied Sciences, 15(2), 958. https://doi.org/10.3390/app15020958

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres

Abstract

1. Introduction

2. Materials and Methods

2.1. Anthropomorphic Voxelized Phantom and Common Dosimetric Methods

2.2. Monte Carlo Simulations

2.3. Evaluation of Voxel-Based Dosimetric Approaches

3. Results

3.1. VSV Kernel for Soft and Lung Tissue for ¹⁶⁶Ho

3.2. ¹⁶⁶Ho Monte Carlo Simulations vs. “Classical” Dosimetric Approaches

3.3. ¹⁶⁶Ho Versus ⁹⁰Y Monte Carlo Simulations

3.4. Impact on Clinical Decision-Making: From the $\bar{A D}$ to the EBRT Dosimetric Constraints for DVH

4. Discussion

5. Conclusions

Supplementary Materials

Author Contributions

Funding

Institutional Review Board Statement

Informed Consent Statement

Data Availability Statement

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI

Article Menu

The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: 166Ho Microspheres

Abstract

1. Introduction

2. Materials and Methods

2.1. Anthropomorphic Voxelized Phantom and Common Dosimetric Methods

2.2. Monte Carlo Simulations

2.3. Evaluation of Voxel-Based Dosimetric Approaches

3. Results

3.1. VSV Kernel for Soft and Lung Tissue for 166Ho

3.2. 166Ho Monte Carlo Simulations vs. “Classical” Dosimetric Approaches

3.3. 166Ho Versus 90Y Monte Carlo Simulations

3.4. Impact on Clinical Decision-Making: From the A D ¯ to the EBRT Dosimetric Constraints for DVH

4. Discussion

5. Conclusions

Supplementary Materials

Author Contributions

Funding

Institutional Review Board Statement

Informed Consent Statement

Data Availability Statement

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI

The Essential Role of Monte Carlo Simulations for Lung Dosimetry in Liver Radioembolization—Part B: ¹⁶⁶Ho Microspheres

3.1. VSV Kernel for Soft and Lung Tissue for ¹⁶⁶Ho

3.2. ¹⁶⁶Ho Monte Carlo Simulations vs. “Classical” Dosimetric Approaches

3.3. ¹⁶⁶Ho Versus ⁹⁰Y Monte Carlo Simulations

3.4. Impact on Clinical Decision-Making: From the $\bar{A D}$ to the EBRT Dosimetric Constraints for DVH