The unrepx package provides much of what one needs to analyze basic unreplicated screening experiments such as two-level factorial and fractional factorial designs, in which one has a set of independent effect estimates but no degrees of freedom for error. The analyses are based on these underlying assumptions about the effect estimates:
The first and last of these are the most important. There are several cases in which effects are not independent – among them Plackett-Burman designs and definitive screening designs – for which other software is needed. And without the effect-sparsity assumption, one needs replications, or some other assumption, in order to assess which effects are likely to be active.
This document provides an overview and a few examples. The user should consult the help pages for details and finer points.
estimation of effects without also providing a design matrix. Instead,
the response values must be arranged in standard order, whereby the
first factor varies the fastest, the second the second fastest, and so
help("yates") for details on these functions, as
well as literature references.
As an example, consider the pilot-plant example given in the text by Box, Hunter, and Hunter (2nd ed.), page 177. The experiment involves three factors T, C, and K, each having two levels. Here we estimate the effects of these factors on process yield, by providing the response values in standard order (as presented in the text):
## Loading required package: unrepx
( pilotEff = yates(c(60,72,54,68,52,83,45,80), labels = c("T","C","K")) )
## T C TC K TK CK TCK ## 23.0 -5.0 1.5 1.5 10.0 0.0 0.5 ## attr(,"mean") ## ## 64.25
Note there are 7 effects from the 8 experimental runs; the mean is added as an attribute.
gyates() functions can also
reverse themselves and produce the original response values, given the
effects. They decide which to do based on the length of the data
## --- +-- -+- ++- --+ +-+ -++ +++ ## 60 72 54 68 52 83 45 80
unrepx provides three graphical tools for helping to identify active effects: [half]-normal plots, reference plots, and Pareto plots.
Half-normal and normal plots are available via
(see details in its help page). Here is a half-normal plot of the
effects we just computed, with labels added for all those that exceed
hnplot(pilotEff, ID = 10) title("Half-normal plot of the pilot-plant effects")
By default, a reference line is added (see the section below on statistics for details). The dots are colored blue for positive and red for negative effects. Visually, it appears that the T and TK effects are the most active, as they appear to be outliers. Ordinary normal plots are possible too, but we strongly discourage using them as they can be misleading.
A reference plot is a dotplot of the effects, along with a reference curve. By default, the absolute effects are plotted along with a half-normal curve:
refplot(pilotEff, ID = 10) title("Reference plot of the pilot-plant effects")
As a by-product, a general-purpose
dot.plot function is
also provided. These plots are resizeable in both directions: with a
vertical resize, the dots are re-stacked with the correct spacing; and
with a horizontal resize, the dots may be re-binned.
Pareto plots are produced using
parplot(). The plot is a
bar graph with the bars ordered by the size of the effect. In addition,
cutoff line(s) are added for the margin of error (ME) and the
simultaneous margin of error (SME) of the effects (see the section below
parplot(pilotEff) title("Pareto plot of the pilot-plant effects")
You may have noticed that all of the plots we have seen are enhanced
by lines or curves that help to identify active effects. This is done
via a pseudo-standard-error estimate. The effect-sparsity and
equal-variances assumptions suggests that an estimate of the standard
error of the effects may be obtained by looking at the variability of
the smaller effects. The
PSE() function calculates these
## Zahn_PSE ## 4.011005
PSE(pilotEff, method = "Lenth")
## Lenth_PSE ## 2.25
help("PSE") for a list of the available methods. The
default Zahn method is based on the slope of a line through the origin
fitted to the smallest 2/3 of the points in the half-normal plot.
For any PSE estimate, one can compute a t-like statistic via
dividing the effect estimate by the PSE. In turn, there is an associated
reference distribution that can be estimated by simulating the
t statistics for independent inactive effects. (It can be
accessed if you like via the
ref.dist() function; see its
help page. Ordinarily, though, the needed simulations are done – and
stored – in the background as needed.) A high percentile (typically the
95th) of the reference distribution of |t| provides a critical
value. This critical value, multiplied by the PSE, gives a margin of
error for assessing the effects. This value is used in obtaining
the cutoff line labeled “ME” in the Pareto plot.
If one instead uses instead the distribution of max|t| from each simulated set of effects, the simultaneous margin of error (SME) results. The idea here is that if it happens that all of the effects are inactive, it is unlikely for any of their |t| ratios to exceed the SME. This, if we observe some effects that do, it is pretty convincing evidence that we do indeed have active effects in the mix. And that is the case with our pilot-plant example, as you can see in the Pareto plot above.
A non-graphical analysis showing all of these observed t
ratios, along with P values obtained as tail areas of the
reference distribution, is available via
## effect Zahn_PSE t.ratio p.value simult.pval ## T 23.0 4.011005 5.734 0.0010 0.0065 ## TK 10.0 4.011005 2.493 0.0242 0.1459 ## C -5.0 4.011005 -1.247 0.1669 0.7762 ## K 1.5 4.011005 0.374 0.7380 1.0000 ## TC 1.5 4.011005 0.374 0.7380 1.0000 ## TCK 0.5 4.011005 0.125 0.9129 1.0000 ## CK 0.0 4.011005 0.000 1.0000 1.0000
Several sets of example effects are provided in the package:
help("bikeEff") for details).
Finally, if there is anyone wishing to construct half-normal plots by
hand, as originally done in Daniel’s 1959 paper, the
daniel.paper() function provides suitable graph paper for