# Combining probabilities of nuclear accidents

The recent events in Japan have made me think about the following.

Nuclear Plants are usually designed to limit risk of serious accidents to a 'design basis probability' for example, say, 10E-6/year. This is the criteria for a single plant. However, when there is a population of hundreds of reactors, how do we combine the individual probabilities of a serious accident ? I know I could probably research this myself but having found this site I a sure there is someone that will be able to answer this question quite easily. Thanks

• The nuclear situation in Japan is a Black Swan event. According to N.N. Taleb, Black Swan are events of very low probability but have a very high impact. His assertion is that such probabilities are uncomputable, and any computed probabilities have very little bearing on real life. – Gilead Mar 19 '11 at 15:23
• en.wikipedia.org/wiki/Black_swan_theory – Gilead Mar 19 '11 at 15:30
• Taleb, [ cringe ]. – cardinal Mar 19 '11 at 15:30
• @cardinal, one often wishes the conveyor of such ideas wasn't a guy like Taleb (whose personality can be overbearing). But I wouldn't dismiss the ideas because of the man. – Gilead Mar 19 '11 at 15:37
• I've read each of his books. Though interesting, I'd say, few, if any, ideas are his. He's been quite successful at popularizing them, though. I've also read a bit of the literature that he cites. Some of it I feel he misrepresents for his own purposes. That perturbs me. – cardinal Mar 19 '11 at 15:47

To answer the pure probabilistic question that J Presley presented, using bayer's notation (p=probability of an item failing), the the probability of at least one element failing is 1-P(none fail)= 1-(1-p)^n. This type of calculation is common in system reliability where a bunch of components are linked in parallel, so that the system continues to function if at least one component is functioning.

You can still use this formula even if each plant item has a different failure probability (p_i). The formula would then be 1- (1-p_1)(1-p_2)...(1-p_n).

• Thank you Gail.....that is exactly the solution I wanted. By the way ...is there a general series (power,taylor or other) expansion for (1-p)^n that you know of ? – user3777 Mar 21 '11 at 2:37
• I have officially accepted Galit's (sorry I called you Gail) answer although Bayer's answer came close to the answer I was expecting to my original question. – user3777 Mar 24 '11 at 22:41

Before you set up your analysis, keep in mind the reality of what the current situation involves.

This meltdown was not directly caused by the earthquake or the tsunami. It was because of a lack of back-up power. If they had enough back-up power, regardless of the earthquake/tsunami, they could have kept the cooling water running, and none of the meltdowns would have happened. The plant would probably be back up and running by now.

Japan, for whatever reason, has two electrical frequencies (50 Hz and 60 Hz). And, you can't run a 50 Hz motor at 60 Hz or vice versa. So, whatever frequency the plant was using/providing is the frequency they need to power up. "U.S. type" equipment runs at 60 Hz and "European type" equipment runs at 50 Hz, so in providing an alternative power source, keep that in mind.

Next, that plant is in a fairly remote mountainous area. To supply external power requires a LONG power line from another area (requiring days/weeks to build) or large gasoline/diesel driven generators. Those generators are heavy enough that flying them in with a helicopter is not an option. Trucking them in may also be a problem due to the roads being blocked from the earthquake/tsunami. Bringing them in by ship is an option, but it also takes days/weeks.

The bottom line is, the risk analysis for this plant comes down to a lack of SEVERAL (not just one or two) layers of back-ups. And, because this reactor is an "active design", which means it requires power to stay safe, those layers are not a luxury, they're required.

This is an old plant. A new plant would not be designed this way.

Edit (03/19/2011) ==============================================

J Presley: To answer your question requires a short explanation of terms.

As I said in my comment, to me, this is a matter of "when", not "if", and as a crude model, I suggested the Poisson Distribution/Process. The Poisson Process is a series of events that happen at an average rate over time (or space, or some other measure). These events are independent of each other and random (no patterns). The events happen one at a time (2 or more events don't happen at the exact same time). It is basically a binomial situation ("event" or "no event") where the probability that the event will happen is relatively small. Here are some links:

http://en.wikipedia.org/wiki/Poisson_process

http://en.wikipedia.org/wiki/Poisson_distribution

Next, the data. Here's a list of nuclear accidents since 1952 with the INES Level:

I count 19 accidents, 9 state an INES Level. For those without an INES level, all I can do is assume the level is below Level 1, so I'll assign them Level 0.

So, one way to quantify this is 19 accidents in 59 years (59 = 2011 -1952). That's 19/59 = 0.322 acc/yr. In terms of a century, that's 32.2 accidents per 100 years. Assuming a Poisson Process gives the following graphs.

Originally, I suggested a Lognormal, Gamma, or Exponential Distribution for the severity of the accidents. However, since the INES Levels are given as discrete values, the distribution would need to be discrete. I would suggest either the Geometric or Negative Binomial Distribution. Here are their descriptions:

http://en.wikipedia.org/wiki/Negative_binomial_distribution

http://en.wikipedia.org/wiki/Geometric_distribution

They both fit the data about the same, which is not very well (lots of Level 0's, one Level 1, zero Level 2's, etc).

 Fit for Negative Binomial Distribution

Fitting of the distribution ' nbinom ' by maximum likelihood
Parameters :
estimate Std. Error
size 0.460949  0.2583457
mu   1.894553  0.7137625
Loglikelihood:  -34.57827   AIC:  73.15655   BIC:  75.04543
Correlation matrix:
size           mu
size 1.0000000000 0.0001159958
mu   0.0001159958 1.0000000000

#====================
Fit for Geometric Distribution

Fitting of the distribution ' geom ' by maximum likelihood
Parameters :
estimate Std. Error
prob 0.3454545  0.0641182
Loglikelihood:  -35.4523   AIC:  72.9046   BIC:  73.84904


The Geometric Distribution is a simple one parameter function while the Negative Binomial Distribution is a more flexible two parameter function. I would go for the flexibility, plus the underlying assumptions of how the Negative Binomial Distribution was derived. Below is a graph of the fitted Negative Binomial Distribution.

Below is the code for all this stuff. If anyone finds a problem with my assumptions or coding, don't be afraid to point it out. I checked through the results, but I didn't have enough time to really chew on this.

 library(fitdistrplus)

#Generate the data for the Poisson plots
x <- dpois(0:60, 32.2)
y <- ppois(0:60, 32.2, lower.tail = FALSE)

#Cram the Poisson Graphs into one plot
par(pty="m", plt=c(0.1, 1, 0, 1), omd=c(0.1,0.9,0.1,0.9))
par(mfrow = c(2, 1))

#Plot the Probability Graph
plot(x, type="n", main="", xlab="", ylab="", xaxt="n", yaxt="n")
mtext(side=3, line=1, "Poisson Distribution Averaging 32.2 Nuclear Accidents Per Century", cex=1.1, font=2)
xaxisdat <- seq(0, 60, 10)
pardat <- par()
yaxisdat <- seq(pardat$yaxp[1], pardat$yaxp[2], (pardat$yaxp[2]-pardat$yaxp[1])/pardat$yaxp[3]) axis(2, at=yaxisdat, labels=paste(100*yaxisdat, "%", sep=""), las=2, padj=0.5, cex.axis=0.7, hadj=0.5, tcl=-0.3) mtext("Probability", 2, line=2.3) abline(h=yaxisdat, col="lightgray") abline(v=xaxisdat, col="lightgray") lines(x, type="h", lwd=3, col="blue") #Plot the Cumulative Probability Graph plot(y, type="n", main="", xlab="", ylab="", xaxt="n", yaxt="n") pardat <- par() yaxisdat <- seq(pardat$yaxp[1], pardat$yaxp[2], (pardat$yaxp[2]-pardat$yaxp[1])/pardat$yaxp[3])
mtext("Cumulative Probability", 2, line=2.3)
abline(h=yaxisdat, col="lightgray")
abline(v=xaxisdat, col="lightgray")
lines(y, type="h", lwd=3, col="blue")

mtext("Number of Nuclear Accidents Per Century", 1, line=1)
legend("topright", legend=c("99% Probability - 20 Accidents or More", " 1% Probability - 46 Accidents or More"), bg="white", cex=0.8)

#Calculate the 1% and 99% values
qpois(0.01, 32.2, lower.tail = FALSE)
qpois(0.99, 32.2, lower.tail = FALSE)

#Fit the Severity Data
z <- c(rep(0,10), 1, rep(3,2), rep(4,3), rep(5,2), 7)
zdis <- fitdist(z, "nbinom")
plot(zdis, lwd=3, col="blue")
summary(zdis)


Edit (03/20/2011) ======================================================

J Presley: I'm sorry I couldn't finish this up yesterday. You know how it is on weekends, lots of duties.

The last step in this process is to assemble a simulation using the Poisson Distribution to determine when an event happens, and then the Negative Binomial Distribution to determine the severity of the event. You might run 1000 sets of "century chunks" to generate the 8 probability distributions for Level 0 through Level 7 events. If I get the time, I might run the simulation, but for now, the description will have to do. Maybe someone reading this stuff will run it. After that is done, you'll have a "base case" where all of the events are assumed to be INDEPENDENT.

Obviously, the next step is to relax one or more of the above assumptions. An easy place to start is with the Poisson Distribution. It assumes that all events are 100% independent. You can change that in all sorts of ways. Here are some links to Non-homogeneous Poisson Distributions:

http://www.math.wm.edu/~leemis/icrsa03.pdf

http://filebox.vt.edu/users/pasupath/papers/nonhompoisson_streams.pdf

The same idea goes for the Negative Binomial Distribution. This combination will lead you down all sorts of paths. Here are some examples:

http://surveillance.r-forge.r-project.org/

http://www.m-hikari.com/ijcms-2010/45-48-2010/buligaIJCMS45-48-2010.pdf

http://www.michaeltanphd.com/evtrm.pdf

The bottom line is, you asked a question where the answer depends on how far you want to take it. My guess is, someone, somewhere will be commissioned to generate "an answer" and will be surprised at how long it takes to do the work.

Edit (03/21/2011) ====================================================

I had a chance to slap together the above mentioned simulation. The results are shown below. From the original Poisson Distribution, the simulation provides eight Poisson Distributions, one for each INES Level. As the severity level rises (INES Level Number rises), the number of expected events per century drops. This may be a crude model, but it's a reasonable place to start.

• How big are those generators? I would have guessed either a Skycrane or Mi-26 could haul them in, at the very least in pieces. – cardinal Mar 18 '11 at 19:44
• There are (at least) two reasons for inadequate backup power in my understanding...1. the tidal wave took out the fuel tanks for the standby generators (inadequate tidal wave protection) 2. Inadequate batteries to keep essential equipment running until replacement power was available (likely impractical). Both of these situations are part of a large and complex probabilistic safety analysis of multiple scenarios. However, the bottom line is...the lower your probability criteria is....the more stringent your design will be (cont'd) – user3777 Mar 18 '11 at 23:33
• As a former nuclear reactor designer I am unaware of anyone who ever considered the 'total reactor population of the world' when estimating risk. The last few days have made me wonder whether this should not be the case in future. This is what prompted my question. – user3777 Mar 18 '11 at 23:34
• Why would one use such a strict criteria ? Because the consequences of such (potentially) low frequency events are so large we need to try to eliminate them completely. Again, economics will limit just how much we can do in this regard. – user3777 Mar 18 '11 at 23:39
• @JPresley: If I had to do such a calculation, I think it's more of a "when" than an "if" situation. As a simple model, I would use a Poisson Distribution for the "when", and maybe a Lognormal Distribution (Gamma??, Exponential??) for the magnitude of the problem. That's why several layers of backups/contingencies are necessary. – bill_080 Mar 19 '11 at 1:05

The underlying difficulty behind the question is that situations that have been anticipated, have generally been planned for, with mitigation measures in place. Which means that the situation should not even turn into a serious accident.

The serious accidents stem from unanticipated situations. Which means that you cannot assess probabilities for them - they are your Rumsfeldian unknown unknowns.

The assumption of independence is clearly invalid - Fukushima Daiichi shows that. Nuclear plants can have common-mode failures. (i.e. more than one reactor becoming unavailable at once, due to a common cause).

Although probabilities cannot be quantitatively calculated, we can make some qualitative assertions about common-mode failures.

For example: if the plants are all built to the same design, then they are more likely to have common-mode failures (for example the known problem with pressurizer cracks in EPRs / PWRs)

If the plant sites share geographic commonalities, they are more likely to have common-mode failures: for example, if they all lie on the same earthquake fault line; or if they all rely on similar rivers within a single climatic zone for cooling (when a very dry summer can cause all such plants to be taken offline).

• Agreed -- it is folly to assign a probability to so-called Fourth Quadrant events such as this, or to even think we can predict them. All we can do is to make the system robust to their negative effects through redundancies etc. – Gilead Mar 19 '11 at 15:32
• I don't fully agree. The Tsunami was not unanticipated, the 'level' of the Tsunami was unanticipated. The plant was 'apparently' designed for a 7 meter Tsunami based on historical probabilty information. This was considered acceptable based on some probabilistic arguments by someone. If the criteria was more stringent than a 'less probable' wave hejght would have been required in the design...this is my point... – user3777 Mar 20 '11 at 0:06

As commentators pointed out, this has the very strong independency assumption.

Let the probability that a plant blows up be $p$. Then the probability that a plant does not blow up is $1-p$. Then the probability that $n$ plants do not blow up is $(1-p)^n$. The expected number of plants blown up per year is $np$.

In case you're interested: binomial distribution.

• @bayer, I am not voting this down (though I'm a bit tempted), but the independence assumption strikes me as wholely inappropriate in this circumstance and would lead to absurd inferences! – cardinal Mar 18 '11 at 15:12
• I'm with @cardinal; independent failures is a ridiculous assumption. What if, say, the plants are near each other and in an area of high tectonic activity... – JMS Mar 18 '11 at 15:30
• @cardinal is absolutely right: this is the crux of the matter. Engineers have used these kinds of independence assumptions without considering the possibility that everything could go wrong at once due to a common cause (such as an earthquake). That's (apparently) why multiple backups have failed in some Japanese installations. – whuber Mar 18 '11 at 15:34
• I think in this particular case, assuming independence could mean miscalculating the relevant probabilities by several orders of magnitude. I am not an expert in either nuclear power plant design nor aircraft design and logistics. However, I would hazard a guess that the dependence structure of the risk to Japanese power plants is substantially stronger than for calculating probabilities of aircraft crashes. The nuclear power plants are simultaneously subject to the same risk factors: Off the top of my head: (a) earthquake, (b) tsunami, (c) electrical power grid, (d) common manufacturer,... – cardinal Mar 19 '11 at 14:18
• ...(e) other geographical/geophysical design risks. Just the fact that they are all in close proximity of one another makes them have a common risk factor for simultaneous attack, let's say. The point is that you have to consider all these possible risk factors and for a lot of these, conditional on one of them occuring, the probabilities of multiple plant failures simultaneously jumps to close to one. Plane crashes seem a bit different as they are largely automonous and operate more widely geographically. Now, if air-traffic control over New York (TRACON) went down completely, all of a... – cardinal Mar 19 '11 at 14:24