Checking for model failure and for prior-data conflict with the constrained multinomial model.

Multinomial models can be difficult to use when constraints are placed on the probabilities. An exact model checking procedure for such models is developed based on a uniform prior on the full multinomial model. For inference, a nonuniform prior can be used and a consistency theorem is proved concerning a check for prior-data conflict with the chosen prior. Applications are presented and a new elicitation methodology is developed for multinomial models with ordered probabilities.

Introduction

Suppose we have a sample of n from a multinomial distribution and let denote the cell counts. There are a number of applications where constraints are placed on such that it lies in some proper subset of the full k-dimensional simplex , see Sect. 2. This is the constrained multinomial model. Also, there are applications where there is additional information concerning that takes the form of a prior probability measure on . Our goal here is to check the model, namely, see if the data indicate that the true lies outside , and if the model is found acceptable, then check the prior, namely, see if the data indicate that the true lies in a region of relatively low prior probability in .

The primary way to determine whether or not an ingredient to a statistical analysis is appropriate is to compare it somehow with the data. For example, if the observed data is surprising for every distribution in the model, then it is reasonable to question the appropriateness of the model. By surprising it is meant that the data falls in a region such that each distribution in the model gives a relatively low probability for that data’s occurrence. Once a model is accepted, a chosen prior is questionable if there is an indication that the true value lies in the tails of the prior. For example, if prior measure with density is used for a model parameter with true value , then a prior is placing its mass in the wrong part of the parameter space if is small. While we cannot know this quantity, an assessment of this can be made based on the data, as discussed in Sect. 4, which is in itself an assessment of whether or not the observed data is incompatible with its prior distribution, namely, whether or not a prior-data conflict exists. As discussed in Al Labadi and Evans (2017b), inferences can be quite sensitive to the choice of the prior when prior-data conflict exists.

While a simultaneous check on the model and prior can be considered, as in Box (1980), separate checks are developed here so that we can identify whether an identified problem lies with the model or with the prior. Also, the check on the prior is performed only if the model has not failed its check, as the prior is implicitly dependent on the model. As discussed in Evans and Moshonov (2006, 2007), the check on the prior can sometimes be further decomposed so that individual components of the prior can be checked. So there is an overall logic for checking the ingredients to a Bayesian statistical analysis.

A notable aspect of the model checking approach taken here is that it is Bayesian but does not involve the prior . Rather we consider first the multinomial model with no constraints and place a uniform prior on . Then, based on a distance measure defined on , the prior distribution of the distance is compared with its posterior distribution. In essence, if the posterior has concentrated less about than the prior, then there is evidence against the constrained model and if the posterior has concentrated more about than the prior, then there is evidence in favor of the constrained model. There are several advantages possessed by this approach. First, up to simulation errors, the assessment is exact and does not depend upon asymptotics. Second, the approach works even when is a complicated lower dimensional subset of and it is straightforward to implement provided is relatively easy to compute. Of course, we have the freedom to choose d to address the latter concern. Finally, it is possible to obtain evidence in favor of the constrained multinomial and not just evidence against as is common with most goodness of fit tests. This makes sense as, provided there is i.i.d. sampling, then some multinomial distribution is correct. A check for i.i.d. sampling can be made by generalizations of the runs test, but we always assume that such sampling is correct here. Model checking is discussed in Sect. 3.

The check for the prior is based on the approach taken in Evans and Moshonov (2006) and is discussed in Sect. 4. The consistency of this check was established in Evans and Jang (2011a) under fairly general conditions. These conditions do not hold for the constrained multinomial, however, and so a consistency theorem is proved for this context. A new result is established for putting a prior on ordered probabilities which arise in many important applications. Computational issues are addressed in the examples. A high level of numerical accuracy is generally not required, however, as knowing relevant probabilities to one or two decimal places is sufficient to detect whether or not a problem exists.

Box (1980) is a key paper in Bayesian model checking. One notable aspect of this paper is that the check is based on the prior predictive distribution of the data, namely, the distribution of the data as induced by the sampling model and the prior distribution. There are also approaches to model checking that are based upon the posterior predictive distribution of the data such as Guttman (1967), Rubin (1984) and Gelman et al. (1996). A different approach is taken in Bayarri and Berger (2000) and Castellanos and Bayarri (2007) where the check is based on the conditional distribution of the data given a statistic that is asymptotically minimally sufficient. These checks are similar to frequentist model checks based on the conditional distribution of the data given a minimal sufficient, which is independent of the model parameter. As noted in Robins et al. (2000), there is a problem with posterior predictive checks when based on p-values as these fail to have to have asymptotically uniform distributions when the model is correct and so can be misleading. The Bayarri and Berger (2000) checks do not suffer from this defect. There have been attempts to recalibrate posterior p-values as in Hjort et al. (2006) and van Kollenburg et al. (2017). There is no reason to suppose, however, that model checks need to be a posteriori or that p-values need to be used as part of this. Posterior distributions are used in inferences about parameters as necessitated by the principle of conditional probability which says that prior beliefs must be updated via conditioning on the observed data. But this says nothing about how checks on the sampling model or the prior should be conducted. Also, there are many issues that lead to doubts about whether or not p-values are valid expressions of evidence. In the case of the constrained multinomial it is natural to put a uniform prior on the full multinomial model and then, as discussed in Sect. 3, assess the submodel based upon a valid measure of evidence.

Several different approaches have been developed for checking the prior. Bousquet (2008) bases such a check on a comparison with an improper reference prior. Presanis et al. (2013) develop diagnostics for such checks. Nott et al. (2020) develops an extension to the approach of Evans and Moshonov (2006) that provides a check with an invariance property necessary for contexts where the data is continuous. The data for the constrained multinomial is discrete so this extension is not required here. Tail probabilities are used for checking the prior, but these are used simply as a measure of surprise and not as a measure of evidence as there is no sense in which a prior is true or false, only that it may be inappropriate.

Examples

The following examples of constrained multinomial models are considered. The notation T in the examples refers to the vector of counts and subscripts on T reference the count of a particular cell. Also, the notation T(x) references the vector of observed counts based on a sample x. It is to be noted that T is a minimal sufficient statistic for the full multinomial model and so the full data x is not necessary for the approach taken here.

Example 1

Quantum state estimation.

A qubit is any binary quantum system. For example, the spin of the electron, which can be up or down, or the polarization of a photon, which can be polarized vertically or horizontally. Quantum theory allows for these to be held simultaneously, called superposition, and this is represented by a certain operator called the state, in turn represented by a Hermitian unit-trace matrix. To estimate the state of such a quantum system, in an experiment, a physical apparatus with detectors is employed that results in one and only one of the detectors recording a click and the i-th detector records a click with probability . So a single measurement is a value from a multinomial distribution. In an experiment, n identical copies of the quantum system are used and i.i.d. measurements are obtained. Therefore, with the observed count from the i-th detector denoted as , then multinomial.

Depending on how the measurements are taken, the probabilities must satisfy certain constraints. For example, when , then the symmetric trine model imposes the constraint , with more involved constraints required for the asymmetric case as discussed in Example 4. When , then the cross-hairs model imposes the constraints and , while the tetrahedron model corresponds to the constraint only. When , then the Pauli model imposes the constraints and . More on how these models arise can be found in Shang et al. (2013). These applications produce a rich variety of constrained multinomial models.

Example 2

Hardy-Weinberg equilibrium (HWE).

Consider a trait in a population governed by a gene with l alleles with population proportions . There are then possible genotypes . Let denote the proportion of individuals in the population of genotype for and note that and are indistinguishable. Then, under the assumptions of the Hardy-Weinberg equilibrium,So it is necessary to check if there is evidence for or against (1) holding and this will be based on the observed counts in a sample of n where denotes the number of observations having genotype .

Example 3

Multinomial with ordered probabilities.

Sometimes it is reasonable to suppose that the satisfy an ordering such asSuch a model arises in contexts where systems exhibit aging as in, for example, Briegel et al. (1994). In certain ROC analyses for a medical context, see Zhou et al. (2011), there is a diagnostic variable for a disease which leads to a classification of patients with known disease status into ordered categories. It may be natural to believe that the probability of being in the i-th category increases with i for diseased patients and decreases for nondiseased patients. In dose-response models, as discussed in Chuang-Stein and Agresti (1997), consider I different locations (or treatments) and J levels of exposure to a toxic substance from from low to high with the probability of a reaction for an individual at location i with exposure level j. A monotone increasing effect due to increasing exposure leads to for . For a constrained multinomial satisfying (2), refers to the count for the i-th cell.

The following result is useful with this model.

Lemma 1

Any satisfying (2) is given by, for some ,and any produces a satisfying (2) via (3), with

Proof

Suppose for . Then and it is clear that and so . It is also immediate that satisfies (2).

If satisfies (2), put and noteTherefore, for and . Since , then and if then which by (2) implies which is false. So for and similarly . Finally, .

A particular model satisfying (2) is given by the Zipf-Mandelbrot distribution where for parameters , for and and denoted here by . When this is the uniform distribution and for fixed this converges to the uniform as . For fixed the distribution becomes degenerate on the first cell as . For large k with and , the zeta distribution, with for , serves as an approximation. As discussed in Izsák (2006), there are a variety of applications of this distribution as in word frequency distributions in texts.

Checking the model

In the context of the constrained multinomial model it is possible to use more formal inference methods than are typically used for model checking. The following subsection describes these inferences in some generality and subsequently these are applied to the problem of interest.

Relative belief inferences

The approach to inference as discussed in Evans (2015) will be employed as this is based upon a direct measure of evidence obtained via the principal of evidence as described, for example, in Achinstein (2000). For example, suppose we have a probability model where is a finite set, that interest is in whether not the event has occurred and the event C, with , has been observed to be true. The principle of evidence then says that there is evidence in favor of A being true if , there is evidence against A being true if and no evidence either way if as then A and C are independent. So there can be evidence in favor of A being true if the belief that A is true has increased after being told the information that C is true and evidence against A being true if the belief that A is true has decreased after being told that C is true. If A and C correspond to specifying the value of functions U and V, say and then, letting and denote the marginal and conditional probability functions of U, there is evidence in favor of being true when or equivalently whenever the relative belief ratio . Note that the cutoff value of 1, for determining whether or not there is evidence in favor or against via the relative belief ratio, is not arbitrary, as it is dictated by the principle of evidence.

It would seem that the bigger is than 1, the more evidence there is in favor of , but this needs to be calibrated. In Evans (2015) this is done via quoting the strength of the evidence asFor when and (4) is large, there is a small probability that an alternative value of U has more evidence in its favor than so it is reasonable to regard the value as indicating strong evidence in favor of and when (4) is small, there is only weak evidence in favor of . Similarly, if and (4) is small, then there is strong evidence against , since there is a large probability that an alternative has either more evidence in its favor, or at least less evidence against it, and when (4) is large, there is only weak evidence against . Note that the relative belief ratio is used here only to order the possible values of U with respect to the evidence. Also, alternative measures of the strength of the evidence can be considered. For example, the conditional probability is certainly relevant, particularly when the number of possible values of U is quite small but in general contexts, (4) works quite well.

For an estimation problem, namely, where it is necessary to choose a value of U from the set of possibilities after observing , the natural estimate is to use the relative belief estimate as this maximizes the evidence in favor. The plausible region , namely, the set of values of u for which evidence in favor has been obtained, can be used to assess the accuracy of , via the size of as well as its conditional probability .

While the simple finite context just described gives the basic ideas underlying the approach to inference adopted here, this can be generalized to infinite discrete models in an obvious way and also to continuous models via taking limits of sets shrinking to points as discussed in Evans (2015). In Bayesian contexts with a proper prior, then there is a joint probability model for the model parameter and data x. Denoting the prior probability measure for by with density with respect to a relevant support measure, the prior probability measure for a marginal parameter by with density and the corresponding posterior quantities by , then the relative belief ratio for at is given by with the other quantities defined similarly and with the same interpretations. Various consistency and optimality results for these inferences are summarized in Evans (2015) and see Sect. 5.

Model checking via relative belief

The check on the constraints should not involve the prior since it is not involved in the production of the data. It is permissible, however, to use another prior on for the model checking step. This leads to formal inference methods for model checking using relative belief. The uniform prior on treats all multinomials equivalently and so this prior is used for the model checking step here. Certainly alternative proper priors could be used in the model check, such as the Dirichlet which is Jeffreys’ prior, but a good justification for this should be provided.

Suppose first that . The posterior probability of is where is the Dirichlet probability measure. The relative belief ratio of equals and the strength of this evidence, since there are only two possibilities, can be assessed via the posterior probability . So, if and is low, there is only weak evidence in favor of while, if and is low, there is strong evidence against . Similarly, if and is high, there is strong evidence in favor of while, if and is high, there is only weak evidence against .

It is often the case, however, that , because is of lower dimension, so and it is not possible to proceed as just described. If is very small, then again the preceding approach seems questionable as it cannot be expected that will be large, and so obtain strong evidence in favor of the model, without a large amount of data. So the approach discussed in Al Labadi et al. (2017a) and Al Labadi and Evans (2018) is used in such situations. A distance measure d on is specified so measures the distance is from which will be denoted hereafter as the set is implicit. In all the examples the distance used is where denotes the Kullback-Leibler divergence between the two multinomial distributions specified by and . Note that, see Pawitan (2013), the MLE of in the constrained model, which is the same as the relative belief estimate of , will converge to the point in which minimizes the KL divergence between and the true value, which lies in but may be outside of . As such, KL divergence seems like a natural statistical distance measure. Still alternative distance measures, such as squared Euclidean distance , can be used if it is felt that there is some relevance to a particular application. In general, the choice of any distance measure leads to a check on the model. In our experience, however, the results, in terms of whether evidence in favor of or against a model is obtained and the strength of that evidence, are not highly dependent on this choice.

The prior and posterior distributions of induce prior and posterior distributions for . If the posterior of concentrates more (less) about 0 than the prior, then this provides evidence in favor of (against) . The relative belief ratio given by measures this as there is evidence in favor of when and evidence against when . Typically the limit cannot be computed exactly so is selected satisfying . In practice, there is a such that, if , then can be regarded as true or as an acceptable approximation. The value of can be determined by bounding the absolute error in the probabilities or bounding the relative error in the probabilities as discussed in Al Labadi et al. (2017a). The strength of the evidence can be measured by discretizing the range of the prior distribution of d into and computing the posterior probability .

The consistency of this approach to model checking follows from results in Evans (2015). As and is true, the relative belief ratio converges to its maximum possible value (greater than 1) and the strength goes to 1 while when is false, the relative belief ratio and the strength go to 0. Note that these tests are all exact as they do not rely on asymptotics and evidence in favor of can be obtained.

Computations

For all the examples considered here it is easy to generate values of from the uniform prior on and from the associated posterior. In the case where , this leads directly to estimates of and based upon large Monte Carlo samples from these distributions. In cases where or is very small, the values of are computed for each generated value of which provides large samples from the prior and posterior distributions of this quantity. As described previously, the prior and posterior densities of are then estimated via density histograms based upon a relevant as determined by the application. This leads directly to computing as the ratio of these density estimates, as well as the straightforward computation of the posterior probability giving the strength of the evidence based upon the estimated posterior density of . The Monte Carlo sample sizes used were determined by doing the computations until the results stabilized and adding some additional samples as insurance.

Examples

Consider applications of this to the examples of Sect. 2.

Example 4

Goodness-of-fit for the trine model.

Table 1 contains data from two separate experiments discussed in Len et al. (2017) where two instances of the trine model are relevant and is the number of clicks on the i-th detector.

Table 1

Results from two experiments based on the trine model in Example 4

	n	T_1(x)	T_2(x)	T_3(x)
Symmetric	7076	3416	1912	1748
Asymmetric	6756	6192	316	248

Here with and is an angle associated with the experiment.

For the symmetric trine so and . Under the posterior distribution of is Dirichlet(3417, 1913, 1749) and sampling from this distribution shows that the entire posterior is concentrated within so . So there is evidence in favor of the symmetric trine model and this is very strong evidence since the posterior content of is effectively 1. For the asymmetric trine case so and . Under the posterior distribution of is Dirichlet(6193, 317, 249) and sampling from this distribution shows that the entire posterior is concentrated within so with strength effectively 1. So there is evidence in favor of the asymmetric trine model and again this is very strong evidence. So with both models one can feel quite confident that the true values of the probabilities lie within the respective . The evidence is definitive here because of the large amount of data.

These examples are representative of what one can expect with any of the other models listed in Example 1. Essentially these models are physically correct, at least as dictated by quantum theory, but they may not be correct based on the data collection actually carried out. In other words, if one of the physical devices was not recording clicks as it should, then the model will fail to pass the check and this is highly meaningful to the experimenter.

Example 5

Goodness-of-fit for the Hardy-Weinberg equilibrium model.

The data in Table 2 is from Laird and Lange (2011) and is concerned with CCR-5 receptor deletion which leads to protection against the HIV virus infecting certain T-cells. Let denote the dominant allele associated with deletion and denote the recessive allele. Under HWE and . The probabilities corresponding to HWE comprise a set of measure 0 so the KL distance is used to check this model via the prior and posterior distributions ofThe minimizing value is so simulation from the prior/posterior distribution of is easy.

Table 2

CCR-5 genotypes in a sample of 212 aids free men in Example 5

n	A_1A_1	A_1A_2	A_2A_2
212	4	33	175

It is necessary to specify a value for . Suppose that, if the relative errors in the probabilities are less than , then the HWE model is considered reasonable. So denoting the HWE probabilities by , we want for all (i, j) for HWE to hold. Then when the relative error is small and so it can also be expected that whenever HWE holds and not otherwise. The KL distance is nonnegative so the HWE model is checked by computing . Using , namely, a relative error of 1%, then and the strength equals 1.00, so this is virtually categorical evidence in favor. Similarly, with , then with strength 0.64 so there is some dependence of the strength of the evidence on the relative error bound.

Laird and Lange (2011) use the chi-squared goodness-of-fit test to check the HWE model which requires an asymptotic justification for the p-value and can have inaccurate type I errors as documented in Wigginton et al. (2005) who recommend an exact permutation test. Wakefield (2010) criticizes this exact test as it can be computationally intensive and the p-values need additional calibration, so a Bayes factor approach is developed. The Bayes factor is the ratio of the posterior odds in favor of to the prior odds in favor and, since has prior probability 0 for any continuous prior, it is necessary to specify a prior probability for to be true, a prior distribution concentrated on and a prior distribution on (effectively on as has 0 volume), with the overall prior a discrete mixture. While the Bayes factor can be seen as a valid characterization of evidence, in the sense that it satisfies the principle of evidence also with cut-off 1, the calibration of the Bayes factor is an unresolved issue although ad hoc scales have been proposed. By contrast the approach taken here requires the specification of a single prior, which is the uniform on , the results are exact and the relative belief ratio is calibrated via the strength. The relationship between Bayes factors and relative belief ratios is discussed in Evans (2015) where it is pointed out that the possibility of doing hypothesis assessment without the need to use a mixture prior is a distinct advantage for the relative belief approach.

Example 6

Goodness-of-fit for the multinomial with ordered probabilities.

A numerical example used in Izsák (2006), based on data concerned with fly diversity found in Papp (1992), is considered where the counts are given by so and . The question of interest is whether or not the data could have come from model (2), or even from the submodel. Note that testing the fit of a multinomial with ordered probabilities is quite different than testing the equality of binomial probabilities versus the alternative that these are ordered, as with the Cochran-Armitage test.

Consider first checking (2). For this model which is very small and so estimating with accuracy is difficult. It is to be noted, however, that given the small prior probability of this set, if any of the values generated from the posterior for some feasible sample size fall in , then this will give clear evidence in favor of the model. For example, in a sample of the posterior content was estimated as and this produces a relative belief ratio of but the strength of this evidence in favor is extremely weak. A better approach is to group the cells into subgroups such that the hypothesized monotonicity in the model is maintained. So 9 groups can be formed by combining 2 consecutive cells or 6 groups can be formed by combining 3 consecutive cells. To select which grouping to use for the test, it makes sense to start with the finest grouping such that the posterior content can be accurately estimated but coarser groupings can also be examined. Choosing 9 groups worked as the posterior content of the relevant set was estimated as 0.0396, 0.0402 and 0.0406 based on Monte Carlo samples of sizes and , respectively. The prior content of the relevant set is , the relative belief ratio is 14726 but the posterior content implies only weak evidence in favor. Using groups of size 3, however, the relative belief ratio is 285.6312 and the posterior content is 0.40. So overall it is reasonable to assume (2) holds.

Consider checking the model. The set of distributions has prior probability 0 with respect to so the KL distance measure is used. A technical difficulty is the need to compute for each generated , to obtain samples from the prior and posterior distribution of d. For this a large table of distributions was created, computed for each element and the minimum value found. There is a redundancy in the parameterization as, for example, uniformity is well-approximated by many values of . As in Example 5, a value of was selected so that if the KL distance between two distributions is less than , then this difference is irrelevant. For any , it is only necessary to consider values of such that the KL distance is greater than or equal to which places an upper bound on . Figure 1 is plot of the prior and posterior densities of d based on Monte Carlo samples of size , using and some smoothing. It is clear that the posterior has become much more concentrated about the model than the prior. Furthermore, and the strength equals 1. So there is ample evidence in favor of the model. This agrees somewhat with the finding in Izsák (2006) who conducted a goodness-of-fit test via computing a p-value based on the chi-squared statistic after grouping and found no evidence against the model. With the methodology developed here, however, there is no need to appeal to asymptotics and evidence in favor of the model has been found.

Fig. 1

Plot of prior (- - -) and posterior (—) densities of d for the model in Example 6

Checking the prior

Suppose the model has been accepted and attention now focuses on the prior . A prior is inappropriate when the prior places relatively little mass in a region of the parameter space containing the true value of the model parameter. The true value is not known but it is still possible to say something about this through the data. The approach taken in Evans and Moshonov (2006) is used here and involves comparing the observed value of a statistic with its prior predictive distribution through a certain tail probability which, when small, indicates that a prior-data conflict exists which in turn suggests that the prior is not putting its mass in the right place. For the contexts considered here, with prior having density on , the prior predictive density of the cell counts , evaluated at where and , is given by the probability function . The prior predictive probability measure associated with this probability function is then denoted by and the tail probability isThe tail probability (5) is measuring where the observed T(x) lies with respect to its prior distribution. If (5) is small, then T(x) lies in a region of low prior probability and the data and the prior are in conflict. The check on the prior should depend on the data only through T(x) because the conditional distribution of the data x given T(x), does not depend on and so implies nothing about the prior.

In Evans and Jang (2011a) it is proved that, under fairly general conditions, (5) converges to as . This supports the assertion that a small value of (5) is indicating that the prior is placing its mass in the wrong place. Note, however, that it is the discrepancy between the observed value T(x) and its prior distribution that leads to the conclusion that a prior-data conflict exists. So the goal here is to see if such a conflict exists and the fact that (5) is consistent estimator of suggests that this check is doing the right thing. In other words, our purpose is not to necessarily estimate . One of the conditions for convergence is that be continuous and here is always discrete. It was proved in Evans and Jang (2011a), however, that when a modified version of can be constructed that yields the consistency result in the binomial case and this is described in the “Appendix”. The following much more difficult result is proved in the “Appendix” for the constrained multinomial where denote the sets of continuity and discontinuity points of .

Theorem 1

Let be a prior on that satisfies the following conditions. (A1) The prior density is bounded above, that is, there exists such that for all . (A2) The prior density is continuous almost surely with respect to volume measure, that is, . (A3) The prior probability of each level set of prior density is a null set with respect to the prior, that is, for any . Then (5) converges to whenever as .

The proof also establishes a generalization.

Corollary 1

If satisfies (A1) and (A2) and is continuous at , then .

So, if is large or if is small, then should reflect this even without A3.

The value (5) is invariant under 1-1 transformations of T and under reparametrizations. So, if we reparameterize and use the transformed prior to compute , this has no effect on the limit. The theorem requires to be a continuous probability measure and cannot be constant on a subregion of having positive volume. For example, the theorem does not cover the uniform prior on but the result still holds in this case since (5) and both equal 1.

Given that the prior, like the sampling model, is a subjective choice, it is a natural scientific requirement that it should be subjected to a check against the data as all ingredients to a statistical analysis should be falsifiable and subjected to such a test. As discussed in Al Labadi and Evans (2017b), the existence of a prior-data conflict can lead to inferences that are sensitive to small changes in the prior and so a check on the prior can be viewed as a way of assessing sensitivity to the prior although such sensitivity can also arise when there is no conflict simply because of a small amount of data. When such sensitivity exists, the prior can be thought of as having too much influence on the inferences as opposed to the data. Also, the correct way to choose a prior for a problem is through elicitation where an expert translates background knowledge about the context in question into a prior probability distribution. The existence of a prior-data conflict, then indicates that there either something wrong with the elicitation process or with the knowledge the expert has expressed. Especially when a default proper prior is chosen, as is sometimes done in studies concerning sparsity, it is important that the prior be subjected to such a check to make sure the choice is not contradicted by the data. In any case something valuable has been learned after such a check has been conducted. More broadly speaking checking for prior-data conflict is about the process, as with model checking, whereby one argues that a particular statistical analysis is justified.

A natural question to consider is what is to be done when a prior is found to be in conflict. A similar question can be asked about the model when it fails its checks. For the prior, Evans and Jang (2011b) consider what it means for one prior to be weakly informative with respect to another and the definition involves specifying mathematically that the more weakly informative prior will result in fewer prior-data conflicts and this is expressed via a precise percentage. So before the data is obtained, one can specify a hierarchy of progressively more weakly informative priors starting from a base elicited prior, and if a prior-data conflict is obtained, then move up the hierarchy until this is avoided. Since the hierarchy of priors is specified before seeing the data, the final selected prior is not strongly dependent on the data. The downside of this approach is that it can be complicated to specify the hierarchy and that is certainly the case for constrained multinomials and so this problem is not pursued further here. Often, however, as with the Examples, intuition works to specify a more weakly informative prior, but this is not always the case as demonstrated in Evans and Jang (2011b).

In general the computation of (5) can be difficult but a high degree of precision is not essential for the purpose of checking the prior. For this computation an algorithm to generate is needed and this can often be accomplished by generating , then generating multinomial, which is the case for all the examples considered here. If is available in closed form (Example 8), then it is straightforward to estimate (5) via a large Monte Carlo sample from and checking the inequality in (5) for each generated . If is not available in closed form (Examples 7 and 9), then it is necessary to approximate the integral . As discussed in these examples, importance sampling can be used to estimate the integral in an inner Monte Carlo loop, for each generated as well as for T(x).

In some examples a parameterization of the Dirichlet distribution with all will be used given by the mode , where , and the concentration parameter . As , the set of all Dirichlets with this mode is indexed by . The mean and variance of the i-th coordinate equal and which converge respectively to and 0 as , so the distribution concentrates at the mode. As , the distribution converges to the uniform on .

Consider now the examples of Sect. 2.

Example 7

Checking the prior for the trine.

For a single qubit, an experimenter without any prior knowledge could assign a prior to the qubit state space that is uniform under the Hilbert-Schmidt measure. When a trine measurement is performed on the qubit, this results in the prior given by when and 0 otherwise, where c and C are as in Example 4. The transformation , where gives independent of beta(3/2, 3/2) which provides an algorithm for generating from . When n is modest, generating from and averaging the likelihood can be used to compute the values needed for (5). In this case n is large so this is too inefficient due to the concentration of the likelihood near the MLE over . The posterior under the uniform prior on also concentrates near the MLE and so importance sampling based on sampling from the Dirichlet is used to estimate for each in a sample drawn from . Sampling from is carried out by generating from and then generating multinomial. The values of are compared to to estimate (5).

This procedure was carried out for the entries in Table 1 with values of generated from and each value of estimated using a sample of from the relevant posterior based on . Prior-data conflict in this example corresponds to the true value of lying near the boundary of the respective set . For the symmetric case (5) was estimated as 0.87 and for the asymmetric case (5) was estimated as 0.15. These results were quite stable over different choices of simulation sample sizes. So in neither case is there any indication of prior-data conflict.

Example 8

Eliciting and checking the prior for the HWE model.

Consider first a Dirichlet prior on the probabilities for the alleles . There are variety of approaches to eliciting the values, see Elfadaly and Garthwaite (2013, 2017), but a method discussed in Evans, Guttman and Li (2017) is used here. This is based on prescribing a lower bound or an upper bound on each probability such that these bounds hold simultaneously with at least some prescribed prior probability . If all lower bounds are specified, then the upper bounds are immediate and of course it is required that the lower bounds satisfy which is easily satisfied by weakening some of the bounds. Similar results are available if all upper bounds or a mixture of lower and upper bounds are specified. Such bounds specify a subsimplex of the full simplex of probabilities, with all edges of the same length, and so it takes into account the dependencies among the . It is also required that the , to avoid singularities on the boundaries. So specifying the mode, for example as the centroid of the subsimplex, and then choosing so that the subsimplex contains of the prior probability, fully specifies the prior. It is natural to take large, say , to reflect the fact that one is ‘virtually certain’ that the bounds hold but still allow for the possibility that they don’t.

It also seems reasonable that one could place such bounds on the probabilities for the genotypes. This requires satisfying and . Since this constraint also implies , this is equivalent to the previous approach but with . Again it is possible to have a mixture of upper and lower bounds.

A closed form expression is available for as this is the expectation of a polynomial in with respect to a Dirichlet distribution. Then (5) is computed by generating , generating t from the relevant multinomial, computing and comparing this to , and doing this multiple times.

In the case of the CCR-5 receptor deletion data of Example 5, and so we specify and which implies . With , then so is chosen so that the beta distribution has of its probability in the interval . For illustrative purposes, suppose and so and then gives that the beta(11.35, 25.15) distribution assigns probability 0.99 to (0.1, 0.5). This leads to and (5) equal to 0.005 indicating that a prior-data conflict exists. So these bounds are inappropriate and indeed we would expect a geneticist to be able to state more realistic bounds for either the allele probabilities or the genotype probabilities for a particular population. A more weakly informative prior on will avoid prior-data conflict. For example, the bounds lead to and a beta(9.215, 45.785) prior on , which assigns 0.99 of its probability to (0.01, 0.30), and then (5) equals 0.180 and there is no prior-data conflict in this case. Of course, one could also use a “noninformative” uniform(0, 1) prior on and then (5) equals 0.241, indicating no prior-data conflict. The recommendation, however, is that an elicited informative prior always be used when this is possible and always be checked for prior-data conflict. The use of a default prior suggests that little thought relevant to the application has gone into the choice of the prior. While what seems like a natural noninformative prior exists for in this model, this is not the case in general.

Example 9

Eliciting and checking the prior for the multinomial with ordered probabilities.

It is necessary to provide an elicitation algorithm for a prior on satisfying (2). For this Lemma 1 helps considerably since (3) implies that any prior on induces a prior on satisfying (2). It is natural to choose a Dirichlet prior on but how should the be chosen? This of course depends upon what is known about the . One approach is to specify ordered probabilities and then use a prior with mode at this point. By Lemma 1 this can be accomplished by a Dirichlet prior on with mode at and the concentration can be chosen to reflect the degree of belief in this mode. This requires that the also comprise a probability distribution. The following result characterizes the choices of that are equispaced as this seems like a reasonable choice although there are many other possibilities that can be characterized similarly.

Lemma 2

The probabilities satisfying (2) are equispaced with iff and and in this case for and . Then the prior is Dirichlet.

Proof

By Lemma 1 the give probabilities satisfying (2) for some iff for and and . Since iff and this satisfies iff .

It remains to determine . For this prescribe an interval (l, u) such that it is believed all the true probabilities lie in (l, u) with virtual certainty, namely,where is a large probability like 0.99. Since , the right-hand side of (6) goes to 1 as . Therefore, satisfying (6) is easily found by simulation and the smallest such value of implies the least concentration for the prior.

For the data of Example 6 with , then , and with as bounds, the value is the estimated smallest value satisfying (6) with . When , then and so are possible and the value is the estimated smallest value satisfying (6) with . Now consider checking the prior corresponding to and . Figure 2 is a plot of density histograms for the first four probabilities based on a sample of from the full prior.

Fig. 2

The marginal priors on when and in Example 9

Our approach to computing (5) is via importance sampling. For this particular data set, , an importance sampler on with this mode and values of produces stable estimates of . A problem occurs, however, when computing other values of necessary for estimating (5). Values of are obtained by generating and multinomial. When n is not large relative to dimension, as is the case here, then typically and so the choice of importance sampler is unclear. A solution is to use an importance sampler based on Lemma 1 where the Dirichlet on has its mode at where is the convex combination of the prior mode and that just satisfies being in , the boundary of the region . This always generates points inside and should at least somewhat mimic the integrand over provided the concentration is not chosen too large or too small. Here was chosen such that both smaller and larger values led to smaller estimates of (5). When this was carried out on this example the value of (5) was estimated as 0.36 which indicates there is no prior-data conflict. This makes sense as the naive estimate is not only in but also satisfies the bounds. Suppose instead the data was observed that clearly violates the monotonicity. In this case (5) equals 0 and so prior-data conflict was detected as is correct. Also, with the original data and the prior determined by satisfying (6), then (5) equals 0 again indicating a definite conflict.

A general approach to the computational problem is to reduce dimension by grouping as in Example 6. Intuitively, as the ratio of dimension to sample size decreases, the values of the generated relative frequencies are more likely to be in as confirmed by a simulation experiment. Implementing the importance sampling for the reduced problem, however, requires the computation of the marginal prior density and this is not in closed form. Since this would be required at each generated value from the importance sampler, the computational advantage is largely negated. A modification of this is based on the elicitation algorithm and could be called marginalizing the elicitation. Consider instead being presented with the reduced problem and applying the elicitation algorithm to it. This will not result in a prior that is the marginal of the full prior but checking this prior for conflict with the data is also assessing whether or not the information being used to choose the prior is appropriate. For example, the original elicitation led to the inequality holding with virtual certainty and recall that necessarily . So, if cells are grouped in pairs to maintain the monotonicity and to make the best use of the bounds, then supposing is even, . If is odd, then the last group can consist of by itself and the lower bound doesn’t change. So for even modest k the bounds will not increase by much and clearly this idea can be extended to groups of 3, 4 etc. Lemma 1 can then be used to obtain the prior for the parameters for the grouped problem. Supposing there are m groups and denotes a value generated from the prior predictive for the reduced problem, our recommendation is that this reduction be continued until a reasonable proportion of the values lie in the reduced parameter space. When this is the case even those that lie outside should be close to the relevant set which will improve the quality of the importance sampling. Since the model has been accepted, however, it seems likely that any points outside will still be close to . For values of generated from the prior, the following values of were obtained where is the proportion that lay inside the relevant parameter space: (18, 0.00), (9, 0.04), (6, 0.25), (5, 0.44), (4, 0.61), (3, 0.78) and (2, 0.93). The values 0.41 and 0.55 were obtained for (5) when and , respectively. So one can feel fairly confident that the elicited information is not in conflict with the data.

The model was also considered but a significant problem remains unresolved. It is unclear how to elicit a prior on as the interpretation of these parameters is not obvious. The family also imposes sharper constraints on the probabilities. For example, the maximum probabilities over all for are 1.00, 0.33, 0.26, 0.25, respectively, which contrasts with 1.00, 0.50, 0.33, 0.25 for the general model for ordered probabilities. Since the general model for ordered probabilities is easier to use and interpret, it is recommended over the model.

Performance measurement

It is reasonable to ask how the model checking and checking for prior-data conflict approaches discussed here perform relative to alternatives. For checking for prior-data conflict, the only approach for which consistency results exist is, as far as we know, the one described here. While simulation results can lead to insight, nothing definitive can be concluded by such an approach to making comparisons with other methods for detecting prior-data conflict. The use of (5) has been considered in many contexts, as discussed in Evans and Moshonov (2006), with favorable results both when the prior can be considered as appropriate and otherwise. So there is no hesitancy in recommending the check for prior-data conflict discussed here.

More can said about the model checking method and this is because optimality results have been established for relative belief inferences generally. For example, consider the probabilitieswhere denotes the conditional prior predictive measure for T given that and the probability function is obtained as described previously for based on prior , but now the integration is with respect to the conditional prior . So these conditional prior predictive distributions are now obtained using the prior , appropriately conditioned, that was used for model checking. When the first probability in (7) is the prior probability of obtaining evidence in favor of the model given that the model is false at a specific distance and the second is the prior probability of obtaining evidence against the model given that the model is true. Put and so can be considered as an acceptance region for the model and can be considered as a rejection region for the model. Applying Proposition 4.7.9 of Evans (2015) implies that, among all acceptance regions A satisfying , then minimizes the prior probability of accepting the model given that it is false and, among all rejection regions satisfying , then maximizes the prior probability of rejecting the model when it is false. Actually a fair comparison demands that there be equality and and then it is seen that the relative belief approach to determining evidence is clearly optimal. These conclusions depend on the measure of distance d used, but basing this on KL divergence is a natural choice.

These results are similar to the Neyman-Pearson formulation for hypothesis testing problems. Some intuition for these results can be seen from the Savage-Dickey ratio result which implies , see Evans (2015). Also, Evans and Guo (2019) establishes close links between relative belief inferences and the ability to satisfy frequentist criteria. Of course, the general optimality results depend on the model and prior being correct but here there is no prior-data conflict from the uniform prior and some specific multinomial holds provided the sampling is i.i.d. So these results represent strong support for the model checking approach taken here. Many other good properties have been established for relative belief inferences and these are summarized in Evans (2015).

The consistency of relative belief inferences has been established in Evans (2015). So, if the model is wrong, the relative belief ratio will converge to 0 as well as the posterior probability used to measure the strength, indicating maximally strong evidence against. Similarly, when the model is correct the relative belief ratio converges to the maximal possible value and the strength converges to 1, indicating maximally strong evidence in favor. So the only real issue is whether there is enough data available in a problem to detect model failure or avoid such a finding when it is not appropriate. This in turn is dependent on the dimension of the model and, when the sample size is too low relative to the dimension, then problems can ensue, as in Example 6, and some modifications are in order. Currently there is no general result prescribing appropriate sample sizes for a constrained multinomial.

Conclusions

Constrained multinomial models arise in a number of interesting contexts and pose some unique challenges. A model checking procedure that allows for evidence in favor of as well as evidence against a constrained multinomial has been presented. A significant consistency theorem has been established for the check on the prior. As a particular application, a general model for ordered probabilities has been developed, together with an elicitation algorithm for a prior, and the results applied. Constrained multinomial models from quantum state estimation and genetics have also been discussed as applications.

There are many other models that correspond to constraints placed upon a multinomial. For example, the papers Iliopoulos et al. (2007) and Demirhan (2013) discuss Goodman’s RC model and generalizations of this as well as various constrained loglinear models. The approach discussed in this paper can be applied to check these models although in some cases there may be computational issues that need to be resolved depending on the distance measure d employed. In general, the primary issue in implementing the model check is the need to compute the distance which can be difficult when the region is complex. In such contexts it may be necessary to simplify the problem somewhat by using some lower dimensional characteristic such that the distance of , for an unrestricted value of , to the restricted set of such values can be computed. While this is not the full formal model check described in the paper, this will still serve as a check. Such issues arise with some higher dimensional quantum state models and at present there is no general solution to propose for these. It is to be noted that model checks can take a variety of forms with there being no absolutely correct approach.

The computations associated with checking the prior are, as exemplified by several of the examples, more difficult. The primary issue there is the evaluation of the prior predictive density which can be difficult to compute depending on the choice of the prior. When there is a convenient algorithm for , then implementation is often straightforward although, as exemplified by Example 9 this can be difficult if the dimension of is high relative to the amount of data. As carried out in that example, it can be necessary to make some modifications to the full formal check in such cases as in choosing a characteristic for which the induced prior is computationally tractable and then checking that prior for prior-data conflict.

The R programs associated with the examples discussed here are available at: http://utstat.utoronto.ca/mikevans/software/Constrainedmultinomial/