How to tell the probability of failure if there were no failures?

Question

I was wondering if there is a way to tell the probability of something failing (a product) if we have 100,000 products in the field for 1 year and with no failures? What is the probability that one of the next 10,000 products sold fail?

Answer 1

The probability that a product will fail is surely a function of time and use. We don't have any data on use, and with only one year there are no failures (congratulations!). Thus, this aspect (called the survival function), cannot be estimated from your data.

You can think of failures within one year as draws from a binomial distribution, however. You still have no failures, but this is now a common problem. A simple solution is to use the rule of 3, which is accurate with large $N$ (which you certainly have). Specifically, you can get the upper bound of a one-sided 95% confidence interval (i.e., the lower bound is $0$) on the true probability of failure within one year as $3/N$. In your case, you are 95% confident that the rate is less than $0.00003$.

You also asked how to compute the probability that one or more of the next 10k fails. A quick and simple (albeit extreme) way to extend the above analysis is to just use the upper bound as the underlying probability and use the corresponding binomial CDF to get the probability that there won't be $0$ failures. Using R code, we could do: 1-pbinom(0, size=10000, prob=0.00003), which yields a 0.2591851 chance of seeing one or more failures in the next 10k products. By having used the upper bound, this is not the optimal point estimate of the probability of having at least one failure, rather you can say it is very unlikely that the probability of $\ge 1$ failure is more than $\approx 26\%$ (recognizing that this is a somewhat 'hand-wavy' framing). Another possibility is to use @amoeba's suggestion of the estimate from Laplace's rule of succession. The rule of succession states that the estimated probability of failure is $(F+1)/(N+2)$, where $F$ is the number of failures. In that case, $\hat p = 9.9998\times 10^{-06}$, and the calculation for the predicted probability of $1^+$ failures in the next 10,000 is 1-pbinom(0, size=10000, prob=9.9998e-06), yielding 0.09516122, or $\approx 10\%$.

Answer 2

You can take a bayesian approach. denote the probability of failure by $\Theta$ and think of it as a random variable. A priori, before you see the results of the experiments, you might believe that $\Theta \sim U(0,1)$. If you trust the engineers to make this product reliable, maybe you can take $\Theta \sim U(0,0.1)$ or so. This is up to you. Then, you can use Bayes' theorem to calculate the posterior distribution of $\theta$. Denote $A$ the event that you've observed ($n$ experiments with zero failures).

$$ p(\Theta = \theta | A) = \frac{p (A | \Theta = \theta) p(\Theta = \theta )}{p(A)} = \frac{p (A |\theta) p(\theta )}{\int p (A |\theta) p(\theta )d\theta}. $$ Everything is simple: $\Theta$ is uniform, so $p(\theta)$ is some constant. Since you run $n$ experiments, $p(A | \theta)$ is just the probability of no failures in $n$ bernouli trials with probability of failure $\theta$.

Once you have $p(\theta | A)$ you're gold: you can calculate the probability of any event $B$ by integrateion: $\mathbb{P}(B) = \int p(B |\theta) p(\theta |A) d\theta$

Below, I work through a detailed solution, following the above approach. I'll take a few standard shortcuts.

Let the prior be $U(0,1)$. Then: $$ p(\theta |A)\propto p(A|\theta) \cdot 1 = (1-\theta)^n. $$ The normalization constant $p(A) = \int p(A|\theta)p(\theta) d\theta$ is found to be $B(1,n+1)$ - see wikipedia pages beta function and beta distribution. So, $p(\theta |A) = \frac{(1-\theta)^n}{B(1,n+1)}$, which is a beta distribution with parameters $1, n+1$.

Denote the probability of no failures in $m$ products in the next year by $B$. The probability of at least one failure is $1 -\mathbb{P}( B )$. Then $$ 1- \mathbb{P}(B) =1 - \int (1-\theta)^m\frac{(1-\theta)^n}{B(1,n+1)}d\theta = \frac{B(1,n+m+1)}{B(1,n+1)} $$

which is roughly $0.1$, using $n= 100,000, m = 10,000$. Not very impressive? I took a uniform distribution on the probability of failure. Perhaps you have better prior faith in your engineers.

Answer 3

Rather than computing a probability, why not predict how many products might fail?

Modeling the Observations

There are $n=100000$ products in the field and another $m=10000$ under consideration. Assume their failures are all independent and constant with probability $p$.

We may model this situation by means of a Binomial experiment: out of a box of tickets with an unknown proportion $p$ of "failure" tickets and $1-p$ "success" tickets, draw $m+n=110000$ tickets (with replacement, so that the chance of failure stays the same). Count the failures among the first $n$ tickets--let that be $X$--and count the failures among the remaining $m$ tickets, calling that $Y$.

Framing the Question

In principle, $0\le X \le n$ and $0 \le Y\le m$ could be anything. What we are interested in is the chance that $Y = u$ given that $X+Y=u$ (with $u$ any number in $\{0,1,\ldots, m\}$). Since the failures could occur anywhere among all $n+m$ tickets, with every possible configuration having the same chance, it is found by dividing the number of $u$-subsets of $m$ things by the number of $u$-subsets of all $n+m$ things:

$$p(u;n,m) = \Pr(Y = u\,|\, X+Y=u) = \frac{\binom{m}{u}}{\binom{n+m}{u}} \\= \frac{m(m-1)\cdots(m-u+1)}{(n+m)(n+m-1)\cdots(n+m-u+1)}.$$

Comparable formulas can be used for the calculation when $X=1, 2, \ldots.$

An upper $1-\alpha$ prediction limit (UPL) for the number of failures in those last $m$ tickets, $t_\alpha(X;n,m)$, is given by the smallest $u$ (depending on $X$) for which $p(u;n,m) \le \alpha$.

Interpretation

The UPL should be interpreted in terms of the risk of using $t_\alpha$, as evaluated before either $X$ or $Y$ is observed. In other words, suppose it is one year ago and you are being asked to recommend a procedure to predict the number of failures in the next $m$ products once the first $n$ have been observed. Your client asks

What is the chance that your procedure will underpredict $Y$? I don't mean in the future after you have more data; I mean right now, because I have to make decisions right now and the only chances I will have available to me are the ones that can be computed at this moment."

Your response can be,

Right now the chance is no greater than $\alpha$, but if you plan to use a smaller prediction, the chance will exceed $\alpha$.

Results

For $n=10^5$, $m=10^4$, and $X=0$ we may compute that

$$p(0,n,m)=1;\ p(1,n,m)=\frac{1}{11}\approx 0.091;\ p(2,n,m)=\frac{909}{109999}\approx 0.0083; \ldots$$

Thus, upon having observed $X=0$,

• For up to $1-\alpha=90.9\%$ confidence (that is, when $9.1\%\le \alpha$), predict there is at most $t_\alpha(0;n,m)=1$ failure in the next $10,000$ products.

• For up to $99.2\%$ confidence (that is, when $0.8\%\le \alpha \lt 9.1\%$), predict there are at most $t_\alpha(0;n,m)=2$ failures in the next $10,000$ products.

• Etc.

---

Comments

When and why would this approach apply? Suppose your company makes lots of different products. After observing the performance of $n$ of each one in the field, it likes to produce guarantees, such as "complete no-cost replacement of any failure within one year." By having prediction limits for the number of failures you can control the total costs of having to back those guarantees. Because you make many products, and expect failures to be due to random circumstances beyond your control, the experience of each product will be independent. It makes sense to control your risk in the long run. Every once in a while you might have to pay more claims than expected, but most of the time you will pay fewer. If paying more than announced could be ruinous, you will set $\alpha$ to be extremely small (and you likely would use a more sophisticated failure model, too!). Otherwise, if the costs are minor, then you can live with low confidence (high $\alpha$). These calculations show how to balance confidence and risks.

Note that we don't have to compute the full procedure $t$. We wait until $X$ is observed and then just carry out the calculations for that particular $X$ (here, $X=0$), as shown above. In principle, though, we could have carried out the calculations for all possible values of $X$ at the outset.

A Bayesian approach (described in other answers) is attractive and will work well provided the results do not depend heavily on the prior. Unfortunately, when the failure rate is so low that very few (or no failures) are observed, the results are sensitive to the choice of prior.

Answer 4

The following is a Bayesian answer to "Out of 10,000 new products, how many are expected to fail if all the former 100,000 produced didn't fail?", but you should consider the sensitivity to different priors.

Suppose that $X_1,\dots,X_n$ are conditionally independent and identically distributed, given $\Theta=\theta$, such that $X_1\mid\Theta=\theta\sim\mathrm{Bernoulli}(\theta)$, and use the conjugate prior $\Theta\sim\mathrm{Beta}(a,b)$, with $a,b>0$.

For $m<n$, we have $$ \mathrm{E}\left[\sum_{i=m+1}^n X_i\;\Bigg\vert\; X_1=0,\dots X_m=0 \right] = \sum_{i=m+1}^n \mathrm{E}\left[ X_i\mid X_1=0,\dots X_m=0 \right] \, . $$

For $m+1\leq i\leq n$, we have $$ \begin{align} \mathrm{E}\left[X_i\mid X_1=0,\dots X_m=0\right] &= \Pr(X_i=1\mid X_1=0,\dots X_m=0) \\ &= \int_0^1 \Pr(X_i=1\mid \Theta=\theta) \,f_{\Theta\mid X_1,\dots,X_m}(\theta\mid 0,\dots,0) \,d\theta \\ &= \frac{\Gamma(m+a+b)}{\Gamma(m+a+b+1)} \frac{\Gamma(a+1)}{\Gamma(a)} = \frac{a}{m+a+b}\, , \end{align} $$ in which we used $\Theta\mid X_1=0,\dots,X_m=0\sim \mathrm{Beta}(a,m+b)$.

Plugging in your numbers, with an uniform prior ($a=1,b=1$) you expect a failure rate around $10\%$, while a Jeffreys-like prior ($a=1/2,b=1/2$) gives you a failure rate close to $5\%$.

This predictive expectation doesn't look like a good summary, because the predictive distribution is highly skewed. We can go further and compute the predictive distribution. Since $$ \sum_{i=m+1}^n X_i \;\Bigg\vert\; \Theta=\theta \sim \mathrm{Bin}(n-m+2,\theta) \, , $$ conditioning as we did before we have $$ \begin{align} \Pr&\left(\sum_{i=m+1}^n X_i=t \;\Bigg\vert\; X_1=0,\dots X_m=0\right) = \\ &\qquad\qquad\qquad\qquad\binom{n-m+2}{t} \frac{\Gamma(m+a+b)}{\Gamma(a)\Gamma(m+b)} \frac{\Gamma(t+a)\Gamma(n-t+2)}{\Gamma(n+a+2)} \, , \end{align} $$ for $t=0,1,\dots,n-m+2$.

I'll finish it later computing a $95\%$ predictive interval.

Answer 5

Several good answers were provided for this question, but recently I had a chance to review few resources on this topic and so I decided to share the results.

There are multiple possible estimators for zero-failures data. Let's denote $k=0$ as number of failures and $n$ as sample size. Maximum likelihood estimator for probability of failure given this data is

$$ P(K = k) = \frac{k}{n} = 0 \tag{1} $$

Such estimate is rather unsatisfactory since the fact that we observed no failures in our sample hardly proves that they are impossible in general. Out out-of-data knowledge suggests that there is some probability of failure even if non were observed (yet). Having a priori knowledge leads us to using Bayesian methods reviewed by Bailey (1997), Razzaghi (2002), Basu et al (1996), and Ludbrook and Lew (2009).

Among simple estimators "upper bound" estimator that assumes (Bailey, 1997)

that it would not be logical for an estimator for P in the zero-failure case to yield a probability in excess of that predicted by the maximum likelihood estimator in the one-failure case, a reasonable upper bound

defined as

$$ \frac{1}{n} \tag{2} $$

can be mentioned. As reviewed by Ludbrook and Lew (2009), other possibilities are "rule of threes" (cf. here, Wikipedia, or Eypasch et al, 1995)

$$ \frac{3}{n} \tag{3} $$

or other variations:

$$ \frac{3}{n+1} \tag{4} $$

"rule of 3.7" by Newcombe and Altman (or by 3.6):

$$ \frac{3.7}{n} \tag{5} $$

"new rule of four":

$$ \frac{4}{n+4} \tag{6} $$

but as concluded by Ludbrook and Lew (2009) "rule of threes" is "next to useless" and "rule of 3.6" (and 3.7) "have serious limitations – they are grossly inaccurate if the initial sample size is less than 50" and they do not recommend methods (3)-(6), suggesting rather to use proper Bayesian estimators (see below).

Among Bayesian estimators several different can be mentioned. First such estimator suggested by Bailey (1997) is

$$ 1 - 0.5^\frac{1}{n} \tag{7} $$

for estimating median under uniform prior

$$ 1 - 0.5^\frac{1}{n+1} \tag{8} $$

or for estimating mean under such prior

$$ \frac{1}{n+2} \tag{9} $$

yet another approach assuming exponential failure pattern with constant failure rate (Poisson distributions) yields

$$ \frac{1/3}{n} \tag{10} $$

if we use beta prior with parameters $a$ and $b$ we can use formula (see Razzaghi, 2002):

$$ \frac{a}{a+b+n} \tag{11} $$

that under $a = b = 1$ leads to uniform prior (9). Assuming Jeffreys prior with $a = b = 0.5$ it leads to

$$ \frac{1}{2(n+1)} \tag{12} $$

Generally, Bayesian formulas (7)-(12) are recommended. Basu et al (1996) recommends (11) with informative prior, when some a priori knowledge is available. Since no single best method exists I would suggest reviewing the literature prior to your analysis, especially when $n$ is small.

---

Bailey, R.T. (1997). Estimation from zero-failure data. Risk Analysis, 17, 375-380.

Razzaghi, M. (2002). On the estimation of binomial success probability with zero occurrence in sample. Journal of Modern Applied Statistical Methods, 1(2), 41.

Ludbrook, J., & Lew, M. J. (2009). Estimating the risk of rare complications: is the ‘rule of three’good enough?. ANZ journal of surgery, 79(7‐8), 565-570.

Eypasch, E., Lefering, R., Kum, C.K., and Troidl, H. (1995). Probability of adverse events that have not yet occurred: A statistical reminder. BMJ 311(7005): 619–620.

Basu, A.P., Gaylor, D.W., & Chen, J.J. (1996). Estimating the probability of occurrence of tumor for a rare cancer with zero occurrence in a sample. Regulatory Toxicology and Pharmacology, 23(2), 139-144.

Answer 6

Using Laplace's sunrise problem approach, we get the probability that a product would fail within a year $$p=\frac{1}{100000+1}$$. Next, the probability that of $n$ new products none fails within a year is $$(1-p)^n$$ Hence, the probability that at least one product of $n$ will fail in next year is $$1-\left(1-\frac{1}{100001}\right)^{n}$$ For $n=10000$ the value is $P_{10000}\approx 0.095$. In whuber's case $P_{200000}\approx 0.87$, quite high, in fact.

Of course, you should keep updating your data while more products are sold, eventually one will fail.

Answer 7

You really need to go back to the designers of your products. It is a fundamental engineering problem not an observational statistical one. They will have an idea of the failure probability of each component and from that the net failure probabilty of the total assembled product. They can give you the expected number of failures over the whole design life of the product.

A civil engineer designs a bridge to have a design life of 120 years. Each component of the bridge has a slight chance of failure. Each loading has a slight chance of being exceeded. To make the bridge economic to build, total collapse would only occur once in 2400 years which is far longer than the bridge will be maintained for. It is not surprising that the bridge does not fail in year 1, nor year 2 to year 120. That is has not collapsed tells you very little. Its various chances of failure with time can only be estimated by the original designers.

Answer 8

This is similar to a problem I faced when we introduced a new manufacturing process to eliminate a failure in production.

The new system produced no failures so people were asking the same question: how do we predict the failure rate? In your case, because you have stipulated a period over which the failure can occur with no concern for when the failure occurs within that period, the temporal effects have been removed. And it is simply a case of whether something failed or not. With that stipulated - on with my answer.

Intuitively, it seems we need at least one failure to be able to calculate the failure rate. However, this assumption has an implicit mistake within it. We will never calculate the failure rate. That is because we are dealing with a sample. Thus we can only estimate a range of probable failure rates. The way to do this is to find a distribution for the failure rate. The distribution that does the job in this instance is a Beta distribution where the parameters are: α = n + 1 and β = N - n + 1

Note: N is the sample size and n is the number of failures (in your case 0)

For your scenario, the distribution of the of the failure rate is shown below. .

You would then feed that distribution into the respective binomial probability formula to get a distribution for the probability of one unit failing (could be done analytically or using Monte Carlo). I suspect that numbers will be very low.

Note that this process is applicable no matter the number of failures in your fist set.

Answer 9

Median unbiased estimates can be used to estimate sample proportions and (non-singular) 95% CIs in Bernoulli samples with no variability. In a sample with no positive cases, you can estimate the upper bound of a 95% confidence interval with the following formula:

$$ p_{1-\alpha/2} : P(Y=0)/2 + P(Y>y) > 0.975$$

that is we seek a value $p_{1-\alpha/2}$ as the upper bound of the CI so that the Bernoulli process with probability $p=p_{1-\alpha/2}$ gives the above probability inequality. In R this is solved with a NR-like uniroot application.

set.seed(12345)
y <- rbinom(100, 1, 0.01) ## all 0
cil <- 0
mupfun <- function(p) {
  0.5*dbinom(0, 100, p) + 
    pbinom(1, 100, p, lower.tail = F) - 
    0.975
} ## for y=0 successes out of n=100 trials
ciu <- uniroot(mupfun, c(0, 1))$root
c(cil, ciu)

[1] 0.00000000 0.05357998 ## includes the 0.01 actual probability