The Predictive Distribution in Decision Theory " A Case Study

In the classical decision theory framework, the loss is a function of the decision taken 
and the state of nature as represented by a parameter θ . Information about θ can be obtained 
via observation of a random variable X . In some situations however the loss will depend not 
directly on θ but on the observed value of another random variable Y whose distribution depends 
on θ . This adds an extra layer to the decision problem, and may lead to a wider choice of actions. 
In particular there are now two sample sizes to choose, for X and for Y , leading to a range of 
behaviours in the Bayes risk. We illustrate this with a problem arising from the cleanup of sites 
contaminated with radioactive waste. We also discuss some computational approaches.


Introduction
Consider the following consulting problem: the client is involved in the cleanup of sites contaminated with radioactive waste, which involves sending bins of radioac- tive material to a nuclear reprocessing plant.Two such plants are available, one of which is less expensive but which will only accept a bin if the level of radioactivity of the material is below a threshold level.The actual level of radioactivity in the bin is determined by sampling at the reprocessing plant; as far as the client is aware only one such sample is taken.If the measured level exceeds the threshold, the ma- terial is returned to the client and must then be sent to the second, more expensive reprocessing plant which will accept material at any level of contamination.
The client wishes to base the decision of which reprocessing plant to use on a sample or samples taken from each bin on site before the material is dispatched.The cost of sampling is small relative to the difference in reprocessing costs, but is not negligible.How many samples should be taken, and how should this information be used?
The Loss Function in this problem clearly has the form L(al y) IS ify<c S / K + P otherwise L(a2,y) S + K Here al is the decision to send to the less expensive reprocessing plant at a cost S, K is the extra cost of sending to the other plant, and P the extra penalty incurred by sending first to the less expensive plant but having the material rejected.This occurs if the sample taken at the plant has an observed value y greater than the threshold level c.Since the value of S plays no part in the decision, we can without loss of generality take S = 0.
The client expects that there will be considerable variation in levels of radioac- tivity of the material within each bin, much of it being at quite a low level but with some "hotspots" of high radioactivity, suggesting a highly skewed distribution.It seems reasonable then to model the result for a single sample as a random variable X with an Exponential distribution.We write X Exp(0) to denote that X has density Here parametrises the "state of nature" and relates to the average level of ra- dioactivity of the material in the bin.Note however that the mean level is 1/.
The alternative parametrisation of the Exponential distribution is more intuitive, but we keep the present form for mathematical convenience.Because of the high skewness it is intuitively clear that a single sample will not provide reliable infor- mation about the average radioactivity level.It may be advantageous for the client to persuade the reprocessor, through financial inducement or otherwise, to take further samples before accepting or rejecting a bin.Thus there may be two sample sizes to consider, relating to the sampling before and after dispatch.Intuitively one might expect that increasing either sample size would be advantageous to the client, and that for a given cost of sampling the total sample size might be shared equally between the two stages.This turns out not to be the case.
In Section 2 we develop and analyse a Bayesian framework for this problem using a conjugate prior distribution for .I n Section 3 we consider the determination of optimal sample sizes at each stage of sampling.In Section 4 we reconsider the choice of prior and discuss some numerical strategies for incorporating a non-conjugate prior.
The basic principles of statistical decision theory, as used here, are described by DeGroot [2], although our notation is closer to that of Ferguson  [3].In the classical approach the loss incurred by the decision maker is a function of the action taken and the true value of an unknown parameter, information about which can be obtained by sampling.The situation in which the loss depends not on a parameter but on future observations was considered by Roberts [7] in the context of statistical prediction.Aitchison and Dunsmore [1] and Geisser [4] provide an overview and many applications of the predictive approach, some of which involve decision making but not the sample size determination problem considered here.A related problem in determining a single sample size in the classical framework when there are two "adversaries" with different priors, was considered by Lindley and Singpurwalla [5], and an application in environmental monitoring of radiation levels given by Wolfson et al. [8].

Bayesian Analysis
We assume that the uncertainty about/9 can be expressed as a Gamma distribution, F(a, A) with prior density 0(0) r() (2) If we wish to use a non-informative prior we might consider a 1 and A --+ 0, although this may not be sensible here (see Section 4).The main advantage of the Gamma prior is that it is a conjugate for the Exponential distribution, so that the posterior distribution for/9 after observing one or more sample values X will also be Gamma.Specifically, for a single observation X x the joint distribution of (X, 8) has density ,o > o (3) so by inspection the posterior density rOlX(8 x) c Oe -(x+=) and ( z) r(a + I, A + z).
The usual procedure, when the loss is a function of O, would be to choose the action (al or a2) which minimizes the expected loss under this posterior distribution, giving the Bayes Rule

5(x) argmin {Ex [L(a, t9 x)]} (4)
Here however the loss depends not on 0 but on the observed value of a second random variable, say Y, representing the result of the sample taken at the repro- cessing plant.The Bayes criterion is now Ex [L(a, y x)], the expected loss under the predictive distribution for Y given X x.If we assume that both the client and the reprocessor use the same sampling and measurement technique, then Y has the same distribution as X (for a given 0), in which case the predictive distribution The expected loss for a2, the expensive reprocessor, is fixed at S + K and for a Ex [L(al, y lx)] S + (K + x<=l-

K+P
The expected loss incurred by using this decision rule, say (x), can now be found by integrating the expected loss at fixed X, as given above, with respect to the marginal distribution fx of X.This gives the Bayes Risk B(r, ) of the rule with respect to the prior distribution r.Formally we may write to show two different ways of calculating the Bayes Risk corresponding to two different forms of iterated expectation.It is more convenient here to use the latter.
The marginal distribution for X has density x(z) Note that A is a scale parameter for the marginal distribution of X (and of Y) so that the problem is invariant to transformations (A, c) -> (kA, kc) for k > 0. This transformation corresponds to a change in the unit of measurement of radioactivity.
Similarly the decision made depends only on the ratio K/P not on the individual values.
Suppose then that with suitable units we take a 3, A 10, c 5, K = 10, P 15.The prior distribution for the mean level of radioactivity 1/0 is shown in Figure 1, and the marginal distribution for X (and Y) in Figure 2. We find that the Bayes Rule is ai if x < 9.423 5(x) (11) a otherwise It is clear from Figure 2 that at will be chosen most (93%) of the time, even though the mean level of radioactivity is often above the critical level c.This occurs because of the extreme skewness of the sampling distribution for Y 0 which makes the "gamble" of using the cheaper reprocessor worthwhile even when the average level of radioactivity in a bin is quite high.In the next section we consider the changes which occur when repeated sampling is used at both ends of the process (i.e. for X and for Y).Using the framework established in the previous section, we now suppose that the client bases his decision on samples Xi, X2,..., X taken from the bin, and assume that these are iid Exp().It is convenient now to work with the total X X1 / X2 /... / Xn which is a sufficient statistic for and is distributed as F(n, ).Similarly the total Y of m samples taken by the cheaper reprocessor, assuming the bin is sent there, will be F(m, ).The decision to accept or reject the bin is then based on the mean of the m samples, so that in the Loss Function c is replaced by mc.
Proceeding as before, we now find that the posterior distribution for X is F(n / c, A / x).The predictive distribution for Y X, following the method of Equations ( 5)-( 7), then has density r(m + n + ) () + x)n+ay m-1 ]rl(Y ]) r(m)r(n + ) (: + z + y)++ If we make the substitution u _.._e__ we find that the predictive density for A+x+y U r given X = x, has the form of a Beta distribution, and we write: for the predictive distribution of the transformed variable.
A closed form for the predictive probability PlY < mc X x] is not possible, but the incomplete Beta distribution is easy to calculate numerically (see is the Pg quantile of Be(m, n + ).
To determine the Bayes Risk B(,) for fixed m and n we use the marginal density of X. Proceeding as before we find that Be(n, a) (17) so from Equation (9) B(r, ) KP[X > f] + (g + P)P[X < f and Y > mc] = K (1-IBx, (n, x)) + (K + P)fo' (1-IB/x (m,n +Cz))fx(x)dx Although numerical integration is now required this can be accomplished quite easily using standard routines (see Press et al. [6]), and evaluation over a range of values of m and n gives a criterion for choosing the optimal (from the client's viewpoint) sampling plan.Using the parameter values from Section 2, Table 1 gives the Bayes Risk for values of m and n ranging from 1 to 6; note that these values do not include the cost of obtaining the samples.The optimal choices for m and n will depend on the sampling cost: if for example each sample determination has a cost of 0.1, we add 0.1 (n + m) to each value in the table and find that the optimum is n m 5.The advantage of not including the sampling cost explicitly in the table is that we can observe the behaviour of the Bayes Risk as n and m are varied.Notice that for n 1 the Bayes Risk initially increases as m increases from 1 to 2, the increased accuracy of determination by the reprocessor being disadvantageous to the client, but thereafter an increase in m results in a lower expected cost.For n however a higher value will always decrease the Bayes Risk, as one would expect: more information for the client should always result in a better decision.Table 2 shows the cutoff point for the Bayes lule, expressed in relation to the sample mean, i.e. if x/n is greater than the tabulated value, the bin should be G.JONES sent to the expensive reprocessor.As the number of samples increases, the cutoff converges quite quickly to the critical value c.Table 3 shows the proportion of bins which would be sent to the expensive reprocessor for each sampling plan.
Several different behaviours are possible, depending on the parameter values.In some cases the Bayes risk may decrease very slowly, or even increase, as m increases from 1; in other cases it decreases quite markedly.It is important therefore to get accurate information about costs and the prior before deciding whether it is worthwhile obtaining extra samples, and whether the extra effort should be devoted to X or Y or both equally.

Non-conjugate Prior
The prior Gamma distribution for employed in Sections 2 and 3 was chosen mainly for mathematical convenience.We now re-examine its appropriateness and how a wider class of priors might be incorporated into the analysis.
To obtain a reasonable prior distribution from the client, he must be invited to speculate on the likelihood of a range of values of .T his may be difficult since 8 is itself not a particularly meaningful parameter.A far more natural parametrization of the problem would be to use 1/ which represents the mean level of radioactivity in the bin; this is something about which the client might reasonably be expected to speculate.We could still proceed by showing the client graphs of the density of 1/ for various choices of a and A, as in Figure 1, but even so we are restricting ourselves to a class of distributions, the Inverted Gamma, which might be thought inappropriate.These distributions are very long-tailed, having less than a-1 finite moments.
Suppose that instead we decide to use a general prior distribution specified for 1/.We can still denote the implied prior for by r0 (), but the integrals needed to evaluate the marginal distribution for X and the predictive distribution for Y will not now involve simple special functions like the Gamma and Beta.It has become commonplace is such situations to employ some form of Monte Carlo integration.
There are essentially three stages to the calculation: Evaluate the risk for fixed cutoff and fixed sample sizes n,m.
Choose to minimize the risk for fixed n,m.
Choose n,m to minimize the Bayes Risk.
If we denote the rule with cutoff by , i.e. al if x < a2 otherwise (18) then we need to evaluate R(r, 5) KP[X > ] + (K + P)P[X < and Y > mc] (19) One approach would be to sample from the joint distribution of (, X, Y).Pro- vided that the prior for 1/ is reasonably easy to simulate, we invert a randomly drawn value to get , then draw X and Y from their conditional distributions F(n, ) and F(m, ) respectively.The conditional independence of X and Y given means that we do not need iterative methods such as the Gibbs sampler.Given a sample (, X, Y), i 1,... ,N we can approximate the risk for given by N R(r, ) _ g2:{;r)<} + (]C + P){>> and :>>$]} (20 where 27 is the indicator function. Because of the dimensionality problem this method requires a huge sample size to achieve even reasonable accuracy, and repeated computation for varying becomes very inefficient.A better approach is to simulate for only and to calculate directly Pe[X > ] and Pe[Y > mc].These probabilities are incomplete gamma functions and can be calculated quite efficiently (see Press et al. [6]), giving N R(r, ) _ g (1 IG,(n)) + (g + P)IGe, (n)(1 IGm, (m)) (21) where IG..(k) denoted the incomplete gamma function IG (k)  uk-e-Udu Note that only the X probabilities depend on (, so the Y probabilities for each 0i may be stored end re-used.This method requires a much smaller sample of 0 values to achieve reasonable accuracy, and is therefore more eNcient han use of the full multivariate joint distribution. Using the prior and parameter values from Section it was found ghat N 10,000 gave sufficient accuracy (2 dp) and a reasonable computation time (about 90s).Now however we are no longer restricted go a small class of priors.The calculagion was repeated using a r(4,1) prior for 1/0.This is shown in Figure for comparison with Figure 1; it is similar but much less long-tailed.The estimated Bayes Risk and cutoff value using this prior axe given in Tables 4 and 5.We now find that with a cost of sampling of 0.1 the best option is n = m = 1.

Discussion
In the above analysis we have considered the problem only from the client's point of view, assuming that he can pay the reprocessor to take extra samples, as well as deciding to take more samples himself, if this is to his advantage.We have also assumed that the critical value c used by the reprocessor is kept fixed for different sample sizes.If we now consider the reprocessor's point of view, he clearly does not want to accept material which has too high a level of radioactivity.We assume that he requires the mean level for each bin to be less than c, but that he does not allow for sampling variability in making his test.Were he to do so, he would want to adjust the critical value depending on the sample size m.It would also be to his advantage to take extra samples.Rather than use decision theory merely to improve the decision-making of one side in the process, as we have done here, it would be more appropriate to use an agreed decision theory framework as a negotiating tool in establishing an optimal sampling scheme which would be of benefit to both parties.
It has already been noted that the optimal solution for the problem we have considered seems to be quite sensitive to the parameter values and prior information.This shows the need for estimated costs to be as accurate as possible, and for prior data to be incorporated in choosing the prior distribution, possibly through an empirical Bayes approach.This sensitivity is probably due in part to the use of long-tailed distributions.There is a considerable range of behavior for different parameter values and distributions.In particular the fact that increasing the sample size for Y may either increase or decrease the risk is interesting, and is the subject of further investigation.

Figure 1 .Figure .
Figure 1.Prior inverted gamma density for mean level of radioactivity (1/8) with a 3 and A=IO