2
$\begingroup$

I am currently implementing double sampling. The goal is to accept or reject a lot by only sampling a few instances from it. In a two-stage testing procedure (double sampling) one first draws a sample of size $n_1$ and compares the number $K_1$ of items “of interest” (e.g., “non compliant”) with two integers: $c_1$, $c_2$. If $K_1 \le c_1$, we accept the null hypotheses ($H_0$) that the batch is of acceptable quality; and if $X > c_2$, we accept the alternative hypotheses ($H_1$) that the batch quality is unacceptable. However, if $c_1 < X \le c_2$, we draw a second sample of size $n_2$. We stop in case we find more than $c_2$. You can play with some plans in here.

double sampling flow chart

My problem is now finding the average sample size (ASN) with partial curtailment, that is we look at the first sample as a whole and then stop looking at the first sample once it is clear that we reject it.

I found this paper that computes it as

$ASN = n_1 + \sum_{i=c_1 + 1}^{c2} P(n_1, i) [ n_2 [ P_1(n_2, c_2 - i) + \frac{c_2 - i + 1}{p} P_2(n + 1, c_2 - i + 2) ]]$

where $P(n_1, i)$ is the probability of observing exactly $j$ defectives in a sample of size $n_1$, $P_1(n_2, c_2 − i)$ is the probability of observing $c_2 − i$ or fewer defectives in a sample of size $n_2$, and $P_2(n_2 + 1, c_2 − i + 2)$ is the probability of observing $c_2 − i + 2$ defectives in a sample of size $n_2 + 1$.

From when I read it, then I would think $P$, $P_1$ and $P_2$ are all PMFs, e.g. of Bernoulli $p$, but when I implement it, then I get very different curves for most of the plans compared to here.

I also had a look at

The Average Sample Number for Truncated Single and Double Attributes Acceptance Sampling Plans Author(s): C. C. Craig Source: Technometrics, Vol. 10, No. 4 (Nov., 1968), pp. 685-692

and at the end it looks like $P_1$ and $P_2$ are CDFs, but I do not manage to get the same results as their tables. Can you help me understand what I need to put in as $P$, $P_1$ and $P_2$?

$\endgroup$
2
  • $\begingroup$ The factor of $1/p$ looks fishy, because as $p$ decreases, the ASN ought to decrease, but that factor makes it explode. Your source doesn't provide a derivation--it only gives a reference but omits the bibliography. As far as I can tell, $P, P_1,$ and $P_2$ are all the same function and are explicitly referred to as Binomial probabilities, of the kinds appearing in formulas (5) through (8). $\endgroup$
    – whuber
    Commented Aug 9, 2022 at 22:05
  • $\begingroup$ @whuber Thank you for the comment! The second source has a derivation, but I find the P_2 most difficult to understand as they use a negative binomial or so to "simplify" from what I understood. Regarding $p$, I think that it does not explode because the terms before go close to zero, as having a low chance for defects would lead to having the second step not occur, leaving us with only sampling $n_1$ instances. $\endgroup$
    – jcklie
    Commented Aug 9, 2022 at 22:25

1 Answer 1

3
$\begingroup$

This comes down to an analysis of "curtailment" (which is a simplified version of sequential sampling). Let's focus on that first. Because I think your reference is incorrect, I will be careful and show the details.

Curtailment occurs during the second stage of double sampling of a process that has an (assumed) probability $q \gt 0$ at each step of yielding a "defective" result. In this stage we plan to sample for as many as $n$ steps, but we will stop sampling at any step where the cumulative number of defectives reaches a threshold $h \gt 1.$ What is the expected sample size for this stage?

It helps to be formal because we need to track two random variables:

  • Let $X(i)$ be the number of defectives observed on or before step $i.$
  • Let $N$ be the number of steps taken to first observe $h$ defectives with unlimited sampling. It is related to $X(i)$ by $N = \min\{i\mid X(i) \ge h\}.$
  • Write $N\wedge n = \min(N, n)$ for the number of samples limited by $n.$

Even if we were to terminate sampling early ($N \lt n$), it would not change the results of our calculations if we went ahead and obtained a full sample of size $n.$ Breaking up the possible outcomes according to the values of $X(n),$ we may express the expected sample size as a weighted average of the expectations for each of those values:

$$E[N\wedge n] = \sum_{k = 0}^n E[N \wedge n\mid X(n) = k] \Pr(X(n) = k).\tag{*}$$

The probabilities are Binomial: $\Pr(X(n)=k)$ is the chance of seeing exactly $k$ defectives in any sample of size $n.$ A formula is

$$\Pr(X(n) = k) = \binom{n}{k} q^k(1-q)^{n-k} = \frac{n!}{k!(n-k)!} q^k(1-q)^{n-k}.$$

We must figure out the conditional expectations in $(*),$ which for brevity I will write as

$$g(n,h,k) = E[N \wedge n\mid X(n) = k].$$

That is, knowing there are $k$ defectives among the next $n$ observations, when will we stop sampling? Clearly, when $k \lt h$ we will obtain all $n$ samples. Otherwise, when $k\ge h,$ we will stop sampling by step $n$ (and likely sooner than that).

Consider the first observation among those $n.$ There are two possibilities: it is defective or not. The key idea is that because there are $k$ defectives among the next $n$ observations, the chance that the next observation is defective is $k/n.$ (The information given by the numbers $(n,k)$ does not determine where in the sequence $1,2,\ldots, n$ the defectives might occur: all subsets of $k$ of those positions are equally likely.) Consequently

$$g(n,h,k) = 1 + \frac{i}{n} g(n-1,h-1,k-1) + \frac{n-i}{n} g(n-1, h, k).$$

The initial $1$ counts the first observation and the next two terms use the law of conditional expectations to find the expected number of additional observations: $g(n-1,h-1,k-1)$ is the expected sample size with $n-1$ observations left to go; the threshold to terminate sampling has been reduced to $h-1$ because we just saw a defective; and there remain $k-1$ defectives among these observations. $g(n-1, h, k)$ is the expected sample size with $n-1$ observations left to go and we haven't yet seen a defective.

The (unique) solution to this recursion (consistent with obvious starting conditions) is

$$g(n, h, k) = \frac{h(n+1)}{k+1}.\tag{**}$$

To see why, note that the formula is correct in all edge cases; e.g., when $k=n,$ every one of the next $n$ observations is a defective and so we stop after seeing the first $h = h(n+1)/(n+1)$ of them. Then all we need to check is the recursion, finding

$$\begin{aligned} g(n, h, k) &= \frac{h(n+1)}{k+1}\\ &= 1 + \frac{k}{n}\frac{(h-1)n}{k} + \frac{n-k}{n} \frac{hn}{k+1}\\ &= 1 + \frac{k}{n} g(n-1,h-1,k-1) + \frac{n-k}{n} g(n-1, h, k). \end{aligned}$$

with straightforward algebra. Writing $F(\ ;n,q)$ for the Binomial$(n,q)$ distribution function we may now express $(*)$ in the closed form

$$\begin{aligned} E[N\wedge n] &= n\Pr(X(n)\lt h) + \sum_{k=h}^n \frac{h(n+1)}{k+1} \binom{n}{k}q^k(1-q)^{n-k}\\ &= nF(h-1;n,q) + h\sum_{k=h}^n \binom{n+1}{k+1}q^k(1-q)^{n-k}\\ &= nF(h-1;n,q) + h\sum_{k=h+1}^{n+1} \binom{n+1}{k}q^{k-1}(1-q)^{n+1-k}\\ &= nF(h-1;n,q) + \frac{h}{q}\sum_{k=h+1}^{n+1} \binom{n+1}{k}q^{k}(1-q)^{n+1-k}\\ &= nF(h-1;n,q) + \frac{h}{q}\left(1 - F(h;n+1,q)\right). \end{aligned}.$$

The first step breaks the expectation into the cases $X(n)\lt h,$ where the full sample will be obtained, and the remaining cases $X(n)\ge h,$ where early termination is possible. A little algebra and a change in the index from $k$ to $k+1$ yields the result.


We are ready to compute the expected sample size for double sampling with curtailment. The initial sample of $n_1$ observations will be completed, contributing $n_1$ to the expectation. If the number of defectives in that sample is $c_1$ or less or if it exceeds $c_2,$ sampling stops with a size of $n_1.$ The lower (red) and upper (purple) paths in the figure illustrate these scenarios.

Figure showing four sampling paths

Otherwise, additional sampling is needed. Sampling terminates either when the cumulative number of defectives exceeds $c_2$ (as in the teal path) or the total sample size reaches $n_1 + n_2$ (green path). For paths like the teal, $N$ (marked with a solid vertical line) can have values between $n_1$ and $n_1 + n_2.$

The expected extra effort is (as before) computed by partitioning the possibilities according to the value of $X(n_1),$ which ranges from $c_1+1$ through $c_2.$ Let this value be $j.$ The threshold to terminate this second stage of sampling with a maximum sample size of $n=n_2$ is $h = c_2+1 - j$ because $j$ defectives have already been observed. Applying the formula $(**)$ gives

$$\begin{aligned} &E[N] \\ &= n_1 + \sum_{j = c_1 + 1}^{c_2} E[N-n_1 \mid X(n_1) = j]\Pr(X(n_1) = j)\\ &= n_1 + \sum_{j = c_1 + 1}^{c_2} \left[nF(c_2-j;n_1,q) + \frac{h}{q}\left(1 - F(c_2+1 - j;n+1,q)\right)\right]\binom{n_1}{j}q^j(1-q)^{n_1-j}. \end{aligned}\tag{***}$$

Although this looks a lot like formula $(13)$ in your reference, it is not the same, because the $F(\cdots)$ values are not the Binomial probabilities described there: they are cumulative probabilities and cannot be equated with or converted into individual probabilities. (This makes me suspect the Technometrics article, which I haven't seen, is likely correct.)

Here is a plot of $E[N]$ as a function of $q$ for the sampling plan in the previous illustration ($c_1=4,$ $c_2 =14,$ $n_1=30,$ and $n_2=40$).

Figure 2

This is a general pattern: the sampling effort will be greatest when the rate of defectives is in a "sweet spot" enabling the sample path to stay out of the gray areas. Otherwise it's highly likely sampling will end after the first round of size $n_1.$ Even so, the expectation is substantially less than the maximum planned sample size of $n_1+n_2.$


To settle the question of whose formula (if either) is correct, we need a third approach. Let's simulate some samples in a way that does not depend on much if any analysis. Small errors will be most apparent in small problems, so here is an example with $n_1=10,$ $n_2=10,$ $c_1=3,$ $c_2=7,$ and $q=0.5.$

#
# The formula.
#
f <- function(n, h, q) n * pbinom(h-1, n, q) + h * (1 - pbinom(h, n+1, q)) / q
a <- Vectorize(function(c, n, q) {
  j <- seq(c[1]+1, c[2])
  n[1] + sum(dbinom(j, n[1], q) * f(n[2], c[2]+1-j, q))
}, "q")
#
# Simulation.
#
cc <- c(3, 7)
nn <- c(10,10)
q <- 0.5
set.seed(17)
N <- replicate(5e3, {
  k1 <- rbinom(1, nn[1], q)         # Take the first sample
  if(k1 <= cc[1] || k1 > cc[2]) nn[1] else {
    pop <- rbinom(nn[2], 1, q)      # Take the second sample, one at a time
    x <- cumsum(pop) + k1           # Compute the sample path
    x[nn[2]] <- cc[2] + 1           # Force termination at the end
    nn[1] + match(TRUE, x > cc[2])  # Find where sampling can first terminate
  }
})
signif(c(Simulation = mean(N), Formula = a(cc, nn, q)), 4)

The output prints the average sample length in 5,000 simulated samples along with what this formula says:

Simulation    Formula 
     14.00      14.04 

They are statistically indistinguishable. Repeated experiments of this nature bear out the correctness of formula $(***).$

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.