| Title: | Bayes Screening and Model Discrimination |
|---|---|
| Description: | Bayes screening and model discrimination follow-up designs. |
| Authors: | Ernesto Barrios based on Daniel Meyer's code. |
| Maintainer: | Ernesto Barrios <[email protected]> |
| License: | GPL (>= 3) |
| Version: | 2023.920 |
| Built: | 2024-09-13 08:04:17 UTC |
| Source: | CRAN |
Bayes screening and model discrimination follow-up designs
| Package: | BsMD |
| Type: | Package |
| Version: | 2023.920 |
| Date: | 2023-09-14 |
| License: | GPL version 2 or later |
The packages allows you to perform the calculations and analyses described in Mayer, Stainberg and Box paper in Technometrics, 1996.
Author: Ernesto Barrios based on Daniel Meyer's code. Maintainer: Ernesto Barrios <[email protected]>
Box and Mayer, 1986; Box and Mayer, 1993; Mayer, Steinberg and Box, 1996.
data(BM86.data)data(BM86.data)
Design factors and responses used in the examples of Box and Meyer (1986)
data(BM86.data)data(BM86.data)
A data frame with 16 observations on the following 19 variables.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Contrast factor.
numeric vector. Log drill advance response.
numeric vector. Tensile strength response.
numeric vector. Shrinkage response.
numeric vector. Yield of isatin response.
Box, G. E. P and R. D. Meyer (1986). "An Analysis of Unreplicated Fractional Factorials". Technometrics. Vol. 28. No. 1. pp. 11–18.
library(BsMD) data(BM86.data,package="BsMD") print(BM86.data)library(BsMD) data(BM86.data,package="BsMD") print(BM86.data)
12-run Plackett-Burman design from the $2^5$ reactor example from Box, Hunter and Hunter (1977).
data(BM93.e1.data)data(BM93.e1.data)
A data frame with 12 observations on the following 7 variables.
a numeric vector. Run number from a $2^5$ factorial design in standard order.
a numeric vector. Feed rate factor.
a numeric vector. Catalyst factor.
a numeric vector. Agitation factor.
a numeric vector. Temperature factor.
a numeric vector. Concentration factor.
a numeric vector. Percent reacted response.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (1978). Statistics for Experimenters. Wiley.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
library(BsMD) data(BM93.e1.data,package="BsMD") print(BM93.e1.data)library(BsMD) data(BM93.e1.data,package="BsMD") print(BM93.e1.data)
12-run Plackett-Burman design for the study of fatigue life of weld repaired castings.
data(BM93.e2.data)data(BM93.e2.data)
A data frame with 12 observations on the following 8 variables.
a numeric vector. Initial structure factor.
a numeric vector. Bead size factor.
a numeric vector. Pressure treat factor.
a numeric vector. Heat treat factor.
a numeric vector. Cooling rate factor.
a numeric vector. Polish factor.
a numeric vector. Final treat factor.
a numeric vector. Natural log of fatigue life response.
Hunter, G. B., Hodi, F. S., and Eager, T. W. (1982). "High-Cycle Fatigue of Weld Repaired Cast Ti-6A1-4V". Metallurgical Transactions 13A, pp. 1589–1594.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
library(BsMD) data(BM93.e2.data,package="BsMD") print(BM93.e2.data)library(BsMD) data(BM93.e2.data,package="BsMD") print(BM93.e2.data)
Fractional factorial design in the injection molding example from
Box, Hunter and Hunter (1978).
data(BM93.e3.data)data(BM93.e3.data)
A data frame with 20 observations on the following 10 variables.
a numeric vector
a numeric vector. Mold temperature factor.
a numeric vector. Moisture content factor.
a numeric vector. Holding Pressure factor.
a numeric vector. Cavity thickness factor.
a numeric vector. Booster pressure factor.
a numeric vector. Cycle time factor.
a numeric vector. Gate size factor.
a numeric vector. Screw speed factor.
a numeric vector. Shrinkage response.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (1978). Statistics for Experimenters. Wiley.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (2004). Statistics for Experimenters II. Wiley.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
library(BsMD) data(BM93.e3.data,package="BsMD") print(BM93.e3.data)library(BsMD) data(BM93.e3.data,package="BsMD") print(BM93.e3.data)
Marginal factor posterior probabilities and model posterior probabilities from designed screening experiments are calculated according to Box and Meyer's Bayesian procedure.
BsProb(X, y, blk, mFac, mInt = 2, p = 0.25, g = 2, ng = 1, nMod = 10)BsProb(X, y, blk, mFac, mInt = 2, p = 0.25, g = 2, ng = 1, nMod = 10)
X |
Matrix. The design matrix. |
y |
vector. The response vector. |
blk |
integer. Number of blocking factors (>=0). These factors are
accommodated in the first columns of matrix |
mFac |
integer. Maximum number of factors included in the models. |
mInt |
integer <= 3. Maximum order of interactions considered in the models. |
p |
numeric. Prior probability assigned to active factors. |
g |
vector. Variance inflation factor(s) |
ng |
integer <=20. Number of different variance inflation factors ( |
nMod |
integer <=100. Number of models to keep with the highest posterior probability. |
Factor and model posterior probabilities are computed by Box and Meyer's Bayesian
procedure. The design factors are accommodated in the matrix X after
blk columns of the blocking factors. So, ncol(X)-blk design factors
are considered. If g, the variance inflation factor (VIF) ,
is a vector of length 1, the same VIF is used for factor main effects and interactions.
If the length of g is 2 and ng is 1, g[1] is used for factor
main effects and g[2] for the interaction effects. If ng greater than 1,
then ng values of VIFs between g[1] and g[2] are used for
calculations with the same value for main effects and interactions.
The function calls the FORTRAN subroutine ‘bm’ and captures summary results.
The complete output of the FORTRAN code is save in the ‘BsPrint.out’
file in the working directory. The output is a list of class BsProb for which
print, plot and summary methods are available.
A list with all output parameters of the FORTRAN subroutine ‘bm’. The names of the list components are such that they match the original FORTRAN code. Small letters used for capturing program's output.
X |
matrix. The design matrix. |
Y |
vector. The response vector. |
N |
integer. The number of runs. |
COLS |
integer. The number of design factors. |
BLKS |
integer. The number of blocking factors accommodated in the first
columns of matrix |
MXFAC |
integer. Maximum number of factors considered in the models. |
MXINT |
integer. Maximum interaction order considered in the models. |
PI |
numeric. Prior probability assigned to the active factors. |
INDGAM |
integer. If 0, the same variance inflation factor ( |
INDG2 |
integer. If 1, the variance inflation factor |
NGAM |
integer. Number of different VIFs used for computations. |
GAMMA |
vector. Vector of variance inflation factors of length 1 or 2. |
NTOP |
integer. Number of models with the highest posterior probability |
.
mdcnt |
integer. Total number of models evaluated. |
ptop |
vector. Vector of probabilities of the top |
sigtop |
vector. Vector of sigma-squared of the top |
nftop |
integer. Number of factors in each of the |
jtop |
matrix. Matrix of the number of factors and their labels
of the top |
del |
numeric. Interval width of the |
sprob |
vector. Vector of posterior probabilities. If |
pgam |
vector. Vector of values of the unscaled posterior density of |
prob |
matrix. Matrix of marginal factor posterior probabilities for each of
the different values of |
ind |
integer. Indicator variable. |
The function is a wrapper to call the FORTRAN subroutine ‘bm’, modification of Daniel Meyer's original program, ‘mbcqp5.f’, for the application of Bayesian design and analysis of fractional factorial experiments, part of the mdopt bundle, available at StatLib.
R. Daniel Meyer. Adapted for R by Ernesto Barrios.
Box, G. E. P and R. D. Meyer (1986). "An Analysis for Unreplicated Fractional Factorials". Technometrics. Vol. 28. No. 1. pp. 11–18.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
print.BsProb, print.BsProb, summary.BsProb.
library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) plot(drillAdvance.BsProb) summary(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) plot(drillAdvance.BsProbG, code = FALSE, prt = TRUE)library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) plot(drillAdvance.BsProb) summary(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) plot(drillAdvance.BsProbG, code = FALSE, prt = TRUE)
Normal plot of effects from a two level factorial experiment.
DanielPlot(fit, code = FALSE, faclab = NULL, block = FALSE, datax = TRUE, half = FALSE, pch = "*", cex.fac = par("cex.lab"), cex.lab = par("cex.lab"), cex.pch = par("cex.axis"), ...)DanielPlot(fit, code = FALSE, faclab = NULL, block = FALSE, datax = TRUE, half = FALSE, pch = "*", cex.fac = par("cex.lab"), cex.lab = par("cex.lab"), cex.pch = par("cex.axis"), ...)
fit |
object of class |
code |
logical. If |
faclab |
list. If |
block |
logical. If |
datax |
logical. If |
half |
logical. If |
pch |
numeric or character. Points character. |
cex.fac |
numeric. Factors' labels character size. |
cex.lab |
numeric. Labels character size. |
cex.pch |
numeric. Points character size. |
... |
graphical parameters passed to |
The two levels design are assumed -1 and 1. Factor effects assumed 2*coef(obj)
((Intercept) removed) are displayed in a qqnorm plot with the effects in
the x-axis by default. If half=TRUE the half-normal plots of effects is
plotted as the normal quantiles of 0.5*(rank(abs(effects))-0.5)/length(effects)+1
versus abs(effects).
The function returns invisible data frame with columns: x, y
and no, for the coordinates and the enumeration of plotted points.
Names of the factor effects (coefficients) are the row names of the data frame.
Ernesto Barrios.
C. Daniel (1976). Application of Statistics to Industrial Experimentation. Wiley.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (1978). Statistics for Experimenters. Wiley.
### Injection Molding Experiment. Box et al. 1978. library(BsMD) # Data data(BM86.data,package="BsMD") # Design matrix and response print(BM86.data) # from Box and Meyer (1986) # Model Fitting. Box and Meyer (1986) example 3. injectionMolding.lm <- lm(y3 ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 + X14 + X15, data = BM86.data) print(coef(injectionMolding.lm)) # Model coefficients # Daniel Plots par(mfrow=c(1,3),oma=c(0,0,1,0),pty="s") DanielPlot(injectionMolding.lm, half = TRUE, main = "Half-Normal Plot") DanielPlot(injectionMolding.lm, main = "Normal Plot of Effects") DanielPlot(injectionMolding.lm, faclab = list(idx = c(12,4,13), lab = c(" -H"," VG"," -B")), main = "Active Contrasts")### Injection Molding Experiment. Box et al. 1978. library(BsMD) # Data data(BM86.data,package="BsMD") # Design matrix and response print(BM86.data) # from Box and Meyer (1986) # Model Fitting. Box and Meyer (1986) example 3. injectionMolding.lm <- lm(y3 ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 + X14 + X15, data = BM86.data) print(coef(injectionMolding.lm)) # Model coefficients # Daniel Plots par(mfrow=c(1,3),oma=c(0,0,1,0),pty="s") DanielPlot(injectionMolding.lm, half = TRUE, main = "Half-Normal Plot") DanielPlot(injectionMolding.lm, main = "Normal Plot of Effects") DanielPlot(injectionMolding.lm, faclab = list(idx = c(12,4,13), lab = c(" -H"," VG"," -B")), main = "Active Contrasts")
Plot of the factor effects with significance levels based on robust estimation of contrast standard errors.
LenthPlot(obj, alpha = 0.05, plt = TRUE, limits = TRUE, xlab = "factors", ylab = "effects", faclab = NULL, cex.fac = par("cex.lab"), cex.axis=par("cex.axis"), adj = 1, ...)LenthPlot(obj, alpha = 0.05, plt = TRUE, limits = TRUE, xlab = "factors", ylab = "effects", faclab = NULL, cex.fac = par("cex.lab"), cex.axis=par("cex.axis"), adj = 1, ...)
obj |
object of class |
alpha |
numeric. Significance level used for the margin of error (ME) and simultaneous margin of error (SME). See Lenth(1989). |
plt |
logical. If |
limits |
logical. If |
xlab |
character string. Used to label the x-axis. "factors" as default. |
ylab |
character string. Used to label the y-axis. "effects" as default. |
faclab |
list with components |
cex.fac |
numeric. Character size used for the factor labels. |
cex.axis |
numeric. Character size used for the axis. |
adj |
numeric between 0 and 1. Determines where to place the
"ME" (margin of error) and the "SME" (simultaneous margin of error) labels
(character size of 0.9* |
... |
extra parameters passed to |
If obj is of class lm, 2*coef(obj) is used as factor
effect with the intercept term removed. Otherwise, obj should be a
vector with the factor effects. Robust estimate of the contrasts standard
error is used to calculate marginal (ME) and simultaneous margin
of error (SME) for the provided significance (1 - alpha) level.
See Lenth(1989). Spikes are used to display the factor effects.
If faclab is NULL, factors are labelled with the effects or
coefficient names. Otherwise, those faclab\$idx factors are labelled
as faclab\$lab. The rest of the factors are blanked.
The function is called mainly for its side effect. It returns a vector with the value of alpha used, the estimated PSE, ME and SME.
Ernesto Barrios. Extension provided by Kjetil Kjernsmo (2013).
Lenth, R. V. (1989). "Quick and Easy Analysis of Unreplicated Factorials". Technometrics Vol. 31, No. 4. pp. 469–473.
DanielPlot, BsProb and plot.BsProb
### Tensile Strength Experiment. Taguchi and Wu. 1980 library(BsMD) # Data data(BM86.data,package="BsMD") # Design matrix and responses print(BM86.data) # from Box and Meyer (1986) # Model Fitting. Box and Meyer (1986) example 2. tensileStrength.lm <- lm(y2 ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 + X14 + X15, data = BM86.data) print(coef(tensileStrength.lm)) # Model coefficients par(mfrow=c(1,2),pty="s") DanielPlot(tensileStrength.lm, main = "Daniel Plot") LenthPlot(tensileStrength.lm, main = "Lenth's Plot")### Tensile Strength Experiment. Taguchi and Wu. 1980 library(BsMD) # Data data(BM86.data,package="BsMD") # Design matrix and responses print(BM86.data) # from Box and Meyer (1986) # Model Fitting. Box and Meyer (1986) example 2. tensileStrength.lm <- lm(y2 ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 + X11 + X12 + X13 + X14 + X15, data = BM86.data) print(coef(tensileStrength.lm)) # Model coefficients par(mfrow=c(1,2),pty="s") DanielPlot(tensileStrength.lm, main = "Daniel Plot") LenthPlot(tensileStrength.lm, main = "Lenth's Plot")
Best follow-up experiments based on the MD criterion are suggested to discriminate between competing models.
MD(X, y, nFac, nBlk = 0, mInt = 3, g = 2, nMod, p, s2, nf, facs, nFDes = 4, Xcand, mIter = 20, nStart = 5, startDes = NULL, top = 20, eps = 1e-05)MD(X, y, nFac, nBlk = 0, mInt = 3, g = 2, nMod, p, s2, nf, facs, nFDes = 4, Xcand, mIter = 20, nStart = 5, startDes = NULL, top = 20, eps = 1e-05)
X |
matrix. Design matrix of the initial experiment. |
y |
vector. Response vector of the initial experiment. |
nFac |
integer. Number of factors in the initial experiment. |
nBlk |
integer >=1. The number of blocking factors in the initial experiment.
They are accommodated in the first columns of matrix |
mInt |
integer. Maximum order of the interactions in the models. |
g |
vector. Variance inflation factor for main effects ( |
nMod |
integer. Number of competing models. |
p |
vector. Posterior probabilities of the competing models. |
s2 |
vector. Competing model variances. |
nf |
vector. Factors considered in each of the models. |
facs |
matrix. Matrix [ |
nFDes |
integer. Number of runs to consider in the follow-up experiment. |
Xcand |
matrix. Candidate runs to be chosen for the follow-up design. |
mIter |
integer. If 0, then user-entered designs |
nStart |
integer. Number of starting designs. |
startDes |
matrix. Matrix |
top |
integer. Highest MD follow-up designs recorded. |
eps |
numeric. A small number (1e-5 by default) used for computations. |
The MD criterion, proposed by Meyer, Steinberg and Box is used to discriminate
among competing models. Random starting runs chosen from Xcand are used
for the Wynn search of best MD follow-up designs. nStart starting points are
tried in the search limited to mIter iterations. If mIter=0 then
startDes user-provided designs are used. Posterior probabilities and
variances of the competing models are obtained from BsProb.
The function calls the FORTRAN subroutine ‘md’ and captures
summary results. The complete output of the FORTRAN code is save in
the ‘MDPrint.out’ file in the working directory.
A list with all input and output parameters of the FORTRAN
subroutine MD. Most of the variable names kept to match FORTRAN code.
NSTART |
Number of starting designs. |
NRUNS |
Number of runs used in follow-up designs. |
ITMAX |
Maximum number of iterations for each Wynn search. |
INITDES |
Number of starting points. |
NO |
Numbers of runs already completed before follow-up. |
IND |
Indicator; 0 indicates the user supplied starting designs. |
X |
Matrix for initial data ( |
Y |
Response values from initial experiment ( |
GAMMA |
Variance inflation factor. |
GAM2 |
If |
BL |
Number of blocks (>=1) accommodated in first columns of |
.
COLS |
Number of factors. |
N |
Number of candidate runs. |
Xcand |
Matrix of candidate runs. ( |
NM |
Number of models considered. |
P |
Models posterior probability. |
SIGMA2 |
Models variances. |
NF |
Number of factors per model. |
MNF |
Maximum number of factor in models. ( |
JFAC |
Matrix with the factor numbers for each of the models. |
CUT |
Maximum interaction order considered. |
MBEST |
If |
NTOP |
Number of the top best designs. |
TOPD |
The D value for the best |
TOPDES |
Top |
ESP |
"Small number" provided to the ‘md’ FORTRAN subroutine. 1e-5 by default. |
flag |
Indicator = 1, if the ‘md’ subroutine finished properly, -1 otherwise. |
The function is a wrapper to call the FORTAN subroutine ‘md’, modification of Daniel Meyer's original program, ‘MD.f’, part of the mdopt bundle for Bayesian model discrimination of multifactor experiments.
R. Daniel Meyer. Adapted for R by Ernesto Barrios.
Meyer, R. D., Steinberg, D. M. and Box, G. E. P. (1996). "Follow-Up Designs to Resolve Confounding in Multifactor Experiments (with discussion)". Technometrics, Vol. 38, No. 4, pp. 303–332.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
### Injection Molding Experiment. Meyer et al. 1996, example 2. library(BsMD) data(BM93.e3.data,package="BsMD") X <- as.matrix(BM93.e3.data[1:16,c(1,2,4,6,9)]) y <- BM93.e3.data[1:16,10] p <- c(0.2356,0.2356,0.2356,0.2356,0.0566) s2 <- c(0.5815,0.5815,0.5815,0.5815,0.4412) nf <- c(3,3,3,3,4) facs <- matrix(c(2,1,1,1,1,3,3,2,2,2,4,4,3,4,3,0,0,0,0,4),nrow=5, dimnames=list(1:5,c("f1","f2","f3","f4"))) nFDes <- 4 Xcand <- matrix(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1, -1,-1,-1,-1,1,1,1,1,-1,-1,-1,-1,1,1,1,1, -1,-1,1,1,-1,-1,1,1,-1,-1,1,1,-1,-1,1,1, -1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1, -1,1,1,-1,1,-1,-1,1,1,-1,-1,1,-1,1,1,-1), nrow=16,dimnames=list(1:16,c("blk","f1","f2","f3","f4")) ) injectionMolding.MD <- MD(X = X, y = y, nFac = 4, nBlk = 1, mInt = 3, g = 2, nMod = 5, p = p, s2 = s2, nf = nf, facs = facs, nFDes = 4, Xcand = Xcand, mIter = 20, nStart = 25, top = 10) summary(injectionMolding.MD) ### Reactor Experiment. Meyer et al. 1996, example 3. par(mfrow=c(1,2),pty="s") data(Reactor.data,package="BsMD") # Posterior probabilities based on first 8 runs X <- as.matrix(cbind(blk = rep(-1,8), Reactor.data[c(25,2,19,12,13,22,7,32), 1:5])) y <- Reactor.data[c(25,2,19,12,13,22,7,32), 6] reactor8.BsProb <- BsProb(X = X, y = y, blk = 1, mFac = 5, mInt = 3, p =0.25, g =0.40, ng = 1, nMod = 32) plot(reactor8.BsProb,prt=TRUE,,main="(8 runs)") # MD optimal 4-run design p <- reactor8.BsProb$ptop s2 <- reactor8.BsProb$sigtop nf <- reactor8.BsProb$nftop facs <- reactor8.BsProb$jtop nFDes <- 4 Xcand <- as.matrix(cbind(blk = rep(+1,32), Reactor.data[,1:5])) reactor.MD <- MD(X = X, y = y, nFac = 5, nBlk = 1, mInt = 3, g =0.40, nMod = 32, p = p,s2 = s2, nf = nf, facs = facs, nFDes = 4, Xcand = Xcand, mIter = 20, nStart = 25, top = 5) summary(reactor.MD) # Posterior probabilities based on all 12 runs X <- rbind(X, Xcand[c(4,10,11,26), ]) y <- c(y, Reactor.data[c(4,10,11,26),6]) reactor12.BsProb <- BsProb(X = X, y = y, blk = 1, mFac = 5, mInt = 3, p = 0.25, g =1.20,ng = 1, nMod = 5) plot(reactor12.BsProb,prt=TRUE,main="(12 runs)")### Injection Molding Experiment. Meyer et al. 1996, example 2. library(BsMD) data(BM93.e3.data,package="BsMD") X <- as.matrix(BM93.e3.data[1:16,c(1,2,4,6,9)]) y <- BM93.e3.data[1:16,10] p <- c(0.2356,0.2356,0.2356,0.2356,0.0566) s2 <- c(0.5815,0.5815,0.5815,0.5815,0.4412) nf <- c(3,3,3,3,4) facs <- matrix(c(2,1,1,1,1,3,3,2,2,2,4,4,3,4,3,0,0,0,0,4),nrow=5, dimnames=list(1:5,c("f1","f2","f3","f4"))) nFDes <- 4 Xcand <- matrix(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1, -1,-1,-1,-1,1,1,1,1,-1,-1,-1,-1,1,1,1,1, -1,-1,1,1,-1,-1,1,1,-1,-1,1,1,-1,-1,1,1, -1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1, -1,1,1,-1,1,-1,-1,1,1,-1,-1,1,-1,1,1,-1), nrow=16,dimnames=list(1:16,c("blk","f1","f2","f3","f4")) ) injectionMolding.MD <- MD(X = X, y = y, nFac = 4, nBlk = 1, mInt = 3, g = 2, nMod = 5, p = p, s2 = s2, nf = nf, facs = facs, nFDes = 4, Xcand = Xcand, mIter = 20, nStart = 25, top = 10) summary(injectionMolding.MD) ### Reactor Experiment. Meyer et al. 1996, example 3. par(mfrow=c(1,2),pty="s") data(Reactor.data,package="BsMD") # Posterior probabilities based on first 8 runs X <- as.matrix(cbind(blk = rep(-1,8), Reactor.data[c(25,2,19,12,13,22,7,32), 1:5])) y <- Reactor.data[c(25,2,19,12,13,22,7,32), 6] reactor8.BsProb <- BsProb(X = X, y = y, blk = 1, mFac = 5, mInt = 3, p =0.25, g =0.40, ng = 1, nMod = 32) plot(reactor8.BsProb,prt=TRUE,,main="(8 runs)") # MD optimal 4-run design p <- reactor8.BsProb$ptop s2 <- reactor8.BsProb$sigtop nf <- reactor8.BsProb$nftop facs <- reactor8.BsProb$jtop nFDes <- 4 Xcand <- as.matrix(cbind(blk = rep(+1,32), Reactor.data[,1:5])) reactor.MD <- MD(X = X, y = y, nFac = 5, nBlk = 1, mInt = 3, g =0.40, nMod = 32, p = p,s2 = s2, nf = nf, facs = facs, nFDes = 4, Xcand = Xcand, mIter = 20, nStart = 25, top = 5) summary(reactor.MD) # Posterior probabilities based on all 12 runs X <- rbind(X, Xcand[c(4,10,11,26), ]) y <- c(y, Reactor.data[c(4,10,11,26),6]) reactor12.BsProb <- BsProb(X = X, y = y, blk = 1, mFac = 5, mInt = 3, p = 0.25, g =1.20,ng = 1, nMod = 5) plot(reactor12.BsProb,prt=TRUE,main="(12 runs)")
12-run Plackett-Burman design matrix.
data(PB12Des)data(PB12Des)
A data frame with 12 observations on the following 11 variables.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
numeric vectors. Contrast factor.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (2004). Statistics for Experimenters II. Wiley.
library(BsMD) data(PB12Des,package="BsMD") str(PB12Des) X <- as.matrix(PB12Des) print(t(X)%*%X)library(BsMD) data(PB12Des,package="BsMD") str(PB12Des) X <- as.matrix(PB12Des) print(t(X)%*%X)
Method function for plotting marginal factor posterior probabilities for Bayesian screening.
## S3 method for class 'BsProb' plot(x, code = TRUE, prt = FALSE, cex.axis=par("cex.axis"), ...)## S3 method for class 'BsProb' plot(x, code = TRUE, prt = FALSE, cex.axis=par("cex.axis"), ...)
x |
list. List of class |
code |
logical. If |
prt |
logical. If |
cex.axis |
Magnification used for the axis annotation.
See |
... |
additional graphical parameters passed to |
A spike plot, similar to barplots, is produced with a spike for each factor.
Marginal posterior probabilities are used for the vertical axis.
If code=TRUE, X1, X2, ... are used to label the factors
otherwise the original factor names are used.
If prt=TRUE, the print.BsProb function is called
and the posterior probabilities are displayed.
When BsProb is called for more than one value of gamma (g),
the spikes for each factor probability are overlapped to show the
resulting range of each marginal probability.
The function is called for its side effects. It returns an invisible
NULL.
Ernesto Barrios.
Box, G. E. P and R. D. Meyer (1986). "An Analysis for Unreplicated Fractional Factorials". Technometrics. Vol. 28. No. 1. pp. 11–18.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
BsProb, print.BsProb, summary.BsProb.
library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) plot(drillAdvance.BsProb) summary(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) plot(drillAdvance.BsProbG, code = FALSE, prt = TRUE)library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) plot(drillAdvance.BsProb) summary(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) plot(drillAdvance.BsProbG, code = FALSE, prt = TRUE)
Printing method for lists of class BsProb. Prints the posterior
probabilities of factors and models from the Bayesian screening procedure.
## S3 method for class 'BsProb' print(x, X = TRUE, resp = TRUE, factors = TRUE, models = TRUE, nMod = 10, digits = 3, plt = FALSE, verbose = FALSE, ...)## S3 method for class 'BsProb' print(x, X = TRUE, resp = TRUE, factors = TRUE, models = TRUE, nMod = 10, digits = 3, plt = FALSE, verbose = FALSE, ...)
x |
list. Object of |
X |
logical. If |
resp |
logical. If |
factors |
logical. Marginal posterior probabilities are printed if |
models |
logical. If |
nMod |
integer. Number of the top ranked models to print. |
digits |
integer. Significant digits to use for printing. |
plt |
logical. Factor marginal probabilities are plotted if |
verbose |
logical. If |
... |
additional arguments passed to |
The function prints out marginal factors and models posterior probabilities. Returns invisible list with the components:
calc |
numeric vector with general calculation information. |
probabilities |
Data frame with the marginal posterior factor probabilities. |
models |
Data frame with model the posterior probabilities. |
Ernesto Barrios.
Box, G. E. P and R. D. Meyer (1986). "An Analysis for Unreplicated Fractional Factorials". Technometrics. Vol. 28. No. 1. pp. 11–18.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
BsProb, summary.BsProb, plot.BsProb.
library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) print(drillAdvance.BsProb) plot(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) print(drillAdvance.BsProbG, X = FALSE, resp = FALSE) plot(drillAdvance.BsProbG)library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) print(drillAdvance.BsProb) plot(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) print(drillAdvance.BsProbG, X = FALSE, resp = FALSE) plot(drillAdvance.BsProbG)
Printing method for lists of class MD. Displays the
best MD criterion set of runs and their MD for follow-up experiments.
## S3 method for class 'MD' print(x, X = FALSE, resp = FALSE, Xcand = TRUE, models = TRUE, nMod = x$nMod, digits = 3, verbose=FALSE, ...)## S3 method for class 'MD' print(x, X = FALSE, resp = FALSE, Xcand = TRUE, models = TRUE, nMod = x$nMod, digits = 3, verbose=FALSE, ...)
x |
list of class |
X |
logical. If |
resp |
logical If |
Xcand |
logical. Prints the candidate runs if |
models |
logical. Competing models are printed if |
nMod |
integer. Top models to print. |
digits |
integer. Significant digits to use in the print out. |
verbose |
logical. If |
... |
additional arguments passed to |
The function is mainly called for its side effects. Prints out the selected
components of the class MD objects, output of the MD function.
For example the marginal factors and models posterior probabilities and
the top MD follow-up experiments with their corresponding MD statistic.
It returns invisible list with the components:
calc |
Numeric vector with basic calculation information. |
models |
Data frame with the competing models posterior probabilities. |
follow-up |
Data frame with the runs for follow-up experiments and their corresponding MD statistic. |
Ernesto Barrios.
Meyer, R. D., Steinberg, D. M. and Box, G. E. P. (1996). "Follow-Up Designs to Resolve Confounding in Multifactor Experiments (with discussion)". Technometrics, Vol. 38, No. 4, pp. 303–332.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
# Injection Molding Experiment. Meyer et al. 1996. Example 2. # MD for one extra experiment. library(BsMD) data(BM93.e3.data,package="BsMD") X <- as.matrix(BM93.e3.data[1:16,c(1,2,4,6,9)]) y <- BM93.e3.data[1:16,10] nBlk <- 1 nFac <- 4 mInt <- 3 g <- 2 nMod <- 5 p <- c(0.2356,0.2356,0.2356,0.2356,0.0566) s2 <- c(0.5815,0.5815,0.5815,0.5815,0.4412) nf <- c(3,3,3,3,4) facs <- matrix(c(2,1,1,1,1,3,3,2,2,2,4,4,3,4,3,0,0,0,0,4),nrow=5, dimnames=list(1:5,c("f1","f2","f3","f4"))) nFDes <- 1 Xcand <- matrix(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1, -1,-1,-1,-1,1,1,1,1,-1,-1,-1,-1,1,1,1,1, -1,-1,1,1,-1,-1,1,1,-1,-1,1,1,-1,-1,1,1, -1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1, -1,1,1,-1,1,-1,-1,1,1,-1,-1,1,-1,1,1,-1), nrow=16,dimnames=list(1:16,c("blk","f1","f2","f3","f4")) ) mIter <- 0 startDes <- matrix(c(9,11,12,15),nrow=4) top <- 10 injectionMolding.MD <- MD(X=X,y=y,nFac=nFac,nBlk=nBlk,mInt=mInt,g=g, nMod=nMod,p=p,s2=s2,nf=nf,facs=facs, nFDes=nFDes,Xcand=Xcand,mIter=mIter,startDes=startDes,top=top) print(injectionMolding.MD) summary(injectionMolding.MD)# Injection Molding Experiment. Meyer et al. 1996. Example 2. # MD for one extra experiment. library(BsMD) data(BM93.e3.data,package="BsMD") X <- as.matrix(BM93.e3.data[1:16,c(1,2,4,6,9)]) y <- BM93.e3.data[1:16,10] nBlk <- 1 nFac <- 4 mInt <- 3 g <- 2 nMod <- 5 p <- c(0.2356,0.2356,0.2356,0.2356,0.0566) s2 <- c(0.5815,0.5815,0.5815,0.5815,0.4412) nf <- c(3,3,3,3,4) facs <- matrix(c(2,1,1,1,1,3,3,2,2,2,4,4,3,4,3,0,0,0,0,4),nrow=5, dimnames=list(1:5,c("f1","f2","f3","f4"))) nFDes <- 1 Xcand <- matrix(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1, -1,-1,-1,-1,1,1,1,1,-1,-1,-1,-1,1,1,1,1, -1,-1,1,1,-1,-1,1,1,-1,-1,1,1,-1,-1,1,1, -1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1,-1,1, -1,1,1,-1,1,-1,-1,1,1,-1,-1,1,-1,1,1,-1), nrow=16,dimnames=list(1:16,c("blk","f1","f2","f3","f4")) ) mIter <- 0 startDes <- matrix(c(9,11,12,15),nrow=4) top <- 10 injectionMolding.MD <- MD(X=X,y=y,nFac=nFac,nBlk=nBlk,mInt=mInt,g=g, nMod=nMod,p=p,s2=s2,nf=nf,facs=facs, nFDes=nFDes,Xcand=Xcand,mIter=mIter,startDes=startDes,top=top) print(injectionMolding.MD) summary(injectionMolding.MD)
Data of the Reactor Experiment from Box, Hunter and Hunter (1978).
data(Reactor.data)data(Reactor.data)
A data frame with 32 observations on the following 6 variables.
numeric vector. Feed rate factor.
numeric vector. Catalyst factor.
numeric vector. Agitation rate factor.
numeric vector. Temperature factor.
numeric vector. Concentration factor.
numeric vector. Percentage reacted response.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (2004). Statistics for Experimenters II. Wiley.
Box G. E. P, Hunter, W. C. and Hunter, J. S. (1978). Statistics for Experimenters. Wiley.
library(BsMD) data(Reactor.data,package="BsMD") print(Reactor.data)library(BsMD) data(Reactor.data,package="BsMD") print(Reactor.data)
Reduced printing method for class BsProb lists. Prints
posterior probabilities of factors and models from Bayesian screening
procedure.
## S3 method for class 'BsProb' summary(object, nMod = 10, digits = 3, ...)## S3 method for class 'BsProb' summary(object, nMod = 10, digits = 3, ...)
object |
list. |
nMod |
integer. Number of the top ranked models to print. |
digits |
integer. Significant digits to use. |
... |
additional arguments passed to |
The function prints out the marginal factors and models posterior probabilities. Returns invisible list with the components:
calc |
Numeric vector with basic calculation information. |
probabilities |
Data frame with the marginal posterior factor probabilities. |
models |
Data frame with the models posterior probabilities. |
Ernesto Barrios.
Box, G. E. P and R. D. Meyer (1986). "An Analysis for Unreplicated Fractional Factorials". Technometrics. Vol. 28. No. 1. pp. 11–18.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
BsProb, print.BsProb, plot.BsProb.
library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) plot(drillAdvance.BsProb) summary(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) plot(drillAdvance.BsProbG) summary(drillAdvance.BsProbG)library(BsMD) data(BM86.data,package="BsMD") X <- as.matrix(BM86.data[,1:15]) y <- BM86.data["y1"] # Using prior probability of p = 0.20, and k = 10 (gamma = 2.49) drillAdvance.BsProb <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.20, g = 2.49, ng = 1, nMod = 10) plot(drillAdvance.BsProb) summary(drillAdvance.BsProb) # Using prior probability of p = 0.20, and a 5 <= k <= 15 (1.22 <= gamma <= 3.74) drillAdvance.BsProbG <- BsProb(X = X, y = y, blk = 0, mFac = 15, mInt = 1, p = 0.25, g = c(1.22, 3.74), ng = 3, nMod = 10) plot(drillAdvance.BsProbG) summary(drillAdvance.BsProbG)
Reduced printing method for lists of class MD. Displays the
best MD criterion set of runs and their MD for follow-up experiments.
## S3 method for class 'MD' summary(object, digits = 3, verbose=FALSE, ...)## S3 method for class 'MD' summary(object, digits = 3, verbose=FALSE, ...)
object |
list of |
digits |
integer. Significant digits to use in the print out. |
verbose |
logical. If |
... |
additional arguments passed to |
It prints out the marginal factors and models posterior probabilities and the top MD follow-up experiments with their corresponding MD statistic.
Ernesto Barrios.
Meyer, R. D., Steinberg, D. M. and Box, G. E. P. (1996). "Follow-Up Designs to Resolve Confounding in Multifactor Experiments (with discussion)". Technometrics, Vol. 38, No. 4, pp. 303–332.
Box, G. E. P and R. D. Meyer (1993). "Finding the Active Factors in Fractionated Screening Experiments". Journal of Quality Technology. Vol. 25. No. 2. pp. 94–105.
### Reactor Experiment. Meyer et al. 1996, example 3. library(BsMD) data(Reactor.data,package="BsMD") # Posterior probabilities based on first 8 runs X <- as.matrix(cbind(blk = rep(-1,8), Reactor.data[c(25,2,19,12,13,22,7,32), 1:5])) y <- Reactor.data[c(25,2,19,12,13,22,7,32), 6] reactor.BsProb <- BsProb(X = X, y = y, blk = 1, mFac = 5, mInt = 3, p =0.25, g =0.40, ng = 1, nMod = 32) # MD optimal 4-run design p <- reactor.BsProb$ptop s2 <- reactor.BsProb$sigtop nf <- reactor.BsProb$nftop facs <- reactor.BsProb$jtop nFDes <- 4 Xcand <- as.matrix(cbind(blk = rep(+1,32), Reactor.data[,1:5])) reactor.MD <- MD(X = X, y = y, nFac = 5, nBlk = 1, mInt = 3, g =0.40, nMod = 32, p = p,s2 = s2, nf = nf, facs = facs, nFDes = 4, Xcand = Xcand, mIter = 20, nStart = 25, top = 5) print(reactor.MD) summary(reactor.MD)### Reactor Experiment. Meyer et al. 1996, example 3. library(BsMD) data(Reactor.data,package="BsMD") # Posterior probabilities based on first 8 runs X <- as.matrix(cbind(blk = rep(-1,8), Reactor.data[c(25,2,19,12,13,22,7,32), 1:5])) y <- Reactor.data[c(25,2,19,12,13,22,7,32), 6] reactor.BsProb <- BsProb(X = X, y = y, blk = 1, mFac = 5, mInt = 3, p =0.25, g =0.40, ng = 1, nMod = 32) # MD optimal 4-run design p <- reactor.BsProb$ptop s2 <- reactor.BsProb$sigtop nf <- reactor.BsProb$nftop facs <- reactor.BsProb$jtop nFDes <- 4 Xcand <- as.matrix(cbind(blk = rep(+1,32), Reactor.data[,1:5])) reactor.MD <- MD(X = X, y = y, nFac = 5, nBlk = 1, mInt = 3, g =0.40, nMod = 32, p = p,s2 = s2, nf = nf, facs = facs, nFDes = 4, Xcand = Xcand, mIter = 20, nStart = 25, top = 5) print(reactor.MD) summary(reactor.MD)