Compute the value of \(K\) (the multiplier of estimated standard deviation) used
to construct a simultaneous prediction interval based on data from a
normal distribution.
The function
predIntNormSimultaneousK is called by predIntNormSimultaneous.
predIntNormSimultaneousK(n, df = n - 1, n.mean = 1, k = 1, m = 2, r = 1,
rule = "k.of.m", delta.over.sigma = 0, pi.type = "upper", conf.level = 0.95,
K.tol = .Machine$double.eps^0.5, integrate.args.list = NULL)a positive integer greater than 2 indicating the sample size upon which the prediction interval is based.
the degrees of freedom associated with the prediction interval. The default is
df=n-1.
positive integer specifying the sample size associated with the future averages.
The default value is n.mean=1 (i.e., individual observations). Note that all
future averages must be based on the same sample size.
for the \(k\)-of-\(m\) rule (rule="k.of.m"), a positive integer
specifying the minimum number of observations (or averages) out of \(m\)
observations (or averages) (all obtained on one future sampling “occassion”)
the prediction interval should contain with confidence level conf.level.
The default value is k=1. This argument is ignored when the argument
rule is not equal to "k.of.m".
positive integer specifying the maximum number of future observations (or
averages) on one future sampling “occasion”.
The default value is m=2, except when rule="Modified.CA", in which
case this argument is ignored and m is automatically set equal to 4.
positive integer specifying the number of future sampling “occasions”.
The default value is r=1.
character string specifying which rule to use. The possible values are
"k.of.m" (\(k\)-of-\(m\) rule; the default), "CA" (California rule),
and "Modified.CA" (modified California rule).
See the DETAILS section below for more information.
numeric scalar indicating the ratio \(\Delta/\sigma\). The quantity
\(\Delta\) (delta) denotes the difference between the mean of the population
that was sampled to construct the prediction interval, and the mean of the
population that will be sampled to produce the future observations. The quantity
\(\sigma\) (sigma) denotes the population standard deviation for both populations.
See the DETAILS section below for more information. The default value is
delta.over.sigma=0.
character string indicating what kind of prediction interval to compute.
The possible values are pi.type="upper" (the default), and
pi.type="lower".
a scalar between 0 and 1 indicating the confidence level of the prediction interval.
The default value is conf.level=0.95.
numeric scalar indicating the tolerance to use in the nonlinear search algorithm to
compute \(K\). The default value is K.tol=.Machine$double.eps^(1/2).
For many applications, the value of \(K\) needs to be known only to the second
decimal place, in which case setting K.tol=1e-4 will speed up computation a
bit.
A numeric scalar equal to \(K\), the multiplier of estimated standard deviation that is used to construct the simultaneous prediction interval.
What is a Simultaneous Prediction Interval?
A prediction interval for some population is an interval on the real line constructed
so that it will contain \(k\) future observations from that population
with some specified probability \((1-\alpha)100\%\), where
\(0 < \alpha < 1\) and \(k\) is some pre-specified positive integer.
The quantity \((1-\alpha)100\%\) is called
the confidence coefficient or confidence level associated with the prediction
interval. The function predIntNorm computes a standard prediction
interval based on a sample from a normal distribution.
The function predIntNormSimultaneous computes a simultaneous prediction
interval that will contain a certain number of future observations with probability
\((1-\alpha)100\%\) for each of \(r\) future sampling “occasions”,
where \(r\) is some pre-specified positive integer. The quantity \(r\) may
refer to \(r\) distinct future sampling occasions in time, or it may for example
refer to sampling at \(r\) distinct locations on one future sampling occasion,
assuming that the population standard deviation is the same at all of the \(r\)
distinct locations.
The function predIntNormSimultaneous computes a simultaneous prediction
interval based on one of three possible rules:
For the \(k\)-of-\(m\) rule (rule="k.of.m"), at least \(k\) of
the next \(m\) future observations will fall in the prediction
interval with probability \((1-\alpha)100\%\) on each of the \(r\) future
sampling occasions. If obserations are being taken sequentially, for a particular
sampling occasion, up to \(m\) observations may be taken, but once
\(k\) of the observations fall within the prediction interval, sampling can stop.
Note: When \(k=m\) and \(r=1\), the results of predIntNormSimultaneous
are equivalent to the results of predIntNorm.
For the California rule (rule="CA"), with probability
\((1-\alpha)100\%\), for each of the \(r\) future sampling occasions, either
the first observation will fall in the prediction interval, or else all of the next
\(m-1\) observations will fall in the prediction interval. That is, if the first
observation falls in the prediction interval then sampling can stop. Otherwise,
\(m-1\) more observations must be taken.
For the Modified California rule (rule="Modified.CA"), with probability
\((1-\alpha)100\%\), for each of the \(r\) future sampling occasions, either the
first observation will fall in the prediction interval, or else at least 2 out of
the next 3 observations will fall in the prediction interval. That is, if the first
observation falls in the prediction interval then sampling can stop. Otherwise, up
to 3 more observations must be taken.
Simultaneous prediction intervals can be extended to using averages (means) in place
of single observations (USEPA, 2009, Chapter 19). That is, you can create a
simultaneous prediction interval
that will contain a specified number of averages (based on which rule you choose) on
each of \(r\) future sampling occassions, where each each average is based on
\(w\) individual observations. For the functions
predIntNormSimultaneous and predIntNormSimultaneousK,
the argument n.mean corresponds to \(w\).
The Form of a Prediction Interval
Let \(\underline{x} = x_1, x_2, \ldots, x_n\) denote a vector of \(n\)
observations from a normal distribution with parameters
mean=\(\mu\) and sd=\(\sigma\). Also, let \(w\) denote the
sample size associated with the future averages (i.e., n.mean=\(w\)).
When \(w=1\), each average is really just a single observation, so in the rest of
this help file the term “averages” will sometimes replace the phrase
“observations or averages”.
For a normal distribution, the form of a two-sided \((1-\alpha)100\%\)
simultaneous prediction interval is:
$$[\bar{x} - Ks, \bar{x} + Ks] \;\;\;\;\;\; (1)$$
where \(\bar{x}\) denotes the sample mean:
$$\bar{x} = \frac{1}{n} \sum_{i=1}^n x_i \;\;\;\;\;\; (2)$$
\(s\) denotes the sample standard deviation:
$$s^2 = \frac{1}{n-1} \sum_{i=1}^n (x_i - \bar{x})^2 \;\;\;\;\;\; (3)$$
and \(K\) denotes a constant that depends on the sample size \(n\), the
confidence level, the number of future sampling occassions \(r\), and the
sample size associated with the future averages, \(w\). Do not confuse the
constant \(K\) (uppercase K) with the number of future averages \(k\)
(lowercase k) in the \(k\)-of-\(m\) rule. The symbol \(K\) is used here
to be consistent with the notation used for tolerance intervals
(see tolIntNorm).
Similarly, the form of a one-sided lower prediction interval is: $$[\bar{x} - Ks, \infty] \;\;\;\;\;\; (4)$$ and the form of a one-sided upper prediction interval is: $$[-\infty, \bar{x} + Ks] \;\;\;\;\;\; (5)$$
Note: For simultaneous prediction intervals, only lower
(pi.type="lower") and upper
(pi.type="upper") prediction intervals are available.
The derivation of the constant \(K\) is explained below.
The Derivation of K for Future Observations
First we will show the derivation based on future observations (i.e.,
\(w=1\), n.mean=1), and then extend the formulas to future averages.
The Derivation of K for the k-of-m Rule (rule="k.of.m")
For the \(k\)-of-\(m\) rule (rule="k.of.m") with \(w=1\)
(i.e., n.mean=1), at least \(k\) of the next \(m\) future
observations will fall in the prediction interval
with probability \((1-\alpha)100\%\) on each of the \(r\) future sampling
occasions. If observations are being taken sequentially, for a particular
sampling occasion, up to \(m\) observations may be taken, but once \(k\)
of the observations fall within the prediction interval, sampling can stop.
Note: When \(k=m\) and \(r=1\), this kind of simultaneous prediction
interval becomes the same as a standard prediction interval for the next
\(k\) observations (see predIntNorm).
For the case when \(r=1\) future sampling occasion, both Hall and Prairie (1973)
and Fertig and Mann (1977) discuss the derivation of \(K\). Davis and McNichols
(1987) extend the derivation to the case where \(r\) is a positive integer. They
show that for a one-sided upper prediction interval (pi.type="upper"), the
probability \(p\) that at least \(k\) of the next \(m\) future observations
will be contained in the interval given in Equation (5) above, for each of \(r\)
future sampling occasions, is given by:
$$p = \int_0^1 T(\sqrt{n}K; n-1, \sqrt{n}[\Phi^{-1}(v) + \frac{\Delta}{\sigma}]) r[I(v; k, m+1-k)]^{r-1} [\frac{v^{k-1}(1-v)^{m-k}}{B(k, m+1-k)}] dv \;\;\;\;\;\; (6)$$
where \(T(x; \nu, \delta)\) denotes the cdf of the
non-central Student's t-distribution with parameters
df=\(\nu\) and ncp=\(\delta\) evaluated at \(x\);
\(\Phi(x)\) denotes the cdf of the standard normal distribution
evaluated at \(x\); \(I(x; \nu, \omega)\) denotes the cdf of the
beta distribution with parameters shape1=\(\nu\) and
shape2=\(\omega\); and \(B(\nu, \omega)\) denotes the value of the
beta function with parameters a=\(\nu\) and
b=\(\omega\).
The quantity \(\Delta\) (upper case delta) denotes the difference between the
mean of the population that was sampled to construct the prediction interval, and
the mean of the population that will be sampled to produce the future observations.
The quantity \(\sigma\) (sigma) denotes the population standard deviation of both
of these populations. Usually you assume \(\Delta=0\) unless you are interested
in computing the power of the rule to detect a change in means between the
populations (see predIntNormSimultaneousTestPower).
For given values of the confidence level (\(p\)), sample size (\(n\)),
minimum number of future observations to be contained in the interval per
sampling occasion (\(k\)), number of future observations per sampling occasion
(\(m\)), and number of future sampling occasions (\(r\)), Equation (6) can
be solved for \(K\). The function predIntNormSimultaneousK uses the
R function nlminb to solve Equation (6) for \(K\).
When pi.type="lower", the same value of \(K\) is used as when
pi.type="upper", but Equation (4) is used to construct the prediction
interval.
The Derivation of K for the California Rule (rule="CA")
For the California rule (rule="CA"), with probability \((1-\alpha)100\%\),
for each of the \(r\) future sampling occasions, either the first observation will
fall in the prediction interval, or else all of the next \(m-1\) observations will
fall in the prediction interval. That is, if the first observation falls in the
prediction interval then sampling can stop. Otherwise, \(m-1\) more observations
must be taken.
The formula for \(K\) is the same as for the \(k\)-of-\(m\) rule, except that Equation (6) becomes the following (Davis, 1998b): $$p = \int_0^1 T(\sqrt{n}K; n-1, \sqrt{n}[\Phi^{-1}(v) + \frac{\Delta}{\sigma}]) r\{v[1+v^{m-2}(1-v)]\}^{r-1} [1+v^{m-2}(m-1-mv)] dv \;\;\;\;\;\; (7)$$
The Derivation of K for the Modified California Rule (rule="Modified.CA")
For the Modified California rule (rule="Modified.CA"), with probability
\((1-\alpha)100\%\), for each of the \(r\) future sampling occasions, either the
first observation will fall in the prediction interval, or else at least 2 out of
the next 3 observations will fall in the prediction interval. That is, if the first
observation falls in the prediction interval then sampling can stop. Otherwise, up
to 3 more observations must be taken.
The formula for \(K\) is the same as for the \(k\)-of-\(m\) rule, except that Equation (6) becomes the following (Davis, 1998b): $$p = \int_0^1 T(\sqrt{n}K; n-1, \sqrt{n}[\Phi^{-1}(v) + \frac{\Delta}{\sigma}]) r\{v[1+v(3-v[5-2v])]\}^{r-1} \{1+v[6-v(15-8v)]\} dv \;\;\;\;\;\; (8)$$
The Derivation of K for Future Means
For each of the above rules, if we are interested in using averages instead of
single observations, with \(w \ge 1\) (i.e., n.mean\(\ge 1\)), the first
term in the integral in Equations (6)-(8) that involves the cdf of the
non-central Student's t-distribution becomes:
$$T(\sqrt{n}K; n-1, \frac{\sqrt{n}}{\sqrt{w}}[\Phi^{-1}(v) + \frac{\sqrt{w}\Delta}{\sigma}]) \;\;\;\;\;\; (9)$$
Barclay's California Code of Regulations. (1991). Title 22, Section 66264.97 [concerning hazardous waste facilities] and Title 23, Section 2550.7(e)(8) [concerning solid waste facilities]. Barclay's Law Publishers, San Francisco, CA.
Davis, C.B. (1998a). Ground-Water Statistics \& Regulations: Principles, Progress and Problems. Second Edition. Environmetrics \& Statistics Limited, Henderson, NV.
Davis, C.B. (1998b). Personal Communication, September 3, 1998.
Davis, C.B., and R.J. McNichols. (1987). One-sided Intervals for at Least \(p\) of \(m\) Observations from a Normal Population on Each of \(r\) Future Occasions. Technometrics 29, 359--370.
Fertig, K.W., and N.R. Mann. (1977). One-Sided Prediction Intervals for at Least \(p\) Out of \(m\) Future Observations From a Normal Population. Technometrics 19, 167--177.
Gibbons, R.D., D.K. Bhaumik, and S. Aryal. (2009). Statistical Methods for Groundwater Monitoring, Second Edition. John Wiley & Sons, Hoboken.
Hahn, G.J. (1969). Factors for Calculating Two-Sided Prediction Intervals for Samples from a Normal Distribution. Journal of the American Statistical Association 64(327), 878-898.
Hahn, G.J. (1970a). Additional Factors for Calculating Prediction Intervals for Samples from a Normal Distribution. Journal of the American Statistical Association 65(332), 1668-1676.
Hahn, G.J. (1970b). Statistical Intervals for a Normal Population, Part I: Tables, Examples and Applications. Journal of Quality Technology 2(3), 115-125.
Hahn, G.J. (1970c). Statistical Intervals for a Normal Population, Part II: Formulas, Assumptions, Some Derivations. Journal of Quality Technology 2(4), 195-206.
Hahn, G.J., and W.Q. Meeker. (1991). Statistical Intervals: A Guide for Practitioners. John Wiley and Sons, New York.
Hahn, G., and W. Nelson. (1973). A Survey of Prediction Intervals and Their Applications. Journal of Quality Technology 5, 178-188.
Hall, I.J., and R.R. Prairie. (1973). One-Sided Prediction Intervals to Contain at Least \(m\) Out of \(k\) Future Observations. Technometrics 15, 897--914.
Millard, S.P. (1987). Environmental Monitoring, Statistics, and the Law: Room for Improvement (with Comment). The American Statistician 41(4), 249--259.
Millard, S.P., and Neerchal, N.K. (2001). Environmental Statistics with S-PLUS. CRC Press, Boca Raton, Florida.
USEPA. (2009). Statistical Analysis of Groundwater Monitoring Data at RCRA Facilities, Unified Guidance. EPA 530/R-09-007, March 2009. Office of Resource Conservation and Recovery Program Implementation and Information Division. U.S. Environmental Protection Agency, Washington, D.C.
USEPA. (2010). Errata Sheet - March 2009 Unified Guidance. EPA 530/R-09-007a, August 9, 2010. Office of Resource Conservation and Recovery, Program Information and Implementation Division. U.S. Environmental Protection Agency, Washington, D.C.
predIntNormSimultaneous,
predIntNormSimultaneousTestPower,
predIntNorm, predIntNormK,
predIntLnormSimultaneous, tolIntNorm,
Normal, estimate.object, enorm
# NOT RUN {
# Compute the value of K for an upper 95% simultaneous prediction
# interval to contain at least 1 out of the next 3 observations
# given a background sample size of n=8.
predIntNormSimultaneousK(n = 8, k = 1, m = 3)
#[1] 0.5123091
#----------
# Compare the value of K for a 95% 1-of-3 upper prediction interval to
# the value for the California and Modified California rules.
# Note that the value of K for the Modified California rule is between
# the value of K for the 1-of-3 rule and the California rule.
predIntNormSimultaneousK(n = 8, k = 1, m = 3)
#[1] 0.5123091
predIntNormSimultaneousK(n = 8, m = 3, rule = "CA")
#[1] 1.252077
predIntNormSimultaneousK(n = 8, rule = "Modified.CA")
#[1] 0.8380233
#----------
# Show how the value of K for an upper 95% simultaneous prediction
# limit increases as the number of future sampling occasions r increases.
# Here, we'll use the 1-of-3 rule.
predIntNormSimultaneousK(n = 8, k = 1, m = 3)
#[1] 0.5123091
predIntNormSimultaneousK(n = 8, k = 1, m = 3, r = 10)
#[1] 1.363002
#==========
# Example 19-1 of USEPA (2009, p. 19-17) shows how to compute an
# upper simultaneous prediction limit for the 1-of-3 rule for
# r = 2 future sampling occasions. The data for this example are
# stored in EPA.09.Ex.19.1.sulfate.df.
# We will pool data from 4 background wells that were sampled on
# a number of different occasions, giving us a sample size of
# n = 25 to use to construct the prediction limit.
# There are 50 compliance wells and we will monitor 10 different
# constituents at each well at each of the r=2 future sampling
# occasions. To determine the confidence level we require for
# the simultaneous prediction interval, USEPA (2009) recommends
# setting the individual Type I Error level at each well to
# 1 - (1 - SWFPR)^(1 / (Number of Constituents * Number of Wells))
# which translates to setting the confidence limit to
# (1 - SWFPR)^(1 / (Number of Constituents * Number of Wells))
# where SWFPR = site-wide false positive rate. For this example, we
# will set SWFPR = 0.1. Thus, the confidence level is given by:
nc <- 10
nw <- 50
SWFPR <- 0.1
conf.level <- (1 - SWFPR)^(1 / (nc * nw))
conf.level
#[1] 0.9997893
#----------
# Compute the value of K for the upper simultaneous prediction
# limit for the 1-of-3 plan.
predIntNormSimultaneousK(n = 25, k = 1, m = 3, r = 2,
rule = "k.of.m", pi.type = "upper", conf.level = conf.level)
#[1] 2.014365
#==========
# Cleanup
#--------
rm(nc, nw, SWFPR, conf.level)
# }
Run the code above in your browser using DataLab