The asymptotic method of dead-time correction in poissonian distribution

NUCLEAR

INSTRUMENTS

AND

METHODS

78 (t97o) 70-76; © N O R T H - H O L L A N D

PUBLISHING

CO.

THE ASYMPTOTIC M E T H O D OF DEAD-TIME C O R R E C T I O N IN P O I S S O N I A N DISTRIBUTION* F. C A R L O N I t , A. C O R B E R I , M. M A R S E G U E R R A

a n d C. M. P O R C E D D U

C.N.E.N., Centro di Calcolo, Bologna, Italy Received 12 A u g u s t 1969 T h e p r o b l e m o f correcting for dead-time losses a n experimental c o u n t i n g distribution m e a s u r e d with a nonparalyzable device exposed to a poissonian source of particles is here discussed. An a p p r o x i m a t e analytical m e t h o d for estimating the m e a n freq u e n c y of the source a n d the dead-time of the system is given.

The involved statistical and systematic errors are detailed and graphs allowing their rapid determination or the proper design of an experiment are reported. The correctness of the proposed method has been successfully tested against stochastic simulations and experiments.

1. Introduction

mean frequency 2 related to the particle emission rate by the overall interaction efficiency. The distribution is obtained by an experiment consisting of a large number N,,~ 105-107 of successively repeated counting measurements of known length T ~ 10 -3 sec called "gate time" or briefly "gate". During the experiment a recording equipment, usually a multichannel analyzer, accumulates in the k-th channel the number of gates with k counts. Actually, after each gate the system remains inhibited for a known fixed "transfer time r " during which the number of observed counts is properly registered. Nevertheless, being the poissonian distribution homogeneous in time, these transfer times may be suppressed and the measurement considered as it were of total length NT instead of N(T+ z). In this way a sample from the population of the detected events (counts) observable in an interval T is withdrawn. The sample and the population will be respectively called "observed" and "observable" distributions. In an ideal counting device with infinite resolving time, each interaction event is detected and therefore the probability of recording k counts in a gate T is still a Poisson distribution nk(m) whose mean value m = 2 T is estimated by the first moment of the measured distribution. Thus 2 may be obtained being T exactly known. In a real arrangement, because of the finite resolving time, some interactions are not detected and therefore the measured pulse frequency is less than 2. Furthermore, being the lost pulses correlated to the observed ones, the observable distribution is no more poissonian. In order to estimate, in this case too, the true frequency 2, we assume the losses as mainly due to the dead-time of the counting apparatus, which is supposed to behave essentially like a non paralyzable device. That is to say, we postulate the system to remain blocked for a

It is a well known fact that from measured distributions of events, parameters of great physical interest may be inferred. A typical example of such a situation occurs in reactor physics where, for instance, from the observation of the statistical fluctuation of the neutron counts in a fixed time interval, the prompt neutron decay constant may be determined. In this type of experiment, the required distribution may be directly obtained at the output of a slightly modified multichannel-analyzer fed with counter pulses. However, in order to achieve a great accuracy in a short time, the counting rate must reach levels at which the dead-time losses begin to be relevant. Therefore the measured distribution is distorted and must be corrected, otherwise the precision claimed is only apparent. For example the difference between the measured frequency of a long-lived radioactive source and the true one, may quite often exceed by two or three decades the attainable accuracy. This paper deals with such corrections taking as a model a poissonian source emitting particles detected by a non-paralyzable counter. An analytical method for estimating the mean frequency of the source and the deadtime of the system is given together with a detailed analysis of the involved statistical and systematic errors. The correctness of the proposed method has been tested with a stochastic simulation of the chosen model. In addition to this internal check, experimental results showing good agreement between the expected and the observed distributions, are reported. 2. Problem's formulation Let us consider an experimental counting device suitable for distribution measurements, exposed to a poissonian source of particles. The counter sensitive volume is then a source of poissonian distributed 1) interaction events - not necessarily all counted - with

* Work performed u n d e r C o n t r a c t Reattori Veloci, C . N . E . N . t G u e s t researcher.

70

PBS-ADR

Programma

DEAD-TIME CORRECTION

IN P O I S S O N I A N D I S T R I B U T I O N

fixed unknown time 0 (dead time) following the detection of a particle, no matter how many unrecorded interaction events occurred within this time 0. Even with such a simple model, which adds, however, one more parameter (namely 0) to be estimated, the problem becomes much more complicated compared with the above case of infinite resolving time. Nevertheless, proceeding in a similar way, we shall first derive the expression of the observable distribution and then obtain estimates for both 2 and 0 by the method of moments. 3. Rigorous expressions of the moments To determine the distribution we shall proceed heuristically, referring to 2), for a more general treatment. Starting from an arbitrary time origin, when the system is unblocked, let pk+l(t)dt be the probability that the (k + l)-th pulse will be observed between t and t + dt. Its expression may be obtained considering the first k pulses and the last one separately. If this has to be observed between t and t + d t , the system must be unblocked at that time. Therefore the probability of observing the event is simply 2dt and the one observed immediately before should have arrived at the time t - O at the most. Each of the k events observed in (0, t-O), paralyzes the system for a time 0. Due to the complete lack of memory of the homogeneous poissonian processes, these dead-times may all be lumped together in a subinterval of kO < t length during which the system is inhibited. The observed k pulses may thus be treated ignoring the dead-time, as if they had arrived during the remaining portion t - k O of the interval, thus obeying a Poisson distribution with mean value 2 ( t - kO). Therefore:

71

vi(T) = L k'Pk(T) = L k i[ Q k ( T ) - Qk +l(T)] = k=O

k=O

= k=O Z [(k+

=

kM

= Z [(k+l)i-ki]]s(k+l,m-kz)/k!,

(5)

k=0

where k~a is the maximum integral number contained in T/O. The above expressions implicitly define 2 and O in terms of the population moments vi. If a system made of any pair of these equations - e.g. the first two of them - could be solved, rigorous explicit expressions for 2 and 0 would be found. Then equating the involved population moments to the correspondent (experimentally obtained) sample moments, the numerical values of the required parameters could be estimated with an accuracy only dependent on the statistical precision of the measurement, i.e. on the sample size. However, due to the lack of an analytical solution and to the limited convenience of a numerical approach, we have to resort to more practical approximate expressions, therefore introducing a systematic error.

where U(t) is the usual step-function. The probability Qk+ 1(T) of observing k + 1 events or more in (0, T) is:

4. Asymptotic expressions of the moments In the following paragraphs we shall see that the so called 2, 4) "asymptotic" expressions of the moments of the observable distribution - rigorous only for T-+ oo are well suited for the estimation of Z and 0 from the measured sample, giving simple and accurate expressions for both parameters even for gates of the order of few dead-times. To show this we shall give the first two eqs. (5) a different form by taking first their Laplace transform and by performing then their inverses using the residue method. The L-transformed i-th moment is

Qk+ 1 ( T ) =

vi(s)

Pk+ l(t) = U(t-kO)27~k[2(t-kO)],

(1)

Pk+ l ( t ) d / = •

0

(2) where z = 20 and 7(J, t) is the incomplete gamma function defined as 3) = V(T-kO)~(k+l,m-kz)lk!,

),(j,t)=

e - ~ x J - a d x = ( j - l)t 0

~j_ ~(x)dx.

(3)

o

Finally the required probability Pk(T) of observing exactly k events in (0, T) is

Pk( r) = Qk( T ) - Qk +, ( r).

(4)

The i-th moment (i = 1, 2 . . . . ) of this distribution is then

--_

(J.e-S°] k

1 2 L [-(k+l)/-k'] s 2 + s k=o t )-S~"~"~"~"~-s ! '

( i = 1 , 2 .... ),

(6)

from which we have

",(s)--- v(~) -- ,~{sb + ;~(i - e-~°)]} - ' (7) This function has a double pole at s = 0 and a denumerable infinite set of complex conjugate simple poles Sk = (e~__+iflk)/0, (k-----1, 2 . . . . ) lying in the halfplane Re(s) <0. It may be shown that/7 k is the solution of

z + ln z = ln (

flk ) f l k sin/4k

tg flk ' with ( 2 k - 1 ) ~ < ~ k < 2 k ~ ,

(8)

72

F. CARLONI et al.

and fig ~k = z + tgfik

(9)

whose absolute values increase with k. The first five solutions o f eq. (8) have been numerically found for various values of z and are shown in table I together with the c o r r e s p o n d e n t ~k. The inverse t r a n s f o r m o f eq. (7) is given as the sum of an a s y m p t o t i c p a r t due to the residue at the d o u b l e pole v.~(r)=~

+1+z-2\1+20!

+1+20'

(10)

plus a transient p a r t due to the sum of the residues at the simple poles

Vtr(T) =,,_Z, Aksinp,, - + B;cos ;,

e

(11)

where

form coincide with those of the first m o m e n t . Thus, m a k i n g the same kind o f a p p r o x i m a t i o n s , we m a y write

v2(T )

(l+z)3

I- l ~ z

fik(l + Z +2~k) 2

Bk = 2z

2

2

Z

+ 0~

+ \1~!

" (14)

The a b o v e formulas are the before m e n t i o n e d asymptotic expressions for the observable d i s t r i b u t i o n ' s moments. The m o s t i m p o r t a n t results here o b t a i n e d concerns the rapidity of the convergence o f the rigorous expressions t o w a r d s their a s y m p t o t i c values.

5. Systematic and statistical errors in the parameters' estimation The required a p p r o x i m a t e " a s y m p t o t i c " expressions "~as a n d 0as of the p a r a m e t e r s 2 and 0 are easily obtained from eqs. (13) and (14) in terms o f the p o p u l a tion m o m e n t s as

v

Ak=ZZ

(1+).0) ~

,,'

O5)

( k = 1,2 .... ), 2

+ &2]

'

2

ak--fia +ak( 1 +Z) (~# + ,q2) [(I + z + ~k)2 -~- fl~]

(12)

F r o m the d a t a of table 1 we see t h a t the transient p a r t dies-away in few d e a d - t i m e units, so t h a t in practical cases for T > 100, only the a s y m p t o t i c p a r t needs to be considered. M o r e o v e r in this case the c o n s t a n t term in eq. (10) is negligible and we have v(T)

m 2T 1+z - 1+20

(13)

As far as the second m o m e n t is concerned, it is easily recognizable t h a t the poles of its Laplace transTABLE 1 First pairs xA:, fi~. (k = 1, 2, ..., 5) as obtained by eqs. (8) and (9) of the text.

where P2 = v 2 - v 2 is the second central m o m e n t (variance) of the observable distribution. F r o m the above f o r m u l a t i o n of the p r o b l e m , it identically follows /~a = ~ . a ~ T = v

/(~2),

(17)

As a consequence o f the use of practical but a p p r o x imate expressions (13) and (14) instead o f the rigorous ones (5), a systematic error on 2 and 0 is made. To evaluate it, numerical values of m, z a n d T - i.e. of 2 and 0 - are assigned within the d o m a i n of interest. The c o r r e s p o n d e n t p o p u l a t i o n m o m e n t s v and v2 are then exactly calculated by eqs. (5) and introduced into eqs. (15) and (16) to get 2~ and 0-~s. The differences 2 - 2 , ~

z

10-1

10-~

10-:~

~1 ~2 ~3 ~4 ~5

-4.0158 - 4.7500 - 5.1799 -- 5.4784 - 5.7099

- 6.6274 - 7.1392 -- 7.5206 - 7.8043 - 8.0277

- 9.1941 - 9.5471 - 9.8751 -- 10.1374 - 10.3504

fil fiz fi3 f14 fi5

3.9286 10.5815 16.9883 23.3354 29.6582

3.6451 10.3944 16.8597 23.2383 29.5804

3.5059 10.2455 16.7460 23.1492 29.5078

and 0 - ~ , ~ are therefore obtained. A n o t h e r systematic error results from the assumed independence o f successive intervals. In the analytical t r e a t m e n t the system was assumed initially u n b l o c k e d , that is no pulse arrived during the time 0 preceding the beginning o f the gate time. However, this is not strictly true, d e p e n d i n g on the frequency 2 a n d on the r a t i o T/O. W i t h o u t considering this p o i n t in detail, we assume that for 2 0 < 1 and T/O>50, the involved error is negligible*. * Actually the gate is not too crowded with registered pulses if vas0~ T. Making use ofeq. (10) this condition gives: T/O>/~½(20)2.

DEAD-TIME

CORRECTION

IN

POISSONIAN

73

DISTRIBUTION

So far we m e n t i o n e d the systematic errors only. H o w ever, the statistical ones need also to be considered because f r o m the m e a s u r e d distribution, instead o f the required p o p u l a t i o n m o m e n t s v~ a n d / ~ we actually get the c o r r e s p o n d e n t sample m o m e n t s n~ and m v Introducing these values into eqs. (15) a n d (16), as prescribed by lhe m e t h o d o f momentsS), we o b t a i n the estimates

12:51 sloellaslle simulation

result

experimental

,

<~

= ~

,

O~s = T E l

053

-- J ( ~ ) ] .

(16')

Since 2,s a n d 0,~ are sample characteristics based u p o n s a m p l e m o m e n t s and being the sample size N large, we shall assume 6) t h a t 2a~ and 0a~ are n o r m a l l y d i s t r i b u t e d a r o u n d their expected values )'a~ and 0,~. Their s t a n d a r d errors are given by the well k n o w n expression 7) for the variance o f a function of sample moments. F o r 2,s(n, m2) w e have:

~,

,

10

+ l

+

\Tn-n / o

{~-m~ l Itt,(n,m2)+O(N ~),

(19)

2~0

w h e r e D2(~) is the variance o f the r a n d o m variable ¢,

#~ i n , rn2) is the covariance between n and m2 and the

!~

3

~

T

T

[

T

[

i I : r

i

Ii

~

T

stochastic

~ T

U ~

simulation

_ _ _ ~

L ~ ]

I I 1 [ ] _ _ 1

16 3

I

L ~

Ill]

n~ I

11~2

Fig. 2. Relative error of the estimate 0~s, multiplied by ~/N, vs the adimensional parameter zas = 2as0as. Experimental and simulated results are also reported. The standard error of each abscissa is obtained from the variance of the set of the Zas values; the standard error of each ordinate is calculated assuming a normal distribution.

result

experimental

subscript zero indicates t h a t the derivatives are calculated in the p o i n t n = v and m2 = #2. Substituting the u n k n o w n p o p u l a t i o n m o m e n t s v, /~i by their consistent estimates n, ml a n d neglecting the terms o f the o r d e r o f N -~-, the s t a n d a r d error reads

2~

A

D(<)

1 r i-

I_ 2Tx/N L

V3

6ms] i,

(20)

V2

t:

E /

0

[i .i

.

J

~

LL!

' [ ~

~__

i

I

' ~._I

I I I I

t

10

IIlas

F i g . 1. R e l a t i v e e r r o r of" t h e e s t i m a t e

2as,

m u l t i p l i e d b y ~ / N , vs

the adimensional parameter rn~ = ~-a~T. Experimental and simulated results are also reported. The standard error of each abscissa is obtained from the variance of the set of the mas values and it results to be negligible. The standard error of each ordinate is calculated assuming a normal distribution.

where V is the s a m p l e ' s relative variance defined as the ratio between the absolute variance m 2 a n d the m e a n value n. M a k i n g use o f eqs. (5) to get the p o p u l a t i o n m o m e n t s , the s t a n d a r d error o f 2,s was calculated for various values o f m a n d z within the d o m a i n of interest. T h e results are shown in fig. 1 where the relative e r r o r R(mas)=R(2as)=D(,~as)/~as multiplied by ~ / N is p l o t t e d vs ~.s- Since the curve c o r r e s p o n d e n t to the various Lras values practically coincide, at least in the

74

F. CARLONI et al.

T~as

mas

17

10

i

I '°

I

__ 10 7

_ 10~

10g

_10 3

, i 105

~108

practically coincide, at least in the examined range, it turns out that the relative error made in the determination of an u n k n o w n dead-time does not essentially depend upon the gate time. To reduce the total measurement time, short gates are therefore preferable, provided the aforementioned condition T > 50 0 is satisfied, By increasing the number N of examined gates the standard error of 2a~ and 0,~ could be in principle reduced at will*. Such an error should nevertheless always be kept much greater than the correspondent systematic one, otherwise the claimed precision would only be apparent. The connection between these two kinds of error is here established through their ratio -

~--1

_10 7

16~1 _10 s

b

\\

I

b

I

=10 3

10 ~

iF rb

-lO"

I

iII L-- I

Fig, 3, Nomogram for the determination of the ratios between the statistical and the systematic errors of ).,s and 0~s, for given m~ and z~ values. The gate length, in unit of dead-time, can also be read. The dotted line refers to the numerical example of the text. examined range, it turns out that the relative error made in the determination of an u n k n o w n frequency is essentially independent of the counting system's deadtime, i.e. independent of the distorting cause size. F o r O,~(n, m2),following an identical procedure, we have:

D(Oas) N

T._.~F m 4 - m 2 -.[n~,~/N L 4rn 2 + V ( 1 - - zax/V)2 +

m~ ) (l-~x/V) x/(nm2

D().as) D(m.,)

(22)

D(0,s)

(23)

_

and

__10 4

r

_

.

(21)

The relative error R(O,~) multiplied by ~ / N is plotted in fig. 2 vs ~,,~. Since the curves for various ~ . s values

p(0as ) --

10.,-01

These ratios multiplied b y ~ / N were calculated for various m and z values utilizing eqs. (5), (15), (16), (20) and (21). The results are shown in the n o m o g r a m of fig. 3. For given mas and £'as values, the n o m o g r a m allows a rapid evaluation of both functions together with the correspondent value of T/O,s. The values T/O,~ < 50 are merely included for completeness. F r o m the n o m o g r a m we see that in practical instances the p~/N values are so large that even for p > 10 large values of N are still possible. Therefore the standard error may be quite small while the systematic error is negligible in comparison. In conclusion ).,~ and 0,~ are good estimates of the required parameters 2 and 0. All the mentioned computations were performed with an IBM 7094 computer. The codes were always written in double precision in order to avoid numerical instabilities, especially calculating the moments v i by eqs. (5).

6. Numerical examples We shall now illustrate the practical utilization of the reported formulas and graphs. A n experiment lasting 33 sec and consisting of N = l0 s gates of length T = 300 /~sec, and transfer time • = 30 /~sec was performed. F r o m the observed distribution the moments n = 2.923 and m2 = 2.78 were calculated. Eqs. (17) and (18) then yielded the estimates m , ~ = 2 . 9 9 6 and % s = 2 . 5 1 x 1 0 -2. Their standard errors are evaluated from figs. 1 and 2 tentatively • Electronic drifts and other instabilities actually limit the feasible number N of gates.

DEAD-TIME

CORRECTION

IN P O I S S O N I A N

I DETECTORF

! DELAY

CLOCK t - - ~

CHANNEL ANALYZER MULTI

l V'ULb~.-"1 GATE STOP ~CIRCUIT PULSE /]

k

Direct Clock

pulses

k

L

Delayed

I I ! I

I 1

I I

I

',11

Analyzed detector

pulses

75

DISTRIBUTION

- T

X

I

I

[1# °

I

T I

I

Fig. 4. Set-up and time scale of the experiment. assuming that - being N very large - m,s and Za~ may replace their correspondent expected values ~a~ and ~as without the plotted functions changing appreciably. We find ~/(N)R(2,~) = 0.91, v/(N) R(Oa~)= 29, so that* 2,s = 9986 + 2 9 sec- 1, 0a~ = 2.5-t-0.2 sec. (24) The smallness o f the standard errors confirms the feasibility of replacing ~a~ and ~,s with the sample values m,s and z,~. It therefore confirms that N is actually a large number. However, N could be so large to make the systematic error no longer negligible relatively to the standard error. F r o m fig. 3, joining with a straight line the previously obtained m ~ and Zas values - since T/Oas~120 we calculate p(Aa~)> 1000 and p(O,~),,~20. Therefore we conclude that the systematic errors are actually negligible and the values (24) are g o o d estimates of the required parameters 2 and 0. Let us now assume instead the same experimental m o m e n t s of the previous example as obtained from N = 107 gates in 55 rain. The above quoted standard errors must then be reduced by a factor of ten and from fig. 3 it turns out that p(0,s),,~ 2. This is too small, while p(2.~) = 380 is acceptable. The error of 0.~ is therefore apparent, and must be increased, while the one of 2~s is still significant. * The standard errors could have alternatively been calculated from the first four experimental moments by means of eqs. (20) and (21).

7. M o n t e

Carlo and experimental

results

To test the validity of the above analytical treatment, a stochastic simulation of the assumed model was performed. Starting from given values of 2, 0, T, • and N, the poissonian distribution of the interaction events was generated assigning interarrival times t = - 2 -t lnr, where r is a p s e u d o r a n d o m n u m b e r in (0,1). Each event was properly accumulated in the "experimental" distribution only if it arrived during the gate time and was separated from the preceding one by at least a dead-time. After a N(T+ "c)time, when 2 N ( T + v) events have been processed, the "experiment" is over and from expressions (15') and (16') the estimates 2,s and 0,s of the estimands '~as and 0.~ are obtained. Since the estimate's standard error is a useful indicator of how close the estimate is to the u n k n o w n estimand, the correctness of expressions (20) and (21) was checked first. To do this three Monte Carlo experiments were repeated several times. The first was repeated forty times with N = l0 s, while the others were performed twenty times with N = 2 x 104. The mean values and the standard deviations of the set of the estimates 2.s and 0,s were then calculated. The standard deviations of "~/(N) R(Aas) and ~,/(N)R(O,~) were calculated assuming these quantities as normally distributed. The results, in good agreement with the analytical forecast, are shown in figs. 1 and 2. The closeness between the estimates ),,s and 0,~ and the true values 2 and 0 was successively checked through the same set of "experiments". Since these "experiments" were designed utilizing the n o m o g r a m o f fig. 3 so as to keep the systematic error negligible, a

76

F. CARLON1 et al.

gaussian sampling distribution around the true values should result. This was checked on the whole set of 100 available "experiments", and it turned out that 1. the 57%, 94% and 97% of the "~as values differed from 2 by less than once, twice and three times respectively the standard error D(2,s); 2. the 63%, 94% and 100% of the 0~ values differed from 0 by less than once, twice and three times respectively the standard error D(O,~). The positive results of the Monte Carlo tests confirm that the proposed analytical treatment satisfactorily describes the assumed model. Utilizing the experimental apparatus below described, a set of actual experiments was then performed to check the applicability of the model itself. Eleven distributions with N = 105 and one with N = 5 × 105, each repeated fifty times were measured. The relative errors R(A,s) and R(O,~) were calculated as before from the set of the ),,~ and 0a~ values. These experimental results are also reported in figs. 1 and 2 and their good agreement with the theoretical previsions confirms that the counting apparatus actually behaves like a non paralyzable device. The observed distributions were obtained by the experimental set up schematically shown in fig. 4. Pulses from the detector are continuously fed to the gate circuit and are registered or disregarded according to whether the gate is open or closed. The gate circuit's operation - i.e. the succession of the gate and transfer times - is controlled by the clock pulses. These are equally spaced according to a T + z interval and applied both directly and delayed of a ~ interval. Each direct clock pulse closes the gate circuit, thus starting the transfer time. The same pulse, when delayed, opens the gate circuit a time z later, thus ending the transfer time and starting the gate time. During the gate time T each detector pulse operates the channel stepping of the address register into the MCA. During the transfer time z a "plus one" is added to the contents of that channel whose address corresponds to the number of registered detector pulses. Then the address register is reset and the analyzer is ready for a new cycle. After N cycles, the measurement is over and the observed dis-

tribution is directly obtained dividing by N the contents of each channel.

List of Symbols N = Number of examined gates T, z = Gate and transfer-time = Mean frequency of the interaction events 0 = Experimental apparatus' dead-time ,7.a~, 0a~ = Expected "asymptotic" frequency and dead-time 2~, 0~s = Estimates of '~s and 0,s according to the method of moments pk(t)dt = Probability of the k-th count occurring between t and t + dt Qk(t), Pk(t) = Probability that at least and exactly k counts occur in (0, t) v~,/~i = Moment and central moment of order i of the Pk(t) distribution ni, m~ = Sample moment and central sample moment of order i V = Sample relative variance D(~), R(~) = standard and relative error of the sample characteristic p(~) = Absolute value of the ratio between the standard and the systematic error of the biased estimate 4. The authors wish to thank Mr. T. Martinelli for machine codes' programming, Mr. F. Lolli for performing numerical calculations, and Mr. G. Cantoni who assisted in running the experiments.

References l) L. Janossy, Theory and practice o f the evaluation o f measurements (Oxford, 1965) p. 78. 2) I. De Lotto, P. F. Manfredi and P. Principi, Energia Nucl. (Milan) 11, no. 10 (1964) 557. 3) A. Erdely, ed., High transcendental fimctions, vol. 2 (McGraw Hill, New York, 1953). 4) W. Feller, An introduction to probability theory and its applications, vol. 1, 3 r~l ed. (Wiley, New York, 1968) p. 339. 5) H. Cramer, Mathematical methods o f statistics (University of Princeton Press, 19.58) p. 497. 6) Ibid., p. 363. 7) Ibid., p. 354.