Alternative statistical methods for cytogenetic radiation biological dosimetry
CCornell University Library, arXiv.org, 2014
Alternative statistical methods for cytogenetic radiation biological dosimetry
Krzysztof Wojciech Fornalski PGE EJ 1 Sp. z o.o., Technology and Operations Office, ul. Mysia 2, 00-496 Warszawa, Poland contract of specified task with Central Laboratory for Radiological Protection (CLOR), ul. Konwaliowa 7, 03-194 Warszawa, Poland e-mail: [email protected], [email protected] ABSTRACT
The paper presents alternative statistical methods for biological dosimetry, such as the Bayesian and Monte Carlo method. The classical Gaussian and robust Bayesian fit algorithms for the linear, linear-quadratic as well as saturated and critical calibration curves are described. The Bayesian model selection algorithm for those curves is also presented. In addition, five methods of dose estimation for a mixed neutron and gamma irradiation field were described: two classical methods, two Bayesian methods and one Monte Carlo method. Bayesian methods were also enhanced and generalized for situations with many types of mixed radiation. All algorithms were presented in easy-to-use form, which can be applied to any computational programming language. The presented algorithm is universal, although it was originally dedicated to cytogenetic biological dosimetry of victims of a nuclear reactor accident.
KEY WORDS: biological dosimetry; Bayesian; Monte Carlo; nuclear accident; calibration curve; cytogenetic; radiation; biodosimetry INTRODUCTION
The aim of radiation biodosimetry by cytogenetics is to calculate doses and the associated confidence limits to exposed (or suspected exposed) persons after a radiation accident or incident. Calculating the dose absorbed in the human body is based on the observed chromosomal aberration (e.g. dicentrics) frequency in the lymphocytes of peripheral blood sampled from the exposed person (IAEA, 2001). This process requires the use of the fitted coefficients of the calibration dose-response curve that is produced by exposure of human blood in vitro to doses of the appropriate quality of radiation. For accurate assessment of radiation doses and coefficients of the calibration dose response curve, a large number of mathematical and statistical methods need to be employed. A number of authors have suggested that the Bayesian statistics approach may be useful for analysis of cytogenetic data because it increases both the accuracy and quality assurance of radiation dose estimates (Brame & Groer, 2003; Ainsbury et al., 2013a, 2013b). The objective of this paper is to present classical, Bayesian and Monte Carlo statistical methods which can be used in cytogenetic radiation biodosimetry for fitting the proper calibration curve, estimating the coefficients, selection of the dose-response model and estimating the dose and dose components from mixed radiation. The practical applications of those methods were presented in (Pacyniak et al., 2014). .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 CALIBRATION CURVES
In the production of an in vitro calibration curve for radiation dose assessment, the dose-response data obtained for different blood donors for the aberration induction in control and irradiated lymphocytes are collected and fitted to a linear or linear-quadratic model . In this model, two DNA lesions in the two unduplicated chromosomes are required for producing chromosome interchanges, like dicentrics, and these lesions may arise from one or two independent ionization tracks (Kellerer & Rossi, 1974). Dicentrics produced by one track (single ionization) will have a frequency that is proportional to a linear function of dose ( aD ), whereas dicentrics induced by two tracks or photoelectric cascade will have a frequency proportional to the square of the dose ( bD ). In general, the shape of the dose-response relationship of chromosomal aberrations is strictly connected with the type of radiation and the way of ionization as well as with the linear energy transfer, LET. Classical calibration curves
In the popular literature one can find the calibration curve of chromosomal aberrations for neutron irradiation, as (IAEA, 2001; SzΕuiΕska et al., 2005): Y π (π· π ) = π + πΌ π· π (1) as well as for gamma radiation, as (IAEA, 2001; SzΕuiΕska et al., 2005): Y π (π· π ) = π + π½ π· π + πΎ π· π2 (2) Equations (1)-(2) can be written in joined form (when the body is irradiated by mixed n+Ξ³ field) as: Y π+π (π· π+π ) = π + π π· π+π + π π· π+π2 (3) Assuming that both radiation qualities are additive in the production of chromosomal damages, the dose-response relationship of chromosome aberrations may be described by: Y π+π (π· π , π· π ) = π + πΌ π· π + π½ π· π + πΎ π· π2 β‘ π¦ π (4) which is usually called a combined linear-quadratic equation for receiving the frequency of chromosomal aberrations y f after irradiation of mixed dose D n +D g . Parameters Y , Ξ± , Ξ² and Ξ³ are usually found as results of the regression analysis. The parameter y f can be written as a ratio of u/w , where u represents the number of chromosomal aberrations and w β the number of cells. Saturated calibration curves
All the equations (1)-(4) are widely used because of their simplicity and practicality. However, for high doses of several greys and more, linearity and parabolicality are not In general mathematical terminology, the name βlinear-quadraticβ for equations (2)-(4) is incorrect; the correct name is βquadraticβ or βparabolicβ .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 preserved (Sasaki, 2003). It is for this reason that cytogeneticists more experienced in the use of the mathematical methods can use a saturated version of equations (1)-(4). They are more general because equations do not trend to infinity ( lim π·ββ Y β (D) β β ), but the results of low and medium dose calculations are the same. Therefore, the linear function of eq. (1) can be replaced by a quasi-linear function: Y πβ (π· π ) = (π πππ₯ β π ) β ( 1 β π βπΌ π· π ) + π (5) Similarly, the linear-quadratic function of eq. (2) can be replaced by a sigmoid one: Y πβ (π· π ) = (π πππ₯ β π ) Β· ( βπ½ π· π β πΎ π· π2 ) + π (6) The sigmoid function in eq. (6) can also be replaced by an Avrami sigmoid critical function, which is more adequate to radiation induced damages: Y πββ (π· π ) = (π πππ₯ β π ) Β· ( βππ· ππ ) + π (7) The Y max in equations (5)-(7) represents the maximum possible number of chromosomal aberrations per cell (in many cases one can assume simply Y max =1 ), and Y is the natural (non-radiation induced) level of spontaneous aberrations per cell. Therefore, Y can correspond to the experiment average Y β 0.0005 for dicentrics (SzΕuiΕska et al., 2005). Calibration curves for extreme doses
The irradiation by extreme doses, over a dozen greys, is rather an academic or laboratory case. However, in some situations saturated curves presented as eq. (5)-(7) are not accurate. When the dose increases over a certain critical point, the frequency of chromosomal aberrations can decrease due to cell death (Sasaki, 2003). In such a situation it is better to use the curve which is linear (linear-quadratic) in small and medium doses, saturates to critical point and decreases for highest doses, as in: Y π (π· π ) = (π πππ₯ β π ) β πΌ π· π β π βπΌ π· π + π (8) and Y π (π· π ) = (π πππ₯ β π ) Β· (π½ π· π + πΎ π· π2 ) Β· π βπ½ π· π β πΎ π· π2 + π (9) However, equations (5)-(9) can be used in special cases only, such as irradiation by high doses. For lower doses equations (1)-(4) are more common. Especially combined linear-quadratic eq. (4) will be used in the further part of this paper. FITTING METHODS OF CALIBRATION CURVES
In the previous section the presentation covered the classical and modified calibration curves. In practice, the calibration curves are found due to the regression analysis methods of fitting the proper curve to specific experimental data. In the following section the classical Gaussian method is reminded and the robust Bayesian method is introduced. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Gaussian best fit
The method of curve to data fitting (regression analysis) used most often is a least squares method (Wolberg, 2005) connected with the general maximum likelihood method (IAEA, 2001). Just to remind, the maximum likelihood method provides the probability distribution (or probability density function, PDF) of e.g. a normal (Gaussian) distribution for all N experimental data points ( D i , E i ) with vertical uncertainties Ο each , where Y i represents a proposed model (calibration curve): π = β π πππ=1 = β β2π exp [β (π π βπΈ π ) ] ππ=1 (10) The maximization of eq. (10) due to modelβs parameters will give the best fit of proper curve to the data points. However, for simplicity of analytical solutions the general maximum likelihood method allows to use the logarithm of P and the sum instead of a product. Thus, the maximization of P is equivalent to the minimization of the least square function: π = β (π π βπΈ π ) π ππ=1 (11) which is called the least squares method. However, it was assumed that only the symmetrical vertical uncertainties ( Ο ) are taken into account (like presented ones in Fig. 1). In general situation with vertical ( Ο y,i ) and horizontal ( Ο x,i ) uncertainties, the denominator of eq. (11) can be replaced by π π¦,π2 + ( ππ π ππ· π ) π π₯,π2 . For simplicity, throughout the presented paper the horizontal uncertainties are assumed to be insignificant and Ο x,i =0. The least squares method is simple and well accepted worldwide because the best fit to the experimental data points can easily be found. The method identifies the proper values of curve fitting parameters (equations (1)-(9)) very fast and effectively. In the particular case of biological dosimetry the maximum likelihood method can be applied also to the Poisson distribution (Groer & Pereira, 1987; El-Sayyad, 1973). This method was used for neutron dosimetry for higher doses, where the background term Y from eq. (1) can be omitted. The probability distribution of the slope Ξ± (eq. (1)) can be presented as: π(πΌ) β β (πΌ π· π ) π’ π ππ=1 Γ exp (β β π€ π πΌ π· πππ=1 ) (12) where u i represents the number of chromosomal aberrations in i -th sample of w i cells after receiving the dose D i . Finding the maximum of eq. (12), one can easily calculate the best fit of Ξ± parameter. All methods reminded above can be widely found in literature (Wolberg, 2005). However, those methods fail in the case of large scatter of experimental points and/or when at least one outlier point exists. In that case the more universal method is the robust Bayesian regression (fit) analysis described in the textbook by Sivia and Skilling (2006) and applied in authorβs papers (Fornalski & DobrzyΕski, 2009, 2010a, 2010b, 2011; Fornalski et al., 2010). One has to note that in classical Gaussian regression all uncertainties of points, Ο , are the same .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Bayesian best fit
The simple comparison between the robust Bayesian and least squares methods is presented in Fig. 1. One can clearly see that outliers make least squares method very misleading while Bayesian fit copes well and follows the main trend. This results from the fact that each i -th data point can be presented as a probability density function (PDF) composed of a proper Gaussian distribution (so called likelihood function, see eq. (10)) around its expected value as well as the prior function for its probability Ο i : π π = β« π β2π π β (ππβπΈπ)22ππ2 Γ π π π2 ππ πβπ (13) The right-side prior function for Ο i in eq. (13) assumes that the i -th analyzed probability Ο i lies between the original one ( Ο ) and infinity. The procedure makes all outliers insignificant as input to the whole posterior probability distribution P for all N points, where, according to the maximum likelihood method described earlier (Fornalski et al., 2010) one can use a sum instead of a product: π = β π π β π = β ln π π (14) where P i is a result of the integration of eq. (13) for single point i , see eq. (23). Figure 1. The example of robust Bayesian (black solid line) and least squares (grey dashed line) fits to some virtual experimental data with three outliers (outstanding points). More examples in (Fornalski et al., 2010).
After the differentiation of logarithmic probability S over all n fitting parameters Ξ»={Ξ» , Ξ» , β¦, Ξ» n } , one can find the final and general form of a Bayesian fitting equation (Sivia & Skilling, 2006): .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 ππππ = β π π (π π β πΈ π ) ππ π ππ β‘ 0 ππ=1 (15) where weights g i of the points are: π π = π βπΈ π ) [2 β (π π βπΈ π ) π Γ (ππβπΈπ)22π0π2 )β1 ] (16) The equation (15) can be implemented directly into the computational algorithm to find the best Bayesian fit to all N experimental data points ( D i , E i ) with vertical uncertainties Ο each. The Y i means the theoretical shape of the best fit with n fitting parameters, Y i (Ξ» , Ξ» , β¦, Ξ» n ) , for example eq. (1)-(9). The result of the practical application (when eq. (1) for Y i is used) is presented in Fig. 1. Exemplary solution of eq. (15), where Y i is given by eq. (4), is presented below. In this special case of mixed n+Ξ³ radiation field (eq. (4)), the eq. (15) is described by the set of n=4 dependent equations ( Ξ»={Y , Ξ± , Ξ² , Ξ³} ): { β π π (π +πΌ π· π,π + π½ π· π,π +πΎ π· π,π2 βπΈ π )=0 ππ=1 β π π (π +πΌ π· π,π + π½ π· π,π +πΎ π· π,π2 βπΈ π ) π· π,π =0 ππ=1 β π π (π +πΌ π· π,π + π½ π· π,π +πΎ π· π,π2 βπΈ π ) π· π,π =0 ππ=1 β π π (π +πΌ π· π,π + π½ π· π,π +πΎ π· π,π2 βπΈ π ) π· π,π2 =0 ππ=1 (17) The next step is to solve the set of eq. (17) to find all four Ξ» parameters. The easiest way is to use the Cramerβs rule and calculate determinants of matrix: π = det [ β π πππ=1 β π π π· π,πππ=1 β π πππ=1 π· π,π β π πππ=1 π· π,π2 β π πππ=1 π· π,π β π π π· π,π2ππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 π· π,π2 π· π,π β π π π· π,πππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 π· π,π2 β π πππ=1 π· π,π π· π,π2 β π π π· π,π2ππ=1 β π π π· π,π3ππ=1 β π πππ=1 π· π,π3 β π πππ=1 π· π,π4 ] (18a) π π = det [ β π π πΈ πππ=1 β π π π· π,πππ=1 β π πππ=1 πΈ π π· π,π β π πππ=1 π· π,π2 β π πππ=1 π· π,π β π π π· π,π2ππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 π· π,π2 π· π,π β π π πΈ π π· π,πππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 πΈ π π· π,π2 β π πππ=1 π· π,π π· π,π2 β π π π· π,π2ππ=1 β π π π· π,π3ππ=1 β π πππ=1 π· π,π3 β π πππ=1 π· π,π4 ] (18b) π πΌ = det [ β π πππ=1 β π π πΈ πππ=1 β π πππ=1 π· π,π β π πππ=1 πΈ π π· π,π β π πππ=1 π· π,π β π π π· π,π2ππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 π· π,π2 π· π,π β π π π· π,πππ=1 β π π πΈ π π· π,πππ=1 β π πππ=1 π· π,π2 β π πππ=1 πΈ π π· π,π2 β π π π· π,π2ππ=1 β π π π· π,π3ππ=1 β π πππ=1 π· π,π3 β π πππ=1 π· π,π4 ] (18c) π π½ = det [ β π πππ=1 β π π π· π,πππ=1 β π πππ=1 π· π,π β π πππ=1 π· π,π2 β π π πΈ πππ=1 β π π π· π,π2ππ=1 β π πππ=1 πΈ π π· π,π β π πππ=1 π· π,π2 π· π,π β π π π· π,πππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 π· π,π2 β π πππ=1 π· π,π π· π,π2 β π π πΈ π π· π,πππ=1 β π π π· π,π3ππ=1 β π πππ=1 πΈ π π· π,π2 β π πππ=1 π· π,π4 ] (18d) π πΎ = det [ β π πππ=1 β π π π· π,πππ=1 β π πππ=1 π· π,π β π πππ=1 π· π,π2 β π πππ=1 π· π,π β π π πΈ πππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 πΈ π π· π,π β π π π· π,πππ=1 β π π π· π,π π· π,πππ=1 β π πππ=1 π· π,π2 β π πππ=1 π· π,π π· π,π2 β π π π· π,π2ππ=1 β π π πΈ π π· π,πππ=1 β π πππ=1 π· π,π3 β π πππ=1 πΈ π π· π,π2 ] (18e) .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 The request fitting parameters Ξ» can be calculated as Y = W Yo /W , Ξ± = W Ξ± /W , Ξ² = W Ξ² /W and Ξ³ = W Ξ³ /W . However, such parameters have to be iteratively calculated using a computational algorithm , because all Ξ» are also implicit in g i (Ξ») (eq. (16)). Finally, the results are all n=4 fitting parameters Ξ»={Y , Ξ± , Ξ² , Ξ³} of the curve Y n+g (D) given by eq. (4). Uncertainties
The fitting parameters Ξ»={Ξ» , Ξ» , β¦, Ξ» n } found thanks to the eq. (15) have their own uncertainties Ο Ξ» ={Ο , Ο , β¦, Ο n } , which can be estimated using the Hessian matrix, H (Fornalski & DobrzyΕski, 2009, 2010b). Thus, the uncertainty of exemplary Ξ» n parameter can be calculated from n -th variance, which is equal to an h n,n element of main diagonal from invertible Hessian matrix, H -1 : π π = β (β1) π» π,π det π» (19) However, the calculation of eq. (19) is rather difficult (but recommended) for high values of n , so Ο n can be sometimes approximated using the CramΓ©r-Rao theorem (see also eq. (42)): π π β₯ π (20) where π π = π πππ π2 = β β [ΞΎ π (π π β πΈ π ) β π π ] ( ππ π ππ π ) (21) and ΞΎ π = π π6 [2 β (2 + π π2 π + π π4 ) exp ( βπ π2 )] β π π2 (22) The residuals R i = Y i β E i can also be used for all previous equations just for simplicity. Additionally, the P i function is the result of an integral from eq. (13): π π = π π π2 β2π [1 β exp ( βπ π2 )] (23) Using eq. (20)-(23) one can calculate the lower bound of uncertainty Ο n of fitting parameter Ξ» n . However, the eq. (19) should be used to have exact values of uncertainties. sometimes, just for simplicity, the variable D n,i can be presented as a function of D g,i , using eq. (30) The actual value of Ξ» is putting into the g i (Ξ») , the next iteration will give more precise Ξ»β , which is putting into the g i (Ξ»β) etc. CramΓ©r-Rao theorem gives lower bounds of the variances of the estimators, that are just the elements of the diagonal of the inverse of the Fisher information matrix. These lower bounds are asymptotically attained by the variances of maximum likelihood estimators, and the inverse of the Fisher information matrix can be estimated from the inverse of the Hessian matrix, H -1 , of the log-likelihood function evaluated at the maximum likelihood estimates .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Model selection
The Bayesian analysis allows also the possibility of relative estimation of the proposed model reliability . Thanks to that one can check which fitted function, for example from eq. (1)-(9), is the best solution for existed data points (Fornalski & DobrzyΕski, 2010b, 2011). The Bayes theorem connects the probabilities of P(Model|Data) ~ P(Data|Model) , which can be used to estimate the relative reliability of two models, M , in the case of the same data, D . To do so, it is necessary to calculate the reliability function for the model M with the fitting parameter Ξ» , using the marginalization procedure: π(π|π·) β π(π·|π) = β« π(π·, π|π) ππ = β« π(π·|π, π) Γ π(π|π) ππ (24) The
P(D|Ξ»,M) corresponds to the likelihood function, represented by the Gaussian distribution around the expected value Ξ» Β± Ο Ξ» with maximum probability of likelihood function equals P(D|Ξ» ,M) . The prior probability P(Ξ»|M) can be assumed as a uniform distribution
U(Ξ» min ,Ξ» max ) . Because such form of P(Ξ»|M) is independent of Ξ» , the integral (24) can be written as: π(π·|π) = πππ₯ βπ πππ β« π(π·|π , π) π β (πβπ0)22ππ2 ππ π πππ₯ π πππ (25) The result of the integral (25) can be approximated by π(π·|π , π) π π β2π because β the sharp cut-offs at Ξ» min and Ξ» max do not cause a significant truncation of the Gaussian β probability distribution from eq. (25) (Sivia & Skilling, 2006). Because Ξ» corresponds to the parameter found by the robust Bayesian best fit method for model M , the maximum value of likelihood function P(D|Ξ» ,M) can be replaced by the set of P i given by eq. (13) or (23) and the final form of the reliability function can be approximated by (Fornalski & DobrzyΕski, 2011): π(π|π·) β π(π·|π) β β π π Γ π π β2ππ πππ₯ βπ πππ (26) The right-hand term in eq. (26) is called an Ockham factor, which prevents the use of overcomplicated models. Equation (26) corresponds to the situation, where model M have only one ( n=1 ) fitting parameter, Ξ» Β± Ο Ξ» . In the case of zero parameters ( n=0 ) the Ockham factor equals 1 and model M is just a constant value ( Y i =const ). In the case of n fitting parameters Ξ»={Ξ» , Ξ» , β¦, Ξ» n } with their estimated uncertainties Ο Ξ» ={Ο , Ο , β¦, Ο n } , the most general form of eq. (26) can be presented as (Fornalski & DobrzyΕski, 2011): π(π|π·) β β
π,π βπΈ π ) [1 β exp (β (π π,π βπΈ π ) )] ππ=1 Γ β π π β2ππ πππ₯ βπ πππ ππ=1 (27) It is important to remember that N represents the number of experimental points ( D i ,E i ) with vertical uncertainties Ο each, to which model M (the Y(D) relationship) is fitted using n fitting parameters Ξ» Β± Ο Ξ» . The most problematic are the values Ξ» min and Ξ» max , which can be calculated independently or taken arbitrary for all Ξ» . In the case of the arbitrary method, they can be taken as minimum/maximum possible values of the considered parameter Ξ» using the largest span that can be tolerated by the data. In order not to extend the range of Ξ» in the Thus, the classical methods of model selection (like AIC) or other Bayesian ones (like BIC) will be omitted here .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 case of a huge scatter of data, not more than three points are allowed to lie outside the range (Fornalski & DobrzyΕski, 2011). One can also find a non-arbitrary way to calculate both parameters Ξ» min and Ξ» max . Such a method should be independent to Bayesian regression, so one can propose for example Ξ» min = Ξ» β k Ο Ο and Ξ» max = Ξ» + k Ο Ο , where Ο Ο means the standard deviation of Gaussian regression (see eq. (10)) and k the degree of belief (e.g. k β‘ 2 ). The eq. (27) can be applied directly in the computational algorithm (as well as the Excelβs formula) to calculate reliability functions P(M=Y(D)|D) for each model, like eq. (1)-(9). After that one can calculate the relative value of each two models, say A and B , to check which of them is more likely to describe the data: πΎ π΄ = π(π=π΄|π·)π(π=π΅|π·) (28) When W M is greater than 1, model A wins over B . When W M β1 , both models have the same degree of belief. In general, W M can quantify the preference of one model with respect to the other one. In practice, the real values of W M may show that the plausibility of models can differ by orders of magnitude (Fornalski & DobrzyΕski, 2011). DOSE ESTIMATION METHODS
Dose estimation methods in biological dosimetry using e.g. Bayesian statistics can be found in the literature (Brame & Groer, 2003; Ainsbury et al., 2013a, 2013b). However, the methods presented below are more wide and can be treated as a significant extension of the general biodosimetry statistics. Having the fitting parameters of the calibration curves estimated by Gaussian or Bayesian method (previous section), one can use
Y(D) functions to estimate the dose ( D ) and/or the frequency of chromosome aberrations ( Y ) after accidental exposure of people to gamma, neutron or to mixed n+Ξ³ radiation (for example after nuclear reactor accident). Because such mixed radiation is composed of particles having different biological effects, there is a strong need to calculate not only the total absorbed dose but also components of the total dose (here D n and D g ). From the observed chromosome aberration frequency in the lymphocytes irradiated with mixed n+Ξ³ radiation it is not possible to discriminate between those aberrations due to Ξ³-rays and those due to neutrons. However, in the case where a physical estimate of the neutron to gamma absorbed doses ratio is known: π = π· π π· π (29) it is possible to estimate the separate neutron and Ξ³-ray doses by iterative method (IAEA, 2001; SzΕuiΕska et al., 2005), which will not be described here. However, if a physical estimate of neutron and photon components of the absorbed dose is not available, some alternative statistical methods can be applied. In further considerations, the eq. (4) for chromosomal aberration frequency, Y , will be taken into account. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Prior functions
Assuming that calibration curve Y n+g (D n ,D g ) is measured precisely, the most important problem is connected with the precise measurement of neutron to gamma absorbed doses ratio (eg. (29)). The easiest case is when the ratio is known exactly or with an uncertainty Ο Β± Ο Ο . However, using Ο which can vary from zero to infinity is not practical. The ΞΈ parameter normalized to the range of [0,1] is recommended as (Brame & Groer, 2003) π = π· π π· π +π· π = (30) Parameter ΞΈ defined by eq. (30) corresponds to the contribution of gamma dose in the total dose and is more practical in use than Ο . However, the problems for Ο appear also for ΞΈ in the context of precise dose measurement. When ΞΈ (or Ο ) is not precisely known, one can use the PDF to estimate the most probable value of ΞΈ (or Ο ). Such a PDF can be called an prior probability p( ΞΈ ) or p( Ο ) and for the most common cases it can be approximated by Gaussian distribution for ΞΈ (or Ο ) with standard deviation of Ο ΞΈ (or Ο Ο ), for example: π(π) = π exp [β (πβπΜ) π2 ] (31) where πΜ represents the expected value with uncertainty Ο ΞΈ . For normalized values [0,1] one can use analogical Gaussian prior for πΜ Β± π π but transformed into ΞΈ coordinates (Brame & Groer, 2003): π(π) = π π exp [ β12π π2 (( β 1) β πΜ) ] (32) The main assumption for both prior functions of eq. (31) and (32) is that information about the expected value and its uncertainty is necessary (so called informative priors ). However, in some cases such data does not exist and one the only information available is that there is a mixed n+Ξ³ radiation field with unknown ΞΈ . In that case one can assume that there is approximately an even proportion of neutron and gamma doses, which correspond to the uninformative prior (from Beta distribution): π(π) = π β π (33) or π(π) = ππππ π‘ (34) when one can say nothing about the gamma to neutron relation. The prior (34) can also be approximated by the uniform probability distribution within an assumed range, p(ΞΈ)=U(ΞΈ min ,ΞΈ max ) . However, the use of prior (34) is not always possible and rather gives the In some cases the proper informative prior distribution can use Beta or Log-normal distribution When our knowledge about ΞΈ is more and more accurate, informative priors would asymptotically aim towards the Kronecker delta (see http://en.wikipedia.org/wiki/Kronecker_delta) .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 least precise results whenever used. It is the reason why one should always try to assess the potential p( ΞΈ ) shape, which is usually possible. Now, having a sufficient information about ΞΈ , one can try to estimate the n+Ξ³ doses and chromosomal aberration frequencies, especially for dicentrics. Five alternative statistical methods are proposed below. Classical method
In the classical method the contribution of gamma dose in the total dose, ΞΈ , is precisely known and ΞΈ is given as an exact value. Thus no prior function for ΞΈ is needed. The most common approach is an iterative algorithm, precisely described in (IAEA, 2001; SzΕuiΕska et al., 2005). In this approach the doses and chromosomal aberration frequencies are estimated in iterative way from eq. (4) and (29). However, this algorithm is usually laborious and requires many repetitions (iterations) to obtain final results. To automate this method, one can use a system of two equations: eq. (30) and exemplary eq. (4): { π¦ π = π + πΌ π· π + π½ π· π + πΎ π· π2 π = π· π π· π +π· π (35) to calculate both doses, D g and D n , in the function of ΞΈ : {π· π (π) = β(πΌ +π½) +4πΎ(π¦ π β π )β(πΌ +π½)2πΎ π· π (π) = π· π (π) (36) in the special case of Y(D) given by eq. (4). It is assumed, that all constants ( Ξ± , Ξ² , Ξ³ , y f , Y and ΞΈ ) are precisely known due to the experimental data with proper uncertainties (see previous section). Treating ΞΈ as a variable, one can present a set of eq. (36) in Fig. 2. Finally, the uncertainties of dose estimations, Ο D , in the presented method can be calculated using the sum of independent finite increments method: π π· π₯ β‘ β | ππ· π₯ ππ π | π πππ₯ π=1 πΏπ π (37) where x ={ g,n } and j max represents the total number of parameters Ξ΅ j Β± δΡ j . Assuming the Y(D) given by eq. (4), parameters Ξ΅ j ={ Ξ± , Ξ² , Ξ³ , y f , Y , ΞΈ } and δΡ j ={Ο Ξ± ,Ο Ξ² ,Ο Ξ³ ,Ο yf ,Ο Y0 ,Ο ΞΈ } can be found experimentally (see previous section). The eq. (37) will give slightly higher values than in the exact differential method. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Figure 2. The graphical representation of eq. (36) for exemplary experimental data Ξ±=0.832, Ξ²=0.0164, Ξ³=0.0492, y f =1.2 and Y =0.0005 (IAEA, 2001). Enhanced classical method (quasi-Bayesian method)
To enhance the classical method for the possibility, that ΞΈ is given not by the value but by the prior function (PDF), one needs to transform classical relationships into the probability distributions . In practice, the set of eq. (35) needs to be solved to find distributions of ΞΈ in the example case of eq. (4): { π π (π· π ) = π· π π· π + (π¦ π β π βπ½π· π βπΎπ· π2 ) π π (π· π ) = βπ½ β4πΎ( π +πΌπ· π βπ¦ π )βπ½βπ½ β4πΎ( π +πΌπ· π βπ¦ π )βπ½+2πΎπ· π (38) The two different designations of ΞΈ result from the fact, that the parameter is not precisely described by a value but by a prior probability function. Thus, the probability distribution of the dose can be written as: π(π) π β² (π·) = π(π π₯ (π· π₯ )) π π₯β² (π· π₯ ) β‘ ππππ π‘ β π(π· π₯ ) (39) where x={g,n} . For example eq. (39) can be presented with exemplary prior from eq. (33) as a system of two probability distributions: {π(π· π ) β π π (π· π ) β π π (π· π ) π(π· π ) β π π (π· π ) β π π (π· π ) (40) presented in Fig. 3. That way the method can also be called the
Bayesian-frequentist hybrid method. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Finally, the estimation of the dose D x can be found from the maximum of the distribution (39), taken directly from the graph (like exemplary Fig. 3 and eq. (40)) or calculated from the first derivate equation: π π(π· π₯ )ππ· π₯ = 0 (41) The exemplary results shown in Fig. 3 correspond well with the
Classical method calculations of eq. (36) for precisely known ΞΈ =0.5 (equals to the expected value of prior (33)). The uncertainties of dose estimations, Ο Dx , in the presented method can be assessed from the shape of distributions or calculated using the CramΓ©r-Rao theorem: π π· π₯ β₯ π2 ln πππ·π₯2 | (42) where ln(P) is a natural logarithm of P(D x ) due to the maximal likelihood method. Figure 3. The graphical representation of eq. (40) for exemplary experimental data Ξ±=0.832, Ξ²=0.0164, Ξ³=0.0492, y f =1.2 and Y =0.0005 (IAEA, 2001) with prior given by eq. (33). Simplified Bayesian method
The Bayesian reasoning can be reduced to the simple probabilities equation:
πππππΈπ πΌππ ππ ππ΅. = πΏπΌπΎπΈπΏπΌπ»πππ· ππ ππ΅. Γ ππ πΌππ ππ ππ΅. (43) This type of equation was formerly used in the case of eq. (13) and (25). In the case of dose estimation, the selection of prior probability function was discussed for equations (31)-(34). The likelihood probability function has to be found using biophysical basis. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 In general, there are w cells and u events (here: chromosomal aberrations). The expected number of events equals wΒ·y f , so one can use the Poisson statistics for likelihood function as: πΏ = (π€ π¦ π ) π’ π βπ€ π¦π π’! (44) Taking equations (4) and (30) for eq. (44), one can calculate the likelihood functions for both n+Ξ³ doses: { πΏ ( π· π | π ) = [ π€ ( π +πΌ π· π +π½π· π +πΎπ· π2 ) ] π’ π’! Γ π βπ€ ( π +πΌ π· π +π½π· π +πΎπ· π2 ) πΏ ( π· π | π ) = [ π€ ( π +πΌπ· π +π½ π1βπ π· π +πΎ( π1βπ π· π ) ) ] π’ π’! Γ π βπ€ ( π +πΌπ· π +π½ π1βπ π· π +πΎ( π1βπ π· π ) ) (45) Having the likelihood functions (45) and the assumed prior function p( ΞΈ ), one can find the posterior probability β the probability distribution of dose : π(π· π₯ ) = β« πΏ(π· π₯ |π) π(π)ππ (46) where x={g,n} . Similarly to the case of Enhanced classical method , one can find the estimation of dose D x from the maximum of the curve (46) or calculating the first derivate equation, given by the same eq. (41). The result of eq. (41) in the context of eq. (46) gives the searched estimated values of doses. However, for some priors, due to the maximum likelihood method, the natural logarithm of P is much easier in analytical calculations. Figure 4. The application of Bayesian method (eq. (46)) to real experimental data, where w=1000, u=33, Ξ±=0.354, Ξ²=0.0119, Ξ³=0.0557 and Y =0.0005 (Pacyniak et al., 2014). The prior (32) with πΜ = 0.92 was used. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 However, the analytical solution of integral (46) is often difficult (or impossible) to obtain. It is then possible to use computational methods of integration, such as the Monte Carlo integration method or iterative integration method. The latter method is quick, however it can cause some instabilities and fluctuations of the final shape of the
P(D x ) curve due to large values of u variable . Finally, the uncertainty of estimated value of dose, Ο Dx , can be approximated using the eq. (42). Fig. 4 contains the application of presented method to the real experimental data. Additionally, the method presented above can also be simplified with the use of the exact value of ΞΈ instead of prior function, p( ΞΈ ). In that way one should put the known ΞΈ directly into eq. (45) and analyze the distribution of L(D x ) to find its maximum for D x estimation. Full Bayesian method
The complete Bayesian approach of dose estimation was proposed by Dr. R.S. Brame and Prof. P.G. Groer (Brame & Groer, 2003). They assumed that not only ΞΈ parameter is given by a prior function (they used p( ΞΈ ) as eq. (32)), but also Ξ± , Ξ² and Ξ³ are given by a certain probability distribution. In general, prior functions for Ξ± , Ξ² and Ξ³ are given by a proper posterior probability distributions (eq. (14)) for fitted parameters of calibration curve. Another words: priors p( Ξ» ) are the results of robust Bayesian fitting method (see previous section). However, Brame and Groer (2003) assumed for simplicity that those priors can be approximated by the Gamma distribution (the logarithm of Gamma distribution has a simple form for analysis, e.g. for maximum likelihood method): π(π) = π πβ1 π§ π Ξ(π) exp(βπ§π) (47) where Ξ»={Ξ±,Ξ²,Ξ³} and Ξ is a gamma function. Parameters k and z are the shape and scale parameters, respectively. Having such an assumption, one can write the PDF of dose as π(π· π₯ ) β β« β« β« β« πΏ(π· π₯ |πΌ, π½, πΎ, π) π(πΌ) π(π½) π(πΎ) π(π) ππΌ ππ½ ππΎ ππ (48) The likelihood function, L , is also given by the Poisson distribution (44). When parameters Ξ» are given in the classical way with some standard deviation, priors p( Ξ» ) can be introduced by a Gaussian distribution, like eq. (31). One has to note that Brame and Groer (2003) did not use the prior function for Y parameter (see eq. (4)), but generally this should also be included. Finally, the eq. (48) can be calculated using computational or analytical methods (e.g. maximum likelihood method), similarly to eq. (46). Details and results of eq. (48) are presented in the paper by Brame and Groer (2003). Large u! gives zero or infinity in some programs or codes, which can cause wrong results .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Monte Carlo method
The Monte Carlo method is a computational combination of classical iterative method and
Enhanced classical method (quasi-Bayesian method) . This method requires a dedicated computer program and an implementation of the algorithm presented in Fig. 5. In general, the algorithm is composed of two major loops: one over w cells and the second over K iterations. The loop over K is dedicated to obtain average results of K independent simulations. Stable results are obtained when K equals several hundred of Monte Carlo repetitions (thus, proper uncertainties can be estimated as a standard deviation of the resultsβ normal distribution). In the single i -th simulation, the current cell is randomized from all w cells. In such a case, the actual value of ΞΈ is taken. It can be an exact value of ΞΈ when such a value is available or is randomized from the probability distribution of p( ΞΈ ). In the latter case, the actual value of ΞΈ j is taken from the classical Monte Carlo randomization, when the maximal value of p( ΞΈ ) equals p max (for normalized probability p max =1 ). The actual value of ΞΈ is used to select the type of interaction: neutron or gamma. For the damage from neutron interaction, the variable u n increases by 1. For the damage from gamma β u g increases by 1. The DNA damage from radiation is a primary cause which can be used to calculate the absorbed dose from u x damages ( x={g,n} ): π· π₯ = π(π’ π₯ ) β ππππ π‘ Β· π’ π₯ (49) The function f from eq. (49) can be linear, quadratic or generally polynomial or saturated. For simple use of the algorithm, one can assume the linear relationship with the slope const β 0.012 calculated from comparison between Monte Carlo results and real data from (IAEA, 2001). Thus, having the values of D n and D g , one can calculate proper chromosomal aberration frequencies Y n and Y g , e.g. from eq. (1) and (2). The loop over w cells continues when the actual sum Y n +Y g is lower than assumed y f . The whole algorithm is presented graphically in Fig. 5. Owing to the presented Monte Carlo method one can tests all statistical aspects of such a virtual irradiation of cells. For example one can add the information about each damage from irradiation, u x , in the dedicated table T of all w cells. The exact number of damages can then be correlated with dose and chromosomal aberrations in single cells. Also one can obtain many useful distributions of cell parameters, such as Poisson distribution of damage in cells. It provides further statistical studies on virtual irradiation of a group of cells. GENERALIZATION
Methods presented in previous sections were introduced for two types of radiation, especially gamma and neutron. However, generally the problem can be enhanced to the situation, where the potential victim is irradiated by many different types of radiation, as during a severe nuclear accident or space travel. The mentioned computational program was created by the author using C++ language .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Figure 5. The detailed scheme of the Monte Carlo methodβs algorithm. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 The generalized calibration curve, combined Bayesian and Gaussian regression fit as well as generalized Bayesian dose estimation methods are introduced below.
Generalized calibration curve
The calibration curves from eq. (1) and (2) can be written in a more general form, as a polynomial
π(π·) = π + π π· + π π· + . . . + π π π· π (50) In the case of many types of radiation, like eq. (4), the eq. (50) can be presented as π(π· π‘ππ‘ππ ) = π + β β π π,π π· ππππ=1π π=1 (51) where n is a degree of a i -th polynomial and R is a number of radiation types. A more complicated form of generalized Y(D) is expected for saturated and critical equations (5)-(9). However, one can find a single general equation, which is common for all presented before, namely eq. (1)-(9) and (50)-(51):
π(π· π‘ππ‘ππ ) = π + (π πππ₯ β π ) β [β (π π,π(π) + π π,π(π) π· ππ π β β π π,π,π(π) π· πππ,π,ππΎπ=1 ) ππ=0 ] π π=1 (52) For Ξ» (a) =Ξ» (b) =0 and j>0 , eq. (52) becomes polynomial (51). For Ξ» (a) =1 , Ξ» (b) =-1 and n=0 , the generalized eq. (52) becomes a sigmoidal one (5)-(7). For Ξ» (a) = Ξ» (b) =Ξ» (c) and j>0 , eq. (52) becomes eq. (8)-(9). Generalized Gaussian and Bayesian fit
Returning to the regression analysis techniques it is possible to observe that the main assumption of the Gaussian method (eq. (10)) is that all points are treated as correct ones. This methodology causes potential misfits when the outliers exists. On the other hand the Bayesian fit assumes that all points can be potentially treated as outliers. However, generally one can propose the posterior probability function, analogically to eq. (10) and (13), which can connect both methods into a single one (Sivia & Skilling, 2003; Box & Tiao, 1968) π π = π β« π©(πΈ π , π π2 ) π π π2 ππ πβπ + (1 β π) π©(πΈ π , π ) (53) where N is a normal (Gaussian) likelihood distribution and Ο is the probability that data E i is an outlier. It is the reason why the left-hand side of eq. (53) is a Bayesian distribution (same as eq. (13)) and the right-hand the Gaussian one (eq. (10)). This approach is called Mixture of distributions (Box & Tiao, 1968) or
The good-and-bad data model (Sivia & Skilling, 2003). Thus the total posterior distribution for all N points can be found analogically as for eq. (14) and (23) as π = β ln { β2π [π π π π2 (1 β exp ( βπ π2 )) + (1 β π) exp ( βπ π2 )]} ππ=1 (54) .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 where R i = Y i β E i . Using the same reasoning of maximum likelihood method one can find the same solution as eq. (15), but with new weights, instead of g i (16), equal (Fornalski & DobrzyΕski, 2010a) π πβ = β exp( βπ π22π0π2 )β[π π0π2π π2 (2 π0π2π π2 +1)β1+π]β2π π0π4π π4 exp( βπ π22π0π2 )β(π π0π2π π2 β1+π)βπ π0π2π π2 (55) The generalized method presented above is a good alternative for using Bayesian and Gaussian (least squares) fitting in the same time. One can see that for Ο =1 the model became a Bayesian regression, while for Ο =0 the model became a classical Gaussian (least squares) one. However, the mixed model works well just for Ο =0.05 (Ekiz, 2002), because usually outlier points are a minority among all experimental data. Moreover, results obtained by weights g i* (55) are very similar to ones obtained by g i (16), but with uncertainties often reduced by about 30% (Fornalski & DobrzyΕski, 2010a). Generalized Bayesian dose estimation method
For a generalization of presented Bayesian methods for dose estimation it is necessary to assume many parameters ΞΈ i for each i type of radiation: π π = π· π β π· ππ π=1 = π· π π· π‘ππ‘ππ (56) and proper prior functions, p( ΞΈ i ). Priors should be assumed or established experimentally, as in eq. (31)-(34). When choosing, for the benefit of simplicity, the polynomial eq. (51), one can present y f dedicated for exact i -th dose as: π¦ π (π· π ) β‘ π(π· π ) = π + β β π π,π [ π π π ((β π· ππ π=1 )βπ· π )] πππ=1π π=1 (57) Each dose D i was written as a proper part of the total dose using eq. (56) and the reasoning used in eq. (45). Next, the likelihood function based on Poisson distribution (44) and y f (D i ) from eq. (57) can be written as: πΏ(π· π |π , π , β¦ , π ππ , π , π , β¦ , π π ) β [π¦ π (π· π )] π’ Γ π β π€ π¦ π (π· π ) (58) Assuming, that all fitting parameters, Ξ» , are given by their priors (including Ξ» =Y ), the posterior probability distribution for each dose equals π(π· π ) β β« β¦ β« πΏ(π· π |π , π ) β (β π(π π )ππ πππ π=0 ) β (β π(π π )ππ ππ π=1 ) (59) which is the most general form of posterior probability for R types of radiation and (n+1)R number of fitting parameters. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 DISCUSSION AND CONCLUSIONS
The presented paper is composed of four major parts: the first presents several forms of calibration curves, the second the robust Bayesian algorithm of fitting the proper curve to the experimental data, the third discusses five potential methods of dose estimation after mixed n+Ξ³ irradiation and the fourth presents certain generalized forms of methods presented before. The robust Bayesian regression analysis method (Bayesian fit) is useful whenever outlier or scattered data points exists. The proposed method is useful to obtain fitting parameters of all potential curves, from simple linear (1) and quadratic (2) relationships to saturated functions (5)-(7) widely found in real conditions (Pacyniak et al., 2014; Dabrowski & Thompson, 1998). The algorithm allows for the calculation of the proper uncertainties of all parameters found while also providing relative plausibility of alternative models (fitted curves). The proposed algorithm is generally useful to find the best calibration curve for biological dosimetry. It was successfully used in practice on many occasions (Fornalski & DobrzyΕski, 2009, 2010a, 2010b, 2011; Fornalski et al., 2010). The five alternative statistical methods of the estimation of n+Ξ³ doses as well as proper chromosomal aberration frequencies were presented next. One can choose between the classical method, the classical method with probability distributions (called also a quasi-Bayesian method), two Bayesian methods and the Monte Carlo one. All methods have their own advantages as well as many inconveniences. It is not generally possible to indicate one method as the best, since this problem depends on the exact situation and expectations. Classical method is widely known and simple, but uses the exact value of the contribution of gamma dose in the total dose, ΞΈ . When this parameter is known only as a proper probability distribution (prior function), p( ΞΈ ), one can use enhanced classical method (called also quasi-Bayesian one) with probability distributions. Alternatively, it is possible to use the simplified Bayesian method, which allows for (as the only one) the use of the most inexact information about prior, given by eq. (34). The full Bayesian method proposed by Brame and Groer (2003) assumes prior function also for fitting parameters from the calibration curve. The last method, the Monte Carlo one, is an interesting computational alternative, most useful for statistical studies on irradiated groups of virtual cells. All presented methods were used in practice on real data (Pacyniak et al., 2014), see also Figs. 2-4. Other, more general forms of methods presented before were also introduced. They can be particularly useful e.g. when the potential victim is irradiated by many types of radiation at the same time, not only gammas and neutrons. However, for practical cytogenetic dose estimation, certainly in an emergency scenario, the generalized methods could also be useful for separating gradient type inhomogeneous exposures, e.g. irradiation of different fractions of the body. The presented statistical methods can be applied into the computational algorithm and be used in cytogenetic analysis. However, biologists and cytogeneticists less experienced in the use of the mathematical methods may need some help of programmers, mathematicians or physicists. Nevertheless, the methods are usually presented in easy-to-use forms, which can be applied even as Excelβs formulas. The presented paper should be treated as a methodology and statistical guide. The description of detailed application of presented methods to dedicated experimental data was intentionally .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 omitted, except for some examples in Figs. 1-4. The application of all methods to real data is the subject of other studies (Pacyniak et al., 2014). ACKNOWLEDGEMENTS
This work was supported by The National Centre for Research and Development (NCBiR), Poland, under the strategic research project: β
Technologies supporting development of safe nuclear power engineering β, research task no. 6: β
Development of nuclear safety and radiological protection methods for the nuclear power engineeringβs current and future needs β, subject no. 8: β
Implementation of the dicentric assay for biological dosimetry following accidental exposure to fission neutrons β. The author wishes to thank Dr Maria Kowalska from Central Laboratory for Radiological Protection (Poland) for the inspiration, careful review of the manuscript as well as many valuable remarks and discussions about the cytogenetics. Additional thanks for Dr Y. Socol, Prof. P.G. Groer, Dr R.S. Brame, Prof. L. DobrzyΕski, Ms. I. Pacyniak and Ms. K. Buraczewska for additional comments.
REFERENCES
Ainsbury E.A., Vinnikov V.A., Puig P., Higueras M., Maznyk N.A., Lloyd D.C., Rothkamm K. (2013a). Review of Bayesian statistical analysis methods for cytogenetic radiation biodosimetry, with a practical example.
Radiat Prot Dosimetry (on-line), doi: 10.1093/rpd/nct301. Ainsbury E.A., Vinnikov V., Puig P., Maznyk N., Rothkamm K., Lloyd D.C. (2013b). CytoBayesJ: software tools for Bayesian analysis of cytogenetic radiation dosimetry data.
Mutation Research , vol. 756(1-2), pp. 184-191. Box G.E.P., Tiao G. C. (1968). A Bayesian approach to some outlier problems.
Biometrika , vol. 55(1), pp. 119-129. Brame R.S., Groer P.G. (2003). Bayesian methods for chromosome dosimetry following a criticality accident.
Radiation Protection Dosimetry
Proceedings of the International Symposium on Health Effects of Low Doses of Ionizing Radiation: Research into New Millenium ; Ottawa, Canada, pp. 63-74. Ekiz U. (2002). A Bayesian method to detect outliers in multivariate linear regression.
Hacettepe Journal of Mathematics and Statistics , vol. 31, pp. 77-82. El-Sayyad G.M. (1973). Bayesian and Classical Analysis of Poisson Regression.
Journal of the Royal Statistical Society . Series B (Methodological), Vol. 35, No. 3, pp. 445-451. Fornalski K.W., DobrzyΕski L. (2009). Ionizing radiation and the health of nuclear industry workers.
International Journal of Low Radiation , vol. 6, no. 1, pp. 57-78.
Fornalski K.W., DobrzyΕski L. (2010a). Zastosowania twierdzenia Bayesa do analizy niepewnych danych doΕwiadczalnych (in Polish).
PostΔpy Fizyki , vol. 61, no. 5, pp. 178-192. Fornalski K.W., DobrzyΕski L. (2010b). The healthy worker effect and nuclear industry workers.
Dose-Response , vol. 8, no. 2, pp. 125-147. .W. Fornalski: Alternative statistical methods for cytogenetic radiation biological dosimetry. arXiv.org, 2014 Fornalski K.W., DobrzyΕski L (2011). Pooled Bayesian analysis of twenty-eight studies on radon induced lung cancers.
Health Physics , vol. 101, no. 3, pp. 265-273. Fornalski K.W., Parzych G., Pylak M., SatuΕa D., DobrzyΕski L. (2010). Application of Bayesian reasoning and the Maximum Entropy Method to some reconstruction problems.
Acta Physica Polonica A , vol. 117, no. 6, pp. 892-899. Groer P.G., Pereira C.A. (1987). Calibration of a radiation detector: chromosome dosimetry for neutrons. In:
Probability and Bayesian Statistics . Ed R. Viertl (New York: Plenum Publishing), pp. 225-232. IAEA (International Atomic Energy Agency). (2001).
Cytogenetic Analysis for Radiation Dose Assessment. A Manual. Technical Reports Series
No. 405, Vienna. Kellerer A.M., Rossi H.H. (1974). The theory of dual radiation action.
Current Topics in Radiation Research , North-Holland and American Elsevier Publishing Companies, vol. VIII, pp. 85-158. Pacyniak I., Fornalski K.W., Kowalska M. (2014). Employment of Bayesian and Monte Carlo methods for biological dose assessment following accidental overexposures of people to nuclear reactor radiation. In proceeding of:
The Second International Conference on Radiation and Dosimetry in Various Fields of Research (RAD 2014), University of NiΕ‘ (Serbia), 27-30.05.2014, pp. 49-52 (http://rad2014.elfak.rs/php/download2.php?file=../prezentacije/Proceedings%20RAD%202014.pdf). Sasaki, MS. (2003). Chromosomal biodosimetry by unfolding a mixed Poisson distribution: a generalized model.
Int. J. Radiat. Biol. , 79(2), p. 83-97. Sivia D.S., Skilling J. (2006).
Data Analysis. A Bayesian Tutorial (second edition). Oxford University Press. SzΕuiΕska M., Edwards A.A., Lloyd D.C. (2005). Statistical methods for biological dosimetry.
Health Protection Agency / Public Health England report no. HPA-RPD-011