Fluoroscopically guided interventional procedures have been increasing around the world over the past 20 years [1] . These minimally invasive procedures substitute traditional surgery, resulting in reduced risk to the patient. In recent years, technological advancements have led to new development of diagnostic radiology in practice to provide high-quality images with optimization in radiation dose. However, the prolonged and frequent occupational exposure for interventional cardiologists/radiologists put them as one of the highest exposed group among medical staff, which can be high enough to cause concern [1, 2, 3] .
In these exposure scenarios, staff is exposed to the scatter radiation field around the patient. Specifically, in unshielded situation, interventional radiologists are exposed to an inhomogeneous radiation field where the exposure is greatest below the table, less at the operators waist level, and least at the eye level. The operators eyes can be exposed to a significant dose in some situations (large patient, high-dose imaging protocol). In practice, protective shielding is often used such as lead aprons, table-mounted shields for lower body protection, or ceiling suspended shields. The magnitude and distribution of scatter radiation at the interventionalist position depends on various x-ray beam settings and imaging conditions. The radiation dose received by the interventionalist can vary for the same type of procedure and for similar dose to the patient [3] , which shows the complexity of the scatter radiation field in practice. Several variables affect the distribution of doses at the area surrounding the patient. These variables include the angiography beam settings such as beam added filtration, field diameter, positioner rotations, etc. along with non-controllable clinical situations such as patient abdomen thickness, position, and operator position. Evaluating the impact of each parameter is important to understand the occupational exposure in Interventional Radiology (IR).
The basic principle when selecting the x-ray beam settings in radiography is ensuring an appropriate image quality for the purpose of the procedure and, at the same time, the exposure should generate the lowest harm to the patient. From the image generation point of view, sufficient contrast and low noise are needed [4, 5, 6] , which depend on the tissue to be observed (bone and soft tissues require different tube voltage) and the anatomical region to be explored (depending on the thickness of the body part examined tube voltage and tube current will be different). Tube voltage typically determines the maximum energy of the photons generated and, if it is high enough, the presence of characteristic x-rays (e.g. for a typical tungsten anode, k and k characteristic
x-rays appear around 59 and 69 keV respectively). Adding filtration to the beam reduces the number of low-energy photons of the spectrum, thus increasing the mean energy of the resulting x-ray beam. As a principle of dose optimization to the patient, a minimum beam filtration (2.5 mm Al equivalent) is always required to attenuate very low energy x-rays that are absorbed in the patient with minimal diagnostic benefit, but for some examinations additional filtration is also used. Subsequently, changing kVp or adding beam filtration changes the energy distribution of the scatter photons at the position of the medical staff.
On the other hand, the geometrical configuration of the beam, determined by the beam diameter, the beam angular projection, and the source to detector distance, can also influence the penetration and the scattering around the patient.
These geometrical parameters influence mainly the angular distribution of the scatter radiation. Although in interventional radiology (IR) some of these x-ray beam settings (e.g. tube voltage, total filtration) are automatically selected by an automatic exposure control system according to the medical procedure to be carried out and patient thickness, other parameters (angular projection, field size, source to detector distance) are selected by the medical staff. The goal of this paper is to quantify these influences, since the knowledge of them would help to reduce operators doses.
At present, personal dosimetry in IR is typically performed by issuing staff with physical dosemeters that give measurements of operational quantities, H p (d), as an estimate of protection quantities, like the effective dose. Unfor-tunately, large uncertainties with measurements from current dosemeters exist especially in highly inhomogeneous radiation fields in IR. In ICRP publication 75 [7] , it is stated that for personal dosimetry, the over-all uncertainty at the 95% confidence level in the estimation of effective dose around the relevant dose limit may well be a factor of 1.5 in either direction. This is certainly the case for IR workplaces, where the energy and the orientation of the radiation field are rarely known and can vary significantly. The use of protective garments in IR adds a lot of complexity to personal dosimetry because of complete or partial shielding of the organs. For example, practitioners use protective lead aprons as a common practice, which reduces doses significantly for the organs underneath. Currently, several methods have been proposed to estimate the effective dose when lead aprons are used, such as applying a correction factor to the dose measured by a H p (10) dosemeter worn above or below the apron, or the use of a correlation combining the dose from a dosemeter worn above and another one worn below the lead apron. Nevertheless, it remains difficult to provide a conservative estimate with minimized overestimation of the effective dose under all possible exposure conditions [8] .
On the other hand, computational dosimetry has evolved rapidly in the past decades. The use of the Monte-Carlo (MC) method in radiation dosimetry has increased exponentially. This technique allows the simulation of the energy, position and direction of the particles. Hence, absorbed dose and other dosimetric quantities like fluence, kerma, etc. can be calculated directly. Although computer power restricted early applications to simple geometries, like solid prisms or cylinders, the availability of increased computing power nowadays allows the simulation of complex 3-D geometries like a detailed anthropomorphic human body phantom. Schueler et al. (2006) [9] conducted an experiment mimicking various clinical imaging conditions and using a patient-equivalent abdomen physical phantom. The study investigated how operator exposure in interventional radiology is affected by various common fluoroscopic imaging conditions. In particular, the effect of changing the field of view and the effect of adding copper filters; [13] , within the activities of EURADOS Working Group 12, evaluated the exposure condition of the medical staff in few selected interventional radiology procedures using Monte-Carlo simulations. Unlike [10] and [11] , the study evaluated simulated H p (10) and effective dose rather than extremities and eye lens; and comparing the results with measurements within the study and from the study of Olgar et al. [14] . Yet, the study was limited to few cases of changing some projections and beam qualities.
The study of Badal et al. (2013) [15] and Abdelrahman et al. (2020) [16] presented a dose monitoring system based on accurate MC simulations and 3D localization system that can be used to calculate doses for both patient and The main goal of our study is to extend the examination not only to take into account more variables that contribute to the dose, but also to the combinations of these variables over the operating range of fluoroscopy and angiography procedures. In the present study, the examination of different combinations of parameters is focused on the evaluation of the dose levels at the location of the whole body dosemeter. For that purpose, a series of Monte-Carlo simulations were performed. In total, we did the analyses for 5628 different combinations of parameters.
A step further is to use the results of the simulations to build a model by a regression analysis that describes the dose for different combination of parameters within the range of simulated cases. This model could be used to quantify the importance of each variable on the dose. Regression analysis [18] is a powerful statistical method to examine the relationship between one or more dependent variables and independent variable(s).
In interventional radiology procedures, the main elements of the exposure situation are the patient who is exposed to the primary beam, the interventionalist who is mainly exposed to the scatter radiation, the operating table and the image intensifier. Thus, estimation of the personal dose equivalent to the interventionalist requires the use of a model representing the patient lying on the surgical table under an image detector. Besides the geometrical aspects, the angiographic beam can be simulated with an x-ray spectrum located at the position of the x-ray tube and scaled according to the field of view with a conical shape. The generation of the x-ray spectra in tungsten anodes, which is the most common types of anodes in IR set-ups, is performed using validated models found in literature [19, 20] . The characteristics of the x-ray tube taken into account to generate the beam spectra were the following:
• Anode material: Tungsten • Anode angle: 12ž
• Minimum filtration: 2.5 mm Al.
Monte-Carlo simulations were performed using MCNPX v2.7 [21] . MCNPX is a general-purpose Monte-Carlo radiation transport code that tracks nearly all particles at nearly all energies. The approach that was used for the calculation of the staff dose in MCNPX is presented in the following sections. In MCNPX, the reliability of results are assessed based on statistical tests it performs on the tally. One of the most important is the relative error R 1 . R generally must be less than 10% for meaningful results [22] . Thus, the statistical uncertainties in our simulations were limited to relative error R < 10%.
Each studied parameter was set within an operating range found in literature for fluoroscopy and angiography procedures [23] . These typical radiation field characteristics in interventional radiology and cardiology was also supported by statistics from the two hospitals participating in PODIUM (Skåne university hospital in Malmö and St. Jamess hospital in Dublin). The parameters that were studied are the following:
• The x-ray peak tube voltage: 70, 80, 90, 100, 110, and 120 kVp. projections at 15ř, and 30ř as well as combinations of these projections.
• The field size: 20, 25, and 30 cm diameter.
• The source to image intensifier distance (SID): 100, 110, and 120 cm
• The beam quality: without added filtration and with added copper filtration of 0.1, 0.2, 0.3, 0.6, and 0.9 mm thickness.
A representation of a typical geometrical set-up in IR is shown in figure 1 .
In this set-up, the physician/main operator is represented with the Realistic Anthropomorphic Flexible (RAF) phantom [24] . However, for the purpose of the present study and for only calculating the dose at the location of a dosemeter, no phantom was introduced for the operator in the simulations. Also, the patient phantom was simplified as a slab phantom. This simplified geometry was used to optimize the simulation run time, which allowed the extension of the number of simulations that could be achieved. At the same time, these simplifications should not compromise the validity of the model. For that reason, the approximated geometry for the MC simulations was validated against measurements during the PODIUM project as reported in [25] . Figure 2 shows the geometry used in MCNP. The patient was modeled as a prism of 100 cm length, 32 cm width and 23 cm thickness, and composed of the ICRU tissue 2 .
The image intensifier was simulated by a squared layer of 0.3 cm of lead. The patients tabletop and mattress were simulated as a layer of 2.5 mm of aluminum to account for their equivalent attenuation. 
To calculate the absorbed dose that approximate the personal dose equivalent H p (10), the tally F6 was used to score the energy deposition in MeV/g per particle in a volume that is behind 10 mm tissue equivalent volume as described in figure 3. to 120 keV.
Where a is the atomic density (atoms/barn-cm), m is the cell mass (g), E is the energy (MeV), t is the time (s), V is the volume (cm 3 ), is the direction vector, t is the microscopic total cross section (barns), H(E) is the heating response in MeV/g, is the angular flux and Ψ( − → r ,Ω, E, t)) is in (particles/cm 2 /s/MeV/steradian).
The tallied volume is located at a typical position of the operator at (52.5 cm, -20 cm, and 49.5 cm) relative to the iso-center as the origin as shown in figure 2. This position is considered for femoral or radial access procedures; and it was obtained as the most frequent location through several recordings using motion-tracking system as described in [16] . The position is also in agreement with data from literature in the study of Schueler et al. (2006) [9] and Principi et al. [26] .
This calculated absorbed dose at the operator position is given per particle.
Since all calculated doses are given per simulated particles and thus, they need to be normalized to a reference value. This reference value is the kerma area product (KAP) which is experimentally measurable and is reported in the radiation dose structured report (RDSR) defined by DICOM (DICOM PS3.16 2020b)
[27], the most commonly used method of recording and storing dose data. In the simulations, absorbed dose was tallied in a thin layer of air located at the KAP meter location (20 cm from the source), scaled to cover completely the beam depending on the x-ray field size, and rotating together with the source and the image detector.
In particular, the following ratio was calculated:
Absorbed dose at 10mm depth at the dosimeter location Absorbed dose in air at the KAP meter location (2)
The purpose of the regression analysis is to build a machine learning model that is trained using the results of the simulated cases, test its accuracy, and finally obtain the feature importance score for each studied parameter on the on the predictor. Thus, random forest regression was also used for its feature importance function.
The general form of a complete n-degree polynomial regression model in one variable X modelled as:
Where y is the dependent variable, betas (β) are the regression coefficients for different n th powers of the independent variable X and is the residual error.
With the increasing degree of the polynomial, the complexity of the model also increases. Therefore, the value of n must be chosen precisely. If this value is low, then the model will not be able to fit the data properly and if high, the model may overfit the data. Our model was built using the implementation of polynomial regression in Scikit-learn [28] . The function correlates six independent variables (X 1 X 6 ) to the dependent variable (Y ) where: can be calculated as the following:
Where, n is total number of simulated cases.
It is concluded that the 3 rd degree (cubic) polynomial regression yields the optimum accuracy of the model with the least RMSE.
Random forests are among the most popular machine learning methods to model complex dependencies thanks to their relatively good accuracy, robustness and ease of use [29] . The algorithm is an ensemble method that aggregates Mathematically, the random forest algorithm is a type of additive model that predicts by combining decisions from a set of base models; and hence, this group of models can be described as:
Where the final model F is the sum of simple base models f i .
This broad technique of using multiple models to obtain better predictive performance is called model ensembling. In random forests, all the base models are constructed independently using a different subsample of the data as illustrated in figure 5 . Fortunately, many implementations of this ensembling random forest models are available. In this study, the implementation available in Scikit-learn [28] was used. After obtaining the results of the simulations, the dataset for a regression analysis using random forests algorithm can be prepared. The dataset contains the studied parameters as features (X n ) and the corresponding dependent variable (y) as the ratio R calculated by equation 2. This dataset was then split into a training and testing datasets. The training dataset is the data, which is used to build the model, while the testing dataset is the set of data with which the resulted model can be tested. We followed the classical approach of 80/20% training/testing split which is used often in machine learning, and it is based on Pareto principle [30] . Any bias that can results from the random selection of the test dataset should be reflected in the statistical metrics of the resulted model, and in particular in the Mean Squared Error (MSE) 3 . In our case, the process of the split was performed over multiple folds until the model reached a good overall accuracy. In this case, 20% of the results were used for testing the model. The overall accuracy of the model can be calculated with the coefficient of determination 4 (R 2 ). The higher the model accuracy, the better the estimation of the importance coefficients is.
Finally, for regression, the feature importance in Scikit-learn is the variance reduction, which can be calculated with Mean Absolute Error (MAE):
1 n n i=1 |y i − µ| and is the mean given by: 1 n n i=1 y i . The sum of all feature importance scores is 1 such that the higher the score, the more important the feature is.
First, the effect of changing each parameter independently was studied. In the following sections, the magnitude by which each parameter affects the dose to the operator is presented. Moreover, it is also important to study the effect of the parameters when combined which is concluded by the results of the regression analysis. 
The present study examined the effect of six different cases of beam filtration: without added filters, with 0.1mm, 0.2mm, 0.3mm, 0.6mm, and 0.9mm of added copper filters. As the beam filtration is simultaneously changing with kVp in practice, it is useful to show the effect of added filtration on the dose while changing the kVp. Figure 7 (a) shows that the dose increases sub-linearly with the thickness of the copper filter. At lower kVp, the increase in the dose by adding copper filtration is greater than the increase at higher kVp as shown in figure 7(b) . The beam hardening, as a consequence of the additional copper filtration, leads to a higher dose to the operator as the kVp increases. Consequently, the dose can increase by up to a factor 4.3 at 70 kVp by adding copper filters; while at 120kVp, the multiplicative factor goes only up to 2.6. It is important to highlight that in practice the change of the beam quality is also combined with a reduction in the measured KAP, as such when normalizing the results of the Monte-Carlo simulations, the dose for both the patient and the operator can be reduced. In all cases, the dose is also higher for LAO or CRA projections than for RAO or CAU projections. This is due to the shielding effect by either the patient or the image intensifier when using RAO or CAU projections. For this specific location of the dosemeter, the dose decreases to their lowest value at 45ř RAO before it starts to increase again for angles higher than 45ř that is when the image intensifier is no longer shielding the dosemeter.
The study also examined the effect of changing simultaneously both the 
Three field sizes of 20, 25, and 30 cm diameter of the conical beam were studied. These field sizes correspond to circular field areas of 314, 490, 707 cm 2 , respectively projected at the image intensifier. Figure 12 (a) shows that the dose, for this specific configuration, is directly proportional to the field size for all tube voltages, which is found in agreement with the experimental data from Marshal et al. (1992) [31] . For all beam qualities, the doses are higher when a larger field size is used as more scatter radiation is produced when a larger primary beam intersects the patient. The dose increase is greater at higher kVp compared to lower kVp as shown in figure 12 (b). Increasing the field diameter from 20 cm to 30 cm, multiplies the dose by a factor up to 2.5. 
Typically, the source to detector distance (SID) changes for the purpose of magnification and hence it is linked to the field of view as well. Without added collimation and assuming the distance between the source and the isocenter is fixed, the field area decreases when the SID increases. Similar to the change in the field diameter, increasing the distance between the source and the image intensifier increases slightly the dose at the dosimeter position. However, in practice, when the magnification is increased, the radiation dose has to be increased in order to compensate for a smaller area of anatomy. In our simulations, the distance between the source and the patient was fixed which means the distance left between the patient and the detector increased. We studied the effect of increasing the SID on the scatter fraction R and the ratio of R over the projected beam area on the KAP meter (Sv/Gy.cm 2 ). At lower kVps, increasing the source to image intensifier distance from 100 cm to 120 cm has negligible effect on the dose to the operator, while for kVps equal to 120 kV, the scatter fraction R decreases by 15% while the ratio (Sv/Gy.cm 2 ) increases by approximately 20% as shown in figure 13 . Fixing the field of view at 20 cm diameter and increasing the SID, the scatter fraction R slightly decreases. This is due to the decrease in the field area intersecting the patient and the KAP meter.
While, the absorbed dose at the dosimeter location slightly increases as a result of increasing the distance from the patient to the image detector that allows scatterd radiation to reach the dosimeter location. In clinical situations, the change of SID is usually accompanied by changes to the other parameters, but in general, it is always desirable to minimize the distance between the patient and the image intensifier to reduce the dose to the staff.
(a) Non-normalized results (b) Normalized results by 100cm SID Figure 13 : R (Sv/Gy) and R/A(Sv/Gy.cm 2 ) vs SID, for 70 and 120 kVp, with no added copper filtration, 20cm field diameter, and PA projection 3.6. Independent effect of each parameter: summary Table 1 highlights the effect of changing each parameter independently and the ratio with which each parameter can influence the dose. Primary and secondary beam projections were found to have the largest effect on the dose to the operator followed by the effect of changing the x-ray beam quality, while the effect of changing the source to image intensifier distance has the smallest effect. figure 14 where the predicted value is compared to the true value in the test data set. The overall validation score is R 2 =0.92 which is reasonably good as most of the test samples fall within ś 25% error bounds. It means the model can be used to describe the importance of each feature. 
Several studies have examined the influence of some imaging conditions on the occupational exposure in interventional radiology. We evaluated the effect of several parameters extended over the operating range of IR procedures. The study of Schueler et al. [9] evaluated the ratio between air kerma measured at the dosimeter position and the KAP (Gy/Gy.cm 2 ). In their findings, the ratio is proportional to the change in the beam field of view as such increasing the field of view from 20 to 28 cm increases the ratio by more than a factor of 1.5 for specific imaging conditions (80 kV and 3.4 mm Al). In our study and for comparable imaging conditions, increasing the field of view from 20 to 30 cm increases the ratio R (Sv/Gy.cm 2 ) by a factor of 2.5. Moreover, the study showed that by adding 0.2 mm copper filtration, the ratio (Gy/Gy.cm 2 ) increases by 1.6 times, while, in our study, R (Sv/Gy.cm 2 ) was almost doubled for comparable beam parameters. This difference can be attributed to the difference in dosimeter position and the patient thickness in the two models. Our study shows that the dose is directly proportional to both field size and kVp, which is in agreement with the study of Koukorava et al. [10] and the study of Marshal et al. [31] .
Regarding beam projections, we found that the LAO projections (tube closer to the operator) results in higher doses than the RAO projection; while, the cranial projections present significantly higher doses compared to caudal projections.
This finding is in agreement with the study of Koukorava et al. [10] and Carinou et al. [11] . The study of Ferrari et al. 2018 and 2019 [12, 13] evaluated the ratio H p (10)/KAP (Sv/Gy.cm 2 ) for mainly 4 projections ( PA, 25ř LAO, 25ř RAO and 25ř CRA) and 10 different beam qualities. However, the results cannot be directly compared to our results, as the location of the operator in both studies is considerably different. The operator phantom in the model of [13] was placed near the patients abdomen with the hands almost in the x-ray field that led to the dosimeter being located in front of the image intensifier. This location could occur also during jugular/subclavian venous access procedures that deemed to be uncommon in practice [32] . While in our model, the dosimeter was placed at The present study also provides the following model by which the scatter fraction R can be assessed for any given combination of x-ray beam settings within the range studied here:
where T is a term of the cubic equation as given in Appendix A and is the corresponding coefficient as listed in Therefore, this predictive model is valid only within the intervals of the studied parameters and cannot be extrapolated outside those ranges. However, the study here covered a wide range of settings within the operating range in IR; thus, the model is expected to be valid in most of the cases. It is worth mentioning that this absorbed dose is an approximation of the operational quantity Hp(10) as no phantom for the operator was used in the simulations. The contribution of the backscatter component in this set-up is usually in the order of 10%. In the future, the model could be adapted to calculate H p (10) by introducing a correction factor that accounts for the backscatter component at different settings. However, such correction factor would not change the conclusion and the trends reported in this paper, as it is a quasi constant.
The study of the effect of several parameters that affect staff doses at a typical position of the main operator in IR was performed using Monte-Carlo simulations. In this study, we investigated the influence of the angiography parameter to the resulted dose to the operator can be determined using a model ensembling approach such as random forest regression. Accordingly, the results of the random forest regression analysis showed that the contribution of the primary beam projection dominates the change in the operator dose, followed by the contribution of beam secondary projection, field diameter, kVp and added copper filtration. The study provides a predictive model, which can be used to assess the staff dose for any combination of x-ray beam settings within the range studied here. Finally, in the highlights of results, the study could be used for training the staff in interventional radiology to improve their practice by using optimum beam settings, which can help reducing their occupational exposure. Other non-clinical variables such as operator position, dosemeter location, shielding position and orientation, and patient thickness can also affect the dose to the staff in IR. For that reason, another study will follow to study the influence of these parameters.
In this annex, a predictive model using multivariate polynomial regression is described. This model can help the user to assess the ratio R (Sv/Gy) at the given dosemeter location for any combination of x-ray beam settings within the range studied in this paper. Our model was build using the implementation of polynomial regression in Scikit-learn [28] . The function correlates six independent variables (X1 : X6) 